Non-Darcy Flow Experiments of Water Seepage through Rough-Walled Rock Fractures

. The knowledge of flow phenomena in fractured rocks is very important for groundwater-resources management in hydrogeological engineering. The most commonly used tool to approximate the non-Darcy behavior of the flow velocity is the well-known Forchheimer equation, deploying the “inertial” coefficient 𝛽 that can be estimated experimentally. Unfortunately, the factor of roughness is imperfectly considered in the literature. In order to do this, we designed and manufactured a seepage apparatus that can provide different roughness and aperture in the test; the rough fracture surface is established combining JRC and 3D printing technology. A series of hydraulic tests covering various flows were performed. Experimental data suggest that Forchheimer coefficients are to some extent affected by roughness and aperture. At last, favorable semiempirical Forchheimer equation which can consider fracture aperture and roughness was firstly derived. It is believed that such studies will be quite useful in identifying the limits of applicability of the well-known “cubic law,” in further improving theoretical/numerical models associated with fluid flow through a rough fracture.


Introduction
Understanding uid ow behaviors in fractured rock aquifers is of great concern in numerous industrial and scienti c elds, such as water resources management, sustainable urban drainage, contaminant pollution control, hazardous wastes isolation, and geothermal uid or hydrocarbon exploitation.
It has been observed that uid ow behavior through a rock mass depends on the geological-structural and geological-technical (such as degree of fracturing, orientation, persistence, weathering, moisture conditions and seepage aperture and lling, and roughness) [ ]. As we all known, quantity and proximity of fractures pervading a rock mass de ne the degree of fracturing of the rock and, as a consequence, its permeability. Evidently, the higher the persistence of discontinuities, the stronger the permeability. e weathering process determines an increase of porosity. When it comes to a single fracture, aperture and roughness will be the dominant factors. An adequate knowledge of uid ow through a rough-walled fracture is a starting point for a better interpretation of uid ow and solute transport in fractured rock aquifers.
A wide number of groundwater ow models and numerical tools were developed based on the assumption of Darcy's ow through individual fractures. It has been well recognized, however, that the linear Darcy's law is not always adequate to describe the ow behaviors in natural fractures. Typical characteristics of natural rock fractures include rough walls and asperity contact [ , ], and non-Darcy ow may occur as a result of nonnegligible inertial losses. Previous experimental work demonstrated that Darcy's law fails to predict pressure drops in fractures when inertial e ects are relevant before the fully developed turbulence [ -]. In the post-Darcy regimes where inertial e ects are signi cant, two equations are invariably used to describe pressure drop as a function of average velocity: Forchheimer and Ergun equation [ , ]. Other seemingly di erent correlations may be traced back or manipulated to t the basic forms of these two equations [ ].

Geo uids
To account for such inertial losses, a widely used equation referred to as Forchheimer's law was developed. Qian et al. [ ] used well-controlled laboratory experiments to investigate the ow and transport in a fracture under non-Darcy ow conditions, and it is found that the Forchheimer equation ts the experimental V-relationship nearly perfectly, whereas the Darcy equation is inadequate in this respect. Cherubini et al. [ ] investigated the nonlinear ow by analyzing hydraulic tests on an arti cially created fractured rock sample, and the experimental result shows that the relationship between the ow rate and head gradient matches the Forchheimer equation and describes a strong inertial regime well. Javadi et al. [ ] performed both laminar and turbulent ow simulations for a wide range of ow rates in an arti cial three-dimensional fracture.
ey developed a new geometrical model for nonlinear uid ow through rough fractures, which suggested a polynomial expression, like the Forchheimer law, to describe the dependence of pressure drop on discharge. Considering that the mechanisms of non-Darcy two-phase ows in fractures are not well understood and no general model has been presented to describe them, relying on the full cubic law, Radilla et al. [ ] carried out a series of ow experiments and presented a model to solve this preceding issue.
To incorporate the Forchheimer's law in the analytical or numerical solutions for the simulation of non-Darcy ow, the determination of the phenomenological coe cients in Forchheimer equation becomes indispensable. For many practical problems, however, it may be di cult to directly determine the coe cients through experimental tests. In this circumstance, parametric expressions have to be used instead. For porous media, numerous theoretical and empirical expressions of coe cients have been developed, which are typically characterized as a function of the particle diameter and the medium porosity [ ]. For uid ow in rock fractures, however, the study on the parametric expressions for these two coe cients and their physical background is yet insu cient and rarely reported in the literature.
Roughness has a large in uence on uid ow and transport through tight, rough-walled fractures where non-Darcy o wi sp a r t i c u l a r l ye a s yt oo c c u r[ ] .Q i a ne ta l .
[ , ] created the surface roughness by gluing small Plexiglas plates too neo fthetwoverticalwalls,k eep ingtheseco ndvertical wall smooth. Zhang and Nemcik [ ] investigated the ow behavior in fractures, which were created in the laboratory by splitting initially intact samples of a ne grain sandstone block into halves.
More recently, the developed three-dimensional ow models were used to simulate uid ow through various random synthetic rough-walled fractures. e rough-walled fractures were created by combining random elds of aperture and the mean wall topography or midsurface, which quanti es undulation about the fracture plane. Javadi et al.
[ ] performed both laminar and turbulent ow simulations for a wide range of ow rates in an arti cial three-dimensional fracture. e primary motivation of the present study was to experimentally evaluate the Forchheimer equation coecients for non-Darcy ow in Forchheimer equation in roughwalled fractures, where e ects of fracture roughness, fracture aperture, and ow regime were considered and e ects of orientation and fracture networks are not, as the focus is not on bulk ow, but rather on the properties of individual discrete fractures. A series of hydraulic tests were performed under di erent roughness and fractures apertures, with the fracture formed by two planes and the ow taking place in the horizontal direction. e aperture varied from to mm, whereas the "fracture roughness" was created by combining JRC and D printing technology, where the Joint Roughness Coe cient (JRC) means a dimensionless measure of fracture surface roughness ranging from to ( Figure ). Based on the experimental observations, a semiempirical equation of the Forchheimer's coe cients and dependent on hydraulic aperture and roughness was proposed. where and are the Forchheimer coe cients describing pressure losses due to viscous and inertial dissipation mechanisms, respectively [ , ]. Where = / and = , = 2 /12 ( is the aperture of the idealized parallel smooth fracture; is the uid density) is the intrinsic permeability, and is the non-Darcy coe cient or inertial resistance coe cient dependent on the geometrical properties of the medium [ ]. Using dimensional analysis, Schrauf and Evans [ ] rewrote ( ) as

Theories and Experimental Methodology
where and are dimensionless coe cients. By comparing ( ) to ( ), one can obtain =12and = .
. . Sample Preparation. In this study, fracture specimens with di erent joint roughness coe cient were fabricated by D printing technology according to the JRC standard pro le curve proposed by Barton and Choubey [ ], and a high-speed non-Darcy seepage experiment system is established based on this.
e Darcy and non-Darcy seepage experiments of rough fracture with various fracture aperture are carried out in this system. e method of fracture generation is presented in the following three steps.
e JRC pro le curve proposed by Barton is scanned and plotted according to the original size. A er obtaining the digital JRC standard curve, the d fracture model is built by D modeling of CAD so ware. Considering the seal ability and convenience of adjusting the crack width, mm cushion layer is set on both sides of the rough fracture, and mm ow transition zone is added in the inlet and outlet of the smooth plate. e D CAD model is transformed into STL format le and three-dimensional fracture plate model is made based on D printing technology (precision . mm) as shown in Figure . . . . Production of Cement Specimens.
e outer diameter of mm and wall thickness of . mm of a circular PVC pipe were used as a sample mold, with rough fracture plate in the center. Copper pipe was embedded in the mold upstream, central, and downstream for water pressure monitoring. Specimen pouring groove is shown in Figure . e test piece has a diameter of . cm and a total length of cm. e rough fracture portion has a width of cm and a length of cm.
. . . Test Device Assembly. e concrete sample (in Figure (b)) was xed to the test stand, as shown in Figure . During the test, by setting di erent thickness gasket in the smooth edges (cushion area) reserved between the upper and lower plates of the sample to simulate di erent fracture aperture, the fracture aperture was set to , , , , and mm, respectively. e upstream pressure was controlled by a pump and a regulator. e water tank is prepared between the inlet and sample to maintain the stability of the water ow within the specimen and, ultimately, the circulating water supply system which is made up of the storage tank, pumps, test stand, and other equipment ( Figure ). Due to the large ow, delta weir was used as a ow measurement device.
. . Experiment System Test. In order to verify the accuracy of the test equipment and the measurement system, we conduct the fracture ow test under the condition of laminar ow. e results for JRC = ∼ andJRC= ∼ low velocity seepage test in a rough fracture with aperture set as mm is shown in Figure . It can be seen that when the hydraulic gradient is low ( ⩽1 ), there is a slight deviation from Darcy's law between ow and hydraulic gradient at JRC = -. Overall, the ow and hydraulic gradient satisfy a linear relationship, indicating that the fracture ow in this scenario is in the laminar ow stage. e water ow conforms to the linear Darcy's law; the second term of the right side of ( ) can be ignored. e results of the relevant test show the accuracy of the results obtained from another angle.

. . Relation between Flow within Smooth Fracture and Crack
Opening. In order to explore the correlation between velocity, hydraulic gradient, and fracture aperture under high hydraulic gradient, high-speed non-Darcy seepage tests in a smooth fracture are carried out before carrying out rough fracture non-Darcy seepage test.
Figure shows the velocity and hydraulic gradient curve for fracture apertures. e relationship between the velocity and the hydraulic gradient for the large fracture aperture and high velocity can be well described by a power function, and the ow velocity and the hydraulic gradient have deviated from the linear relationship [ -].
rough tting the test data point, the transform Izbash equation ( ) can be used to analyze the results. Izbash's equation describes the relationship between the uid ow velocity and the hydraulic gradient is given a theoretical background arriving from a drag model. e exponent in Izbash's equation [ , ] has been obtained over a wide range of Reynolds' numbers.
where , are the tting coe cients. e result of analysis showed that the correlation coefcient of each curve is more than . ; = . ∼ . , v a l u ei n c r e a s e sa st h ef r a c t u r ea p e r t u r eb e c o m e sl a r g e r . For nonlinear ow yet to approach fully turbulent state, changes between and . , = . f o ra l la p e r t u r e si nt h i s test representing a fully turbulent ow. e owrateisusuallyexpressedintermsofthedimensionless Reynolds number, which quanti es the relative strength of inertia forces as compared to viscous forces. For uid ow through fractures, Reynolds number [ , ] can be de ned as where is the characteristic length (hydraulic radius) of the idealized parallel smooth fracture. velocity of fractures and the hydraulic gradient and fracture aperture is obtained by mathematical derivation.
where is the average velocity of the section, is the Carmen constant, and is the undetermined coe cient. As shown in ( ), by assigning = 0.5,t h ep r o d u c to f the average velocity and the hydraulic gradient 0.5 can be described as a function of the fracture aperture ,assho wn in the following: where 1 , 2 , 3 are the tting coe cient, which can be obtained from the test data.
According to the relationship shown in ( ), a tted curve (as shown in Figure ) could be obtained by sorting out the test data in Figure . e results show that the correlation coe cient is .
As shown in Figure  ef r a c t u r es u r f a c eo f di erent rock cracks is di erent, so whether the correlation of smooth fracture ow can be directly applied to rough fracture is still an open question. erefore, the fracture roughness is introduced, and, based on this, it is of great theoretical value to explore whether there is some relevance in velocity, hydraulic gradient, fracture aperture, and roughness.
In the above-mentioned smooth parallel plate high-speed fracture ow test, the minimum Reynolds number is . Obviously it is in turbulent state. In rough fractures, the hydraulicradiusisclearlylessthanthesmoothfracturedueto the roughness. In order to describe the relationship between ow velocity and hydraulic gradient in complete ow, the Forchheimer equation for Darcy and non-Darcy e ect is still used as benchmark.
Based on the high-speed non-Darcy seepage test, the high-speed non-Darcy seepage tests with fracture aperture of , , , , and mm are carried out to obtain the relationship between the hydraulic gradient and velocity in the fracture. Reynolds numbers are obtained in each case based on ( ). According to the knowledge of the classic pipeline ow, when Reynolds number is greater than , the ow can be considered as turbulence. So the non-Darcy turbulent can be distinguished.
For di erent JRC and fracture apertures, a series of repeated experiments are carried out, respectively, to obtain the ow rate under the corresponding hydraulic gradient. For the same fracture aperture, the greater the JRC is, the greater hydraulic gradient it will require and the stronger the degree of deviation from the linear relationship is.
. . . Analysis Based on Power Function. Basedontheleastsquares method, the tting curve is obtained according to the power function of ( ). e correlation coe cient is . ∼ .
. e tting coe cient and correlation coe cient are summarized in Table . e exponent coe cient is between . and . . Data presented in Table are the results for a mm aperture fracture. e ow state includes the ow from the transition zone to the turbulence zone.
erefore, under di erent roughness conditions, value cannot be uni ed . . A er several ttings, when the hydraulic gradient index = . ∼ . , the correlation coe cient can be maintained above . .

Geo uids
Modeled on ( ), the relationship between / (the ratio of ow rate to hydraulic gradient )a ndJR Cwasp lo t t ed in Figure . It can be observed that there is a negative power function relationship between / and JRC. e fracture width of the test specimen was . m. According to the ow rate calculation equation, the empirical relationship between the unit discharge ,t h eh y d r a u l i c gradient , and the rock joint roughness coe cient JRC can be obtained from the ow calculation equation: where is unit discharge (m 3 /s), is the hydraulic gradient, JRC is the rock joint roughness coe cient, respectively, and , , are the empirical coe cients obtained from test. Empirical coe cient = .
∼ . , = . ∼ . , and = . ∼ . . According to ( ), ( ) is used to establish the empirical relationship between ow rate, hydraulic gradient, and roughness. However, due to the change of hydraulic gradient during the test process, the ow state in the fracture is transformed from transition state to turbulent state. e hydraulic gradient index cannot be used as a certain value. Although the empirical coe cient does not change much, it is impossible to t the equation by using a certain value.

. . . Forchheimer Equation Analysis
. In order to further study seepage in rough fractures under high velocity conditions, the Forchheimer equation is introduced to analyze the experimental data. e results obtained by tting the Forchheimer equation aresummarizedinT able . ecorrelationcoe cient 2 of the experimental data which is tted based on the Forchheimer equation is more than . , and the tting e ect is better. It can be seen from Figure that, as the ow velocity increases, the in uence of inertial force becomes more and more obvious, and the relationship between hydraulic gradient and seepage velocity increasingly deviates from the linear relationship.
In order to analyze the e ect of the fracture aperture on the seepage movement, the relationship between the monomial coe cient and the quadratic coe cient of the ForchheimerequationissummarizedinTable . efracture aperture is plotted, respectively, as shown in Figure .  the fracture ow, respectively. In order to further analyze the relationship between the monomial coe cient, the quadratic coe cient of the Forchheimer equation, and the fracture roughness, the relationship between the monomial coe cient andtheJRCcurveisplottedinFigure .
According to the basic characteristics of ( ), the monomial coe cient is related to the reciprocal of the square of the fracture aperture. Since the in uence of the JRC is taken into account, the actual fracture aperture is changed. erefore, the ordinate is set to × * , in Figures (a) and (b). e a b s c i s s ai ss e tt oJ R C ,s oa st oa n a l y z et h ec o r r e s p o n d i n g relationship between JRC and monomial coe cient and fracture aperture.
As shown in Figures (a) and (b), due to the in uence of roughness, which leads to a changing e ective fracture T : Fitting results of hydraulic gradient and velocity Forchheimer of rough fracture under di erent width and JRC.  of the non-Darcy ow deviating from the linear ow. It can be obtained from ( ) that = .
( ) e non-Darcy in uence coe cient of high-speed non-Darcy uid is ,w h e r e is the quadratic coe cient of Forchheimer obtained by tting the equation and is the acceleration of gravity.
From the point of view of the trend, JRC is linearly related to the quadratic term coe cient for di erent fracture apertures. Under the same fracture aperture condition, with the increase of joint roughness coe cient, the non-Darcy in uence coe cient of rough fracture increases; in the same joint roughness coe cient ssure, the non-Darcy in uence coe cient decreases with the increase of the fracture aperture. It is shown that the extent of the non-Darcy in the high-speed seepage condition is inversely proportional to the fracture aperture, which is proportional to the joint roughness coe cient in the rough rock fracture.
As discussed above, the non-Darcy in uence coe cient d e c r e a s e sw i t ht h ei n c r e a s eo ft h ef r a c t u r ea p e r t u r ef o rt h e same roughness fracture specimen. e unit of the non-Darcy in uence coe cient is m −1 and unit of the fracture aperture is m; it is inferred that the non-Darcy coe cient of the ssure is inversely proportional to the fracture aperture. e curves of the relationship in Figure ( where represents the non-Darcy in uence coe cient of the fracture, and the tting results are summarized in Table . T : e relationship between non-Darcy number and gap width under di erent joint roughness. AscanbeseenfromthedatainT able ,withtheincrease of joint roughness coe cient, coe cient increases. e relationship between and the joint roughness coe cient JRC is plotted in Figure . e relationship between these two can be described as a linear relationship. Based on the leastsquaremethod,the ttingrelationshipis = ⋅ JRC, ( ) where is the linear tting coe cient. e correlation coe cient is . and the linear coe cient is = . ; the empirical relationship between the non-Darcy in uence coe cient and the fracture and the joint roughness coe cient JRC can be obtained as follows: In the equation, is the non-Darcy in uence coe cient, JRC is the joint roughness coe cient of fracture surface, and is the fracture aperture. According to ( ), as the joint roughness increases, the non-Darcy in uence coe cient increases proportionally. It is notable that the focus of this study is on the laboratory scale evaluation of the Forchheimer equation coe cients, which has not yet been attempted to elevate the laboratory results to the eld scale of natural fracturing systems. e ow and hydrodynamic gradient distribution of fractured geological media on a wide scale is a good research topic and has been extensively explored [ , ] and remains to be explored. However, the proposed semiempirical Forchheimer equation in this study also re ects the e ects of aperture and roughnessparametersonhydraulicgradientsand owrates, anditisofimportantreferencevaluefordeterministicfaults or ow headage and crack tip during crack propagation due to engineering disturbances and e ective stress determination, which may provide important guidance for further development or expansion to eld scale ow.

Conclusion
In this study, JRC is used to characterize the joint roughness, and the high-speed non-Darcy seepage test is carried out in fractures of geotechnical soils combining D printing technique. It is mainly used to evaluate the Forchheimer equation coe cients in rough fractures. e empirical relationship between the coe cients A, B of Forchheimer equation and the aperture and roughness of fracture is rst studied, which canbeusedtodescribethehighvelocitynon-Darcy owin rock fractures.
Based on the results of the high-speed non-Darcy seepage test, the obtained pressure gradient and velocity curve show that the quadratic coe cient of the Forchheimer equation can be described as a function of the fracture aperture in the completeturbulencestage.Forchheimer' slawcanbeusedto describe Non-Darcy ow behavior caused by inertia e ect in deformable thick wall fractures.
On the basis of the smooth fracture test, the high-speed non-Darcy seepage was carried out. e in uence of roughness on high-speed non-Darcy ow is analyzed by experiments. e results show that the hydraulic gradient index = . ∼ . ;thecorrelationcoe cientcanbemaintainedat . or more when considering the fracture roughness, according to the power exponent. If the Forchheimer equation is used to t the curve, the correlation coe cient is more than . , and it is clear that the Forchheimer equation has a strong adaptability to the high velocity non-Darcy ow. e in uence of aperture and roughness on Forchheimer coe cient is analyzed. It is observed that the in uence of the fracture aperture on the nonlinear coe cients and is in an accordance with the existing conclusions. e roughness has a signi cant e ect on the nonlinear coe cients and , and the JRC is in linear relationship with the coe cient . However,JRCalsohasane ectonanactualfractureaperture. erefore, the equation requires appropriate correction. JRC and coe cient is linearly related. e analysis shows that whenthefractureapertureissmall,theimpactofJRCismore signi cant. e methods and results may be of interest to industrial g e o l o g i s t s ,h y d r o l o g i s t sa n de n g i n e e r sw h of o c u so nr o c k aquifer systems, and modelers who seek to validate more comprehensive numerical treatments.

Conflicts of Interest
e authors declare that they have no con icts of interest.