Flow and Transport in Porous Media: A Multiscale FocusView this Special Issue
Determining the REV for Fracture Rock Mass Based on Seepage Theory
Seepage problems of the fractured rock mass have always been a heated topic within hydrogeology and engineering geology. The equivalent porous medium model method is the main method in the study of the seepage of the fractured rock mass and its engineering application. The key to the method is to determine a representative elementary volume (REV). The FractureToKarst software, that is, discrete element software, is a main analysis tool in this paper and developed by a number of authors. According to the standard of rock classification established by ISRM, this paper aims to discuss the existence and the size of REV of fractured rock masses with medium tractility and provide a general method to determine the existence of REV. It can be gleaned from the study that the existence condition of fractured rock mass with medium tractility features average fracture spacing smaller than 0.6 m. If average fracture spacing is larger than 0.6 m, there is no existence of REV. The rationality of the model is verified by a case study. The present research provides a method for the simulation of seepage field in fissured rocks.
Seepage problems of the fractured rock mass have always been a heated topic within hydrogeology and engineering geology. Many problems that are closely related to the research of fracture seepages, including the stability of dam foundation and seepage, the stability of bedrock side slopes under the influence of groundwater, fissure deposit, the prevention of mine water inrush, and the leakage and diffusion of nuclear waste. However, the study of the seepage problems of the fractured rock mass is still in its preliminary stage at present. The general method regards fractured media as porous media and uses the permeability tensor of porous media to describe the seepage characteristics of the fractured media. There are different types of structural plane in fractured rock mass, which lead to the complexity of rock mass characteristics. Therefore, the study on rock mass model is always one of the important problems in rock mechanics. The representative elementary volume (REV) is the basis to determine a rock mass mechanics model, and it is necessary to research the REV of fractured rock mass effectively, so that the REV size of the fractured rock mass can be determined.
The concept of the REV was the first introduced in continuum mechanics by Bear , and it is to be used to describe the flow in the porous media. The parameter of interest is both homogeneous and statistically stationary, which will ensure consistency in flow simulation studies. The REV is defined in two situations on unit cell in a periodic microstructure and volume containing a very large set of microscale elements, possessing statistically homogeneous properties. The REV has been discussed by many authors [2–17]. The REV of a fractured rock mass is the smallest volume in during the study of parameter when the hydraulic conductivity is a constant value. One special concern is the evaluation of the REV of the fractured rock masses, due to the fact that fluid flow in fractured rock masses is of high scale effect [2, 18–23]. Previous studies assumed that the anisotropy was achieved by making use of different correlation lengths in the horizontal and vertical directions, and flow barriers were modeled stochastically [24–36].
Snow  concluded the math expressions of single and infinite fracture permeability tensors, assuming that the fracture seepages did not interfere with each other, while the overall permeability tensor of fracture network was the linear superposition of all fractures. As the fracture network was the same as the porous media, the permeability can be expressed by a permeability tensor. A fractured network can be approximately viewed as a porous medium, which was if the equivalent porous media of fractures are existent, then the permeability coefficients can be expressed by a symmetric second-order tensor. Li and Zhang  also conducted research on the REV of fractures. Bear  proved that if or in the polar coordinate system, mapping can form an ellipse (under the condition of three-dimensional ellipsoid), where represents radius vector, represents permeability coefficient in the direction of the flow line, and represents permeability coefficient in the direction of the hydraulic gradient. Long and Witherspoon  proposed that if the representative elementary volume of the fractured rock mass (REV) existed, the following two conditions must be met:(1)The media must be uniform in the area of study; that is, the average permeability coefficient of the area of study changes along with the scope of the study with big change; then it can be determined that this area of study is uniform.(2)In the polar coordinate system, the equivalent permeability coefficient k in each direction within the area of study can be described using an ellipse approximately at this time .
2. Introduction to FractureToKarst Software
This paper uses discrete fracture network (DFN) software, according to the standard of rock classification established by ISRM, and discusses the existence of REV of fractured rock masses with medium tractility. The discrete fracture network (DFN) software FractureToKarst for seepage in fracture rock mass is a kind of software using the Monte Carlo method. It can generate a two-dimensional fracture network of arbitrary shapes, set the common statistical parameters of fracture, filter the fractures within the area of study, and proceed with automatic discrete can be set. The head value and equivalent permeability coefficient of each node can be calculated by using the water balance principle. In the study, set an aspect ratio for the area of study 2 : 1 , calculate an equivalent permeability coefficient every 10° rotation of the area of study, every area of study can get 36 equivalent permeability coefficients, and discuss whether the direction of the equivalent permeability coefficient in polar coordinates can be described as a permeability coefficient ellipse or not, thus determining the existence of REV.
2.1. Mathematical Model
Mathematical model is the water dynamic model of the system. Fracture network is composed of a single, flat, smooth fracture, in a state of the laminar flow in a single fracture and can be described by the cubic law:
Type. is the boundary flux; is the density of water; is the gravitational acceleration; is flow dynamic viscosity coefficient; is fracture width; is two-head difference of the fracture; is the length of the fracture section.
In the fracture network, each node can establish a hydraulic link equation by the water balance principle:
Type. is the node flux; is the boundary flux. The hydrodynamic equations can be constructed by combining formula (1) and formula (2), and the head distribution of the fracture node within the system can be derived by the equations.
2.2. Model Algorithm
Model algorithm is a numerical method. The detailed steps are as follows.
(1) Grid discretization: all the fractures can be divided into the smallest fracture section by the nodes. Each node and connected fractures make a unit balance zone, as shown in Figure 1.
(2) For the hydraulic head, assign an initial value and boundary treatment: the node head is constant in the boundary of fixed water level. The initial water head of the rest of the nodes is 0.99 times the maximum elevation value. The node head of the impermeable boundary is equal to the adjacent nodes water head.
(4) Calculate the permeability coefficient: completing the hydraulic head calculation follows evaluating the equivalent permeability coefficient K in the flow direction using Darcy’s law.
(5) Calculate the permeability coefficient in any direction: keep the shape of the area of study unchanged, make it rotate around its center, and calculate an equivalent permeability coefficient every ~rotation; then we can get 36 equivalent permeability coefficients .
(6) Draw: in the polar coordinate system, it is appropriate to use the equivalent permeability coefficients of 36 directions in Step to draw the diagram, to see whether it can form an ellipse, so we know whether the REV is existent or not.
3. Analytical Solution for Validation of the Software FractureToKarst
Assume that there is a fracture network as shown in Figure 2, the right and left boundaries being the given head border and the upper and lower two boundaries being the impermeable boundary; the start-point and end-point coordinates of the fractures are .
Assume the hydraulic head of two kinds of in and out of boundary values , , respectively. Using the principle of water balance and the cubic flow law, there is
If , there are
Type. ~ is the water head; is the length of fracture section between the two nodes; is the density of water; is the gravitational acceleration; is flow dynamic viscosity coefficient.
Through (4), we can obtain that the values of and the hydraulic head values for each node are
The equivalent permeability coefficient of the flow direction is obtained through Darcy’s law:
Use FractureToKarst to build fracture network, and input parameters, the value of the water head, and the equivalent permeability coefficient are shown in Figure 3.
The calculation result in that the program equals the manual computation result, proving that the program is correct.
4. Simulation of the Fracture Network
According to the standard of rock classification established by ISRM, fractured rock mass with medium tractility refers to rock mass whose trace length is more than 3 m and less than 10 m. Because the permeability of two sets of orthogonal fracture rock masses is closest to being isotropic, the described ellipse is closest to being a circle.
According to the fracture spacing classification of ISRM, the spacing within 20~60 mm is very dense spacing. Within 10 m × 10 m, two sets of orthogonal fractures were generated, the average spacing is 0.06 m, and the aperture is 0.0001 m. The two sets of fracture parameters are shown in Table 1 (two sets of fracture identification for I and II in Table 1). Here, put some fractures in the same direction as a set of fractures. The distribution types of the trace length and the direction are the normal distribution, and the gap width is the logarithmic normal distribution. The right and left boundaries are the given head border and the upper and lower two boundaries are the impermeable boundary.
From Table 1, two sets of orthogonal fractures can be built. The first set has 276 fractures, and the second set has 208 fractures, with a total of 484 generating fractures. The diagram of generated fractures is shown in Figure 4(a). The equivalent permeability coefficient of the study area changes with the scope of the study without changes, so the study area is uniform. Within 10 m × 10 m, select five study areas of 1.0 m × 0.5 m, 2.0 m × 1.0 m, 3.0 m × 1.5 m, 4.0 m × 2.0 m, and 6.0 m × 3.0 m, and the diagram of generated fractures of FractureToKarst is shown in Figures 4(b)–4(f).
(a) All the fractures
(b) 1.0 m × 0.5 m
(c) 2.0 m × 1.0 m
(d) 3.0 m × 1.5 m
(e) 6.0 m × 3.0 m
(f) 7.0 m × 3.5 m
Acquire the permeability coefficients of all the directions in the five regions and make a comprehensive comparison chart of permeability coefficients (shown in Figure 5). Analyzing Figure 5, it can be concluded that every equivalent permeability coefficient is basically stable when the area of study is larger than 1.0 m × 0.5 m. There is no dramatical change for the permeability coefficients in the four regions 2.0 m × 1.0 m, 3.0 m × 1.5 m, 4.0 m × 2.0 m, and 6.0 m × 3.0 m, so more than 1.0 m × 0.5 m in the area of study can be approximately viewed as a homogeneous medium region.
In Figure 5, the meaning of all the symbols as follows.
“△” is the area of study of 1.0 m × 0.5 m; “●” is the area of study of 2.0 m × 1.0 m; “▽” is the area of study of 3.0 m × 1.5 m; “○” is the area of study of 6.0 m × 3.0 m; “⋄” is the area of study of 7.0 m × 3.5 m.
5. The Fitting Calculation of the Uniform Basin
5.1. Calculation of Regional Rotation
As the uniform basin 1.0 m × 0.5 m was determined above, taking the initial angle as the angle between the horizontal direction and the flow direction, rotate the whole area of study clockwise to calculate an equivalent permeability coefficient every 10°; then each area of study has 36 equivalent permeability coefficients, as shown in Table 2.
Table 2 shows that the maximum value of the permeability coefficient is 0.127 m/s when the angle is 0° and 180° between the flow direction and the horizontal orientation (because the angle between the two direction is 180°, so they are the same one flow field); the minimum value of the permeability coefficient is 0.076 m/s when the angle is 40° and 220° between the flow direction and the horizontal orientation. The reason for this phenomenon is that the calculation in Table 2 is derived from the data of the fracture network in Table 1, where 0° (180°) and 90° (270°) are two groups of orthogonal fractures. In Table 1, at four points at angles 0° (180°) and 90° (270°), there will be a local maximum value. The permeability coefficient at 0° (180°) is the local maximum value of fractures in the 0° (180°) direction, and the permeability coefficient at 40° (220°) is the minimum value between the two local maximum values.
5.2. The Polar Coordinate Fitting of the Permeability Coefficient
Assuming that point P is any one point on the ellipse in the polar coordinate system, A, B are, respectively, the two endpoints, and C is the focus of the ellipse. The semi~major axis and the semi~minor axis of the ellipse are assumed as and , respectively. This results in
Assuming that , is a parameter, that is, the angle between the principal axis of the ellipse and the polar coordinate 0° axis.
Solve (7), where the ellipse equation of the permeability coefficient in the polar coordinate system is
The semi~major axis and the semi~minor axis of the fitting ellipse are , , respectively, where radians, and the fitting equation is
6. Determination of the REV
Öhman et al. have conducted a considerable amount of useful research [41–43] on the equivalent medium in the fractured rocks to compare the similarity between the equivalent permeability coefficient of numerical simulation and the ellipse in different area of studies. These two scholars have obtained the minimum values of the square of the simulation value and fitting value, with the following formula:
RMS is the fitting correlation coefficient of the ellipse, is the simulation value, and is the fitting value.
Depending on its similarity to the elliptic curve, RMS can be divided into three categories:(1)When RMS 0.2, the size of the area of study can be used as the REV of the fractured rock mass.(2)When 0.4 > RMS > 0.2, the size of the area of study can be used as the REV of the fractured rock mass under certain conditions.(3)When RMS 0.4, the size of the area of study cannot be used as the REV of the fractured rock mass.
Therefore, after comparing Tables 2 and 3, it can be gleaned that the fitting result is quite ideal, with the correlation coefficient RMS = 0.132064722. The REV of the fractured rock mass of this kind exists, whose size is 1.0 m × 0.5 m.
Further research shows that the REV of very dense spacing (<0.02 m) of the fractured rock mass exists, and it is less than 1.0 m × 0.5 m. The REV of the dense spacing (0.06~0.2 m) of the fractured rock mass is 2.0 m × 1.0 m. The REV of medium spacing (0.2~0.6 m) of the fractured rock mass is 10.0 m × 5.0 m. The REV is not the existence of the fractured rock mass of wide spacing (0.6~2.0 m), with very wide spacing (2.0~6.0 m) and utmost spacing (>0.6 m) because their fractures have no connection and have no hydraulic conductivity, and thus the REV does not exist.
7. Conclusions and Suggestions
In conclusion, the existence of REV is closely related to the fracture conditions. Not all types of fractured rock mass have REV. The more intensive the fractures are, the better the penetration is, the better the permeability of the rocks is, which means the easier it is to become equivalent to porous media. Thus, the following conclusions can be made:(1)The existence condition of the REV for fractured rock mass with medium tractility is that the average spacing of the fractures should be less than or equal to 0.6 m; that is to say, the size of the REV of extremely dense spacing of the fractured rock mass is less than 1.0 m × 0.5 m; the size of the REV of very dense spacing of the fractured rock mass is 1.0 m × 0.5 m; the size of the REV of dense spacing of the fractured rock mass is 2.0 m × 1.0 m; the size of the REV of medium spacing of the fractured rock mass is 10.0 m × 5.0 m.(2)When the average spacing of the fractures is more than 0.6 m, the REV of wide spacing and very wide spacing and extremely wide spacing of the fractured rock mass does not exist.(3)Although this study has obtained certain achievement in the scale effects of the REV of the discrete fracture medium, there are still some shortcomings, including the influence of the different distribution of the geometric elements and the influence of the distribution of the different gap lengths, which all require further discussion.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
This work was financially supported by the Foundation for Teachers from China Earthquake Administration (Grant no. 20150102), by the Foundation for Hebei Province Science and Technology Planning Project (Grant no. 162776436), and by the National Natural Science Foundation of China (Grant no. 41272387).
J. Bear, Dynamics of Fluids in Porous Media, Elsevier, New York, USA, 1972.
M. Hassanizadeh and W. G. Gray, “General conservation equations for multi-phase systems, 1: averaging procedure,” in Flow through porous media, recent developments: a computational mechanics publications, G. F. Pinder, Ed., pp. 1–16, 1983.View at: Google Scholar
H. H. Haldorsen, “Simulator parameter assignment and the problem of scale in reservoir engineering,” in Reservoir characterization, L. W. Lake and H. B. Caroll, Eds., pp. 293–340, Academic Press, Orlando, USA, 1986.View at: Google Scholar
G. H. Shi, “Manifold method of material analysis,” in Transaction of the 9th Army Conference on Applied Mathematics and Computing, pp. 57–76, USA, 1991.View at: Google Scholar
G. H. Shi, “Manifold method,” in Manifold method. Proceeding of first International Forum on Discontinuous Deformation Analysis (DDA) and Simulations of Discontinuous Media. TSI Press, pp. 52–204, Albuquerque, NM, Mexico, USA, 1996.View at: Google Scholar
G. Shi, “Producing joint polygons, cutting joint blocks and finding key blocks for general free surfaces,” Chinese Journal of Rock Mechanics and Engineering, vol. 25, no. 11, pp. 2161–2170, 2006.View at: Google Scholar
A. Hurst, “Sedimentary flow units in hydrocarbon reservoirs: some shortcomings and a case for highresolution permeability data,” in The geological modelling of hydrocarbon reservoirs and outcrop analogues, S. Flint and I. D. Bryant, Eds., pp. 191–204, 1993, Special publication No. 15 of the international association of sedimentologists.View at: Google Scholar
P. S. Ringrose, Pickup G E, J. Jensen L, and M. Forrester, “The Ardross reservoir gridblock analogue: sedimentology, statistical representative and flow upsclaing,” in Reservoir characterization recent advances. Am assoc petrol geol memoir, R. Schatzinger and J. Jordan, Eds., vol. 71, pp. 265–276, The Ardross reservoir gridblock analogue, sedimentology, 1999.View at: Google Scholar
K. Nordahl and P. Ringrose, Identifying the representative elementary volume for permeability in heterolithic deposits using numerical rock models. Math Geosci, vol. 40, 2008.View at: Publisher Site
K.-B. Min, L. Jing, and O. Stephansson, “Determining the equivalent permeability tensor for fractured rock masses using a stochastic REV approach: Method and application to the field data from Sellafield, UK,” Hydrogeology Journal, vol. 12, no. 5, pp. 497–510, 2004.View at: Publisher Site | Google Scholar
P. Blum, R. Mackay, M. S. Riley, and J. L. Knight, “Performance assessment of a nuclear waste repository: upscaling coupled hydro-mechanical properties for far-field transport analysis,” International Journal of Rock Mechanics and Mining Sciences, vol. 42, no. 5-6, pp. 781–792, 2005.View at: Publisher Site | Google Scholar
S. H. Begg, R. R. Carter, and P. Dranfield, “Assigning effective values to simulator gridblock parameters for heterogeneous reservoirs,” SPE Reservoir Engineering (Society of Petroleum Engineers), pp. 455–463, 1989.View at: Google Scholar
A. Henriette, C. G. Jacquin, and P. M. Adler, “The effective permeability of heterogeneous porous media,” Phys Chem Hydrodyn, vol. 11, no. 1, pp. 63–80, 1989.View at: Google Scholar
B. Noetinger and C. Jacquin, “Experimental Tests of a Simple Permeability Composition Formula,” Society of petroleum engineers preprint SPE 22841, 1991.View at: Google Scholar
P. Corbett, S. Anggraeni, and D. Bowen, “The use of the probe permeameter in carbonatesaddressing the problems of permeability support and stationarity,” Log Analyst, vol. 40, no. 5, pp. 316–326, 1999.View at: Google Scholar
V. C. Tidwell and J. L. Wilson, “Heterogeneity, permeability patterns, and permeability upscaling: Physical characterization of a block of Massillon sandstone exhibiting nested scales of heterogeneity,” SPE Reservoir Evaluation and Engineering, vol. 3, no. 4, pp. 283–291, 2000.View at: Publisher Site | Google Scholar
M. D. Jackson, A. H. Muggeridge, S. Yoshida, and H. D. Johnson, “Upscaling permeability measurements within complex heterolithic tidal sandstones,” Math Geol, vol. 35, no. 5, pp. 446–454, 2003.View at: Google Scholar
D. T. Snow, A parallel plate model of fractured permeable media: [dissertation] [dissertation, thesis], Berkeley, Univ. of Calif, [dissertation], 1965.
Q. Yu, D. Chen, and G. Xue, “Hydrodynamics of discontinuous fracture network,” Earth Science-journal of China University of Geosciences, vol. 20, no. 4, pp. 474–478, 1995.View at: Google Scholar
J. Öhman, A. Niemi, and C.-F. Tsang, “Probabilistic estimation of fracture transmissivity from Wellbore hydraulic data accounting for depth-dependent anisotropic rock stress,” International Journal of Rock Mechanics and Mining Sciences, vol. 42, no. 5-6, pp. 793–804, 2005.View at: Publisher Site | Google Scholar