VZJ Journal of Natural Resources and Life Sciences Education
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
QUICK SEARCH:   [advanced]


     

Published online 16 December 2005
Published in Vadose Zone J 5:1-13 (2005)
DOI: 10.2136/vzj2004.0175
© 2005 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Seol, Y.
Right arrow Articles by Ito, K.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Seol, Y.
Right arrow Articles by Ito, K.
GeoRef
Right arrow GeoRef Citation
Agricola
Right arrow Articles by Seol, Y.
Right arrow Articles by Ito, K.
Related Collections
Right arrow Flow
Right arrow Fractured Rock
Right arrow Vadose Zone Processes and Chemical Transport

ORIGINAL RESEARCH

An Evaluation of the Active Fracture Concept in Modeling Unsaturated Flow and Transport in a Fractured Meter-Sized Block of Rock

Yongkoo Seol*, Timothy J. Kneafsey and Kazumasa Ito

Earth Sciences Division, Lawrence Berkeley National Laboratory, One Cyclotron Road, MS 90-1116, Berkeley, CA 94720
* Corresponding author (yseol{at}lbl.gov)

Received 7 December 2004.



   ABSTRACT
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
INTRODUCTION
 MODEL DEVELOPMENT
 NUMERICAL SIMULATIONS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Numerical simulation is an effective and economical method for optimally designing laboratory experiments and deriving practical experimental conditions. We executed a detailed numerical simulation study to examine the active fracture concept (AFC) using a 1-m3-sized block model. The numerical simulations for this study were performed by applying various experimental conditions, including different bottom flow boundaries, varying injection rates, and different fracture–matrix interaction (i.e., by increasing absolute matrix permeability at the fracture–matrix boundary) for a larger fracture interaction under transient or balanced-state flow regimes. Two conceptual block models were developed based on different numerical approaches: a two-dimensional discrete-fracture-network model (DFNM) and a one-dimensional dual-continuum model (DCM). The DFNM was used as a surrogate for a natural block to produce synthetic breakthrough curves of water and tracer concentration under transient or balanced-state conditions. The DCM is the approach typically used for the Yucca Mountain Project because of its computational efficiency. The AFC was incorporated into the DCM to capture heterogeneous flow patterns that occur in unsaturated fractured rocks. The simulation results from the DCM were compared with the results from the DFNM to determine whether the DCM could predict the water flow and tracer transport observed in the DFNM at the scale of the experiment. It was found that implementing the AFC in the DCM improved the prediction of unsaturated flow and that the flow and transport experiments with low injection rates in the DFNM compared better with the AFC implemented DCM at the 1-m3 scale. However, the estimated AFC parameter varied from 0.38 to 1.0 with different flow conditions, suggesting that the AFC parameter was not sufficient to fully capture the complexity of the flow processes in a 1-m3 discrete fracture network.

Abbreviations: AFC, active fracture concept • DCM, dual-continuum model • DFNM, discrete-fracture-network model


   INTRODUCTION
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
MODEL DEVELOPMENT
 NUMERICAL SIMULATIONS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
UNSATURATED FLOW and transport in fractured rock are difficult to predict because of the complexity of fracture network configurations and heterogeneity in fracture–matrix interactions. Meaningful measurements are hard to make because of technological constraints at the scale of interest. One unique phenomenon in unsaturated fractured rock is nonuniform flow along preferential pathways, which often accommodates flow faster than predicted from matrix hydraulic conductivity under saturated conditions (Fabryka-Martin et al., 1996; Pruess et al., 1999a).

To simulate gravity-dominated, nonuniform, preferential water flow within unsaturated media, Liu et al. (1998) implemented the AFC into a dual-permeability concept. In the AFC, the active fractures are defined as the portion of connected fractures actively conducting water within the unsaturated zone, and the fraction of active fractures depends on the water saturation in the fractures. Since only active fractures exchange water and tracers with the matrix, the actual interfacial area between flowing fractures and the matrix is reduced from the total interfacial area calculated from fracture network geometry. The AFC is of special interest because it has been used to represent the effects of nonuniform preferential flow in the fractured unsaturated zone of the proposed nuclear waste repository at Yucca Mountain, Nevada (Liu et al., 2000).

Numerical results suggest the presence of active fractures in the unsaturated zone (Gauthier et al., 1992; Liu et al., 1998, 2002; Pruess et al., 1999a; Zhang et al., 2002). Moreover, laboratory and field experiments provide indirect evidence of preferential flow in heterogeneous porous media and single fractures (Su et al., 1999; Wang et al., 1999). However, there are no controlled experimental results that directly prove their presence and examine the propriety of the AFC in a fracture network. To obtain experimental evidence, an experiment using a 1-m3-sized block of fractured rock has been proposed. However, performing the laboratory experiment will require a significant amount of time and resources, involving a large fractured rock block and numerous preliminary tests to find feasible and appropriate experiments by which to ensure that the experimental objectives are achieved.

Our numerical simulation study was motivated by the need to demonstrate the AFC effectively and economically and, furthermore, to provide practical experimental conditions for the proposed laboratory test using a 1-m3-scale block. A 1-m3-sized block is the largest practical size that can be obtained, transported, and tested in a laboratory experiment. The numerical simulations done for this study focused on determining the AFC parameter defining the fracture–matrix interaction reduction, as well as estimating appropriate experimental conditions that will help make the experimental objectives achievable.

This study adopted the approach outlined by Finsterle (2000). Two conceptual models, a DFNM and a DCM, were developed to simulate flow and tracer transport in a two-dimensional unsaturated fractured block. The DFNM was used as a surrogate block and served to generate synthetic water flow-rate data and tracer breakthrough curves under transient and steady-state conditions. Numerical simulation results from the DCM were compared with the DFNM results to evaluate DCM predictions of flow and transport characteristics. Calibrations of the DCM against the synthetic DFNM data were performed to estimate the AFC parameter for the fracture network.

In this paper we report the development of the two conceptual models and the approach taken in this study to examine the AFC; this is followed by a discussion on calibrated AFC parameters for different flow conditions and recommendations of appropriate experiment conditions.


   MODEL DEVELOPMENT
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MODEL DEVELOPMENT
NUMERICAL SIMULATIONS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Numerical simulations of flow and transport were performed using two conceptual models, the DFNM and the DCM. A DFNM was developed as a surrogate of a natural fractured rock block to synthesize flow and transport data sets. The results of DCM and DCM with the AFC simulations were compared with the synthetic data generated by the DFNM. The DCM implemented with the AFC is the model used for the Yucca Mountain Project because it efficiently and economically represents the model and yields reasonable results at such a large scale of interest. The DFNM selected for this study is two-dimensional because current information available for building a discrete fracture network is limited for a three-dimensional network and also because the three-dimensional discrete network is obviously too large to handle with readily available computational capability.

Discrete Fracture Network Model
To make the DFNM analogous to a natural rock block, the model incorporated statistical information derived from a detailed line survey of fracture data performed at the Exploratory Studies Facility at Yucca Mountain (Mongano et al., 1999), which documented any discontinuities having a minimum trace length of 1 m and intersecting the survey line within 30 cm on either side of an observation line along the drift at the proposed repository site. The detailed fracture survey characterized fracture density, ranges of aperture and trace length, and distribution of fracture orientation.

Fractures with different trace lengths and orientations obtained from the field observations were sorted into groups and counted to calculate the occurrence probability of fractures for each fracture group. For the artificial fracture network generation, fractures were initially placed randomly, and a trace length and orientation were independently assigned based on their occurrence probability. The total number of generated fractures was determined by the fracture line density (the number of fractures per unit length).

The fracture density for the lower lithophysal unit in Yucca Mountain, where a rock block was excavated, was 3.2 fractures per meter, based on a line survey of fractures of length 1 m or greater (Ahlers and Liu, 2000). However, we assumed the same fracture frequency along a line perpendicular to the survey line. Thus, to incorporate intersecting fractures parallel or near-parallel to the survey line, we increased the fracture density to 6.4 fractures per meter.

In addition to increasing the fracture density, we also added small fractures (<1 m) on the basis of small-fracture survey data (CRWMS M&O, 1999) to compensate for the presence of fractures ignored in the detailed line survey (Mongano et al., 1999). In general, larger fractures control global flow through the fracture system, whereas small fractures impact interfacial phenomena between fracture and matrix by providing more interfacial area to the system (Wu et al., 2004).

The two-dimensional DFNM model (1.0 by 1.0 m, with an arbitrary thickness set to 0.01 m) contained 39 fractures of various lengths and orientations (Fig. 1 , Table 1). The aperture was correlated to fracture length as described by Liu and Bodvarsson (2001). Based on the relationship between fracture trace length and filling material (Chiles and de Marsily, 1993), the average aperture, b, was calculated using an empirical equation as a function of trace length, L (m):

[1]
where c and d are empirical constants determined to be 1.008 x 10–4 and 0.317, respectively (Liu et al., 2000). Representative mean and variance values of log(b) for the lower lithophysal unit were –4.01 and 0.05, respectively (Robinson, 2000).



View larger version (139K):
[in this window]
[in a new window]
  		Fig. 1. Discrete fracture network model (DFNM). The thickness of the fracture lines indicates aperture.

 

View this table:
[in this window]
[in a new window]
  		Table 1. Ranges and mean values of fracture data.

 
The DFNM was constructed using the technique by Ito and Seol (2003). The fracture network was transformed into a mesh with 3750 (75 x 50) grid elements having uniform dimensions (0.0133 by 0.02 m) (Fig. 1). Each element in the mesh was assigned to either fracture or matrix. Elements cross-cut by any fractures were allocated to fractures. Matrix elements were assigned a uniform volume, whereas the fracture element volume was calculated using the fracture aperture and the trace length of a fracture segment cross-cutting the element.

Each element was connected to adjacent elements in each direction to allow complete interaction between matrix–matrix, fracture–fracture, and fracture–matrix grid blocks. The methods used for interfacial area calculation and the permeability assignments for connections between elements are depicted in Fig. 2 , where four fracture elements (F1–F4) and two matrix elements (M1 and M2) are shown. Interfacial areas between two matrix elements are the same as the unit areas of the elements (connection a between M1 and M2). Matrix permeability is assigned to the elements. For fracture–fracture element connections along the same fracture (F1 and F2, or F2 and F4), the interfacial area between the two fracture elements is the averaged aperture of the two fracture elements multiplied by the unit thickness of elements (connection b). Fracture permeability is assigned for the water flow between the elements. When two fracture elements of different fractures are connected (F1 and F3, or F3 and F4), the interfacial area is the average of the projected areas of the two fractures on the interface (connection c), and the matrix permeability is assigned to the area. For the interfacial area between fracture and matrix elements (F3 and M1, or F4 and M2), the projected area of a fracture segment onto the corresponding matrix element was used (connection d), and the matrix permeability is assigned to the area. When two or more fractures cross-cut an element, the sum of projection area was the area for the fracture element, and the area was limited by the full interface area. The total fracture–matrix interfacial area was 29.8 m2 m–3.



View larger version (23K):
[in this window]
[in a new window]
  		Fig. 2. Schematic diagram depicting interface area calculation and permeability assignment for connections.

 
We assumed a normal distribution for the fracture aperture within a single fracture, and thus each fracture segment within a fracture element had a different aperture (Table 1). All the fractures were assigned surface roughness factors based on observations by Mongano et al. (1999), and the roughness factor was converted into the standard deviation of the aperture distribution. Based on the cubic law, the original permeability for each segment in an element was modified by a permeability modification factor M, which is the square of the ratio between the randomly sampled aperture and the average aperture. Permeability and capillary pressure for each fracture were modified as follows:

[2]

Changing the aperture of a parallel-plate fracture segment changes the capillary strength as follows:

[3]

In these relations, kf is the permeability of a fracture segment, kave is the average permeability of a fracture, bi is the randomly sampled aperture, bave is the average aperture, M is the permeability modification factor, Pc is the capillary pressure of a fracture segment, and Pc,ave is the average capillary pressure of a fracture. Using the permeability modification factor, heterogeneities in aperture, permeability, and capillary pressure can be incorporated for each individual fracture, as well as for each fracture segment.

Input parameters for fracture and matrix are listed in Table 2 (Ahlers and Liu, 2000). For relative permeability calculations, Corey's curve (Corey, 1954) was selected, whereas the van Genuchten model was used for the capillary pressure function (van Genuchten, 1980).


View this table:
[in this window]
[in a new window]
  		Table 2. Input parameters used for modeling (Ahlers and Liu, 2000).

 
Dual Continuum Model
The one-dimensional DCM was developed to model the same dimensions (1.0 by 1.0 by 0.01 m) as the DFNM, with the mesh generated using the TOUGH2 Meshmaker module for dual permeability (Pruess et al., 1999b, 1996). The DCM consists of a set of parallel planar fractures with matrix between neighboring fracture planes. The fracture spacing was set at 0.06695 m to have the same fracture–matrix interfacial area (29.8 m2 m–3) as the DFNM. Total volumes for the fracture and matrix continua were identical to those in the DFNM. The same input parameters (Table 2) were used except for absolute fracture permeability (kf), which was initially set equal to the equivalent permeability (ke) of the saturated DFNM (see Saturated Flow Simulations below). The van Genuchten relative permeability and capillary pressure functions were used for the constitutive relations of unsaturated flow in the model (van Genuchten, 1980).

Active Fracture Concept
The AFC (Liu et al., 1998) incorporates gravity-dominated preferential liquid flow into the dual-continuum approach for modeling unsaturated flow and transport in fractured rock. It assumes that the number of active fractures within a rock block (control volume) is large enough to allow use of the continuum approach. The fraction of active fractures in a fracture network, fa, is assumed to be a power function of effective water saturation in the connected fractures, Se:

[4]

[5]
where {gamma} is a positive constant depending on properties of the fracture network, Sf is the water saturation of all connected fractures, and Sr is the residual fracture saturation. The fracture–matrix interfacial area reduction factor, afm, is also expressed as a function of effective fracture saturation (Se):

[6]

Capillary pressure and relative permeability functions are modified to account for larger effective saturations in the active fractures than that of the total fracture continuum. Details on the AFC can be found in Liu et al. (1998).

Demonstrating the Active Fracture Concept
The AFC consists of two main features. First, only a portion of the fracture network contributes to water flow. Second, the active interfacial area between fracture and matrix is reduced with the formation of preferential flow pathways. Consequently, evaluation of the AFC can be approached in two different ways, by looking for the portion of fractures in the fracture network contributing to water flow and by estimating interfacial area reduction. The AFC can be demonstrated if active fractures can be observed and quantified in laboratory or field experiments, or if reduced interfacial area for water flow and transport through the preferential pathways can be predicted by numerical modeling. The first AFC demonstration approach requires a large-scale flow domain with numerous fractures; hence, laboratory demonstration for this aspect of the AFC would be difficult to perform with limited time and resources. We focused on the second approach (estimating the reduction in interfacial area between fractures and the matrix).

Demonstrating the second feature of the AFC requires finding the AFC parameter ({gamma}) in Eq. [6]. The {gamma} value is a unique property of a given fractured rock, interpreted as a measure of the "activity" of connected fractures (Liu et al., 1998). To determine the {gamma} value, we calibrated the DCM results against the DFNM data using inverse modeling, fixing all other parameters and solving for the {gamma} value.

Fracture–Matrix Interaction
The AFC attempts to capture the effects of reduced fracture–matrix interaction in a fracture system in which not all fractures carry water. In a system where fracture–matrix interaction is minimal, the system response modeled by the AFC cannot be distinguished from that calculated by the standard DCM. Fracture–matrix interaction is potentially important if (i) matrix permeability is high or (ii) the flow path is long enough for substantial water imbibition. Increasing either matrix permeability or the scale of the system enhances fracture–matrix interaction. Data from large, high-permeability systems are potentially more suitable for demonstrating the AFC, because they better reveal the processes captured by the AFC. In our current DFNM model, the fracture–matrix interface area can be minimal because of the limited fracture–matrix interface area and the lack of data on fine cracks and fractures that would provide a large area for the interaction. To test the impact of enhanced fracture–matrix interaction on the system response captured by the AFC, we simply increased the absolute permeability of matrix at the fracture–matrix connections by a factor of 100. Within a heterogeneous geologic unit, one would not expect a range in fracture–matrix interfacial area to exceed two orders of magnitude.


   NUMERICAL SIMULATIONS
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MODEL DEVELOPMENT
 NUMERICAL SIMULATIONS
RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
The conceptual models, DCM and DFNM, were initially dry and subjected to wetting with varying injection rates that were less than the saturated flow rate, which was determined as the model block's maximum flow rate. The flow rate of the block effluent was monitored until it reached balanced state. Once the balanced state was achieved, a pulse of tracer was released at the top boundary for transport simulations. The same series of simulations was executed for the two conceptual models. Calibrations of the DCM with the AFC against data generated from the DFNM were performed, both to estimate the AFC parameter under various conditions and to demonstrate the ability of the DCM with the AFC to predict the flow and transport occurring in the DFNM. In these calibrations, the AFC parameter was the only parameter estimated.

Atmospheric pressure was applied at the top boundary, and no-flow boundaries were assigned at side boundaries. The bottom boundary was assigned to be either a Dirichlet-type boundary representing open atmosphere (i.e., capillary-barrier boundary) or a free-drainage boundary. The former conditions mimic a laboratory setup in which the rock block has a free bottom, whereas the latter free-drainage condition assumes that the block extends continuously downward, and observations are made at a certain depth. In this case, the capillary pressure immediately below the lower model boundary was set equal to that at the model boundary. Water was released uniformly at various injection rates across the top, and water flow rates were monitored at the bottom.

Flow and transport simulations with various flow rates were performed with the two different bottom flow boundary conditions. Simulations with the larger matrix permeability at the connection between fracture and matrix were performed only with the free-drainage bottom boundary condition.

We used the TOUGH2/iTOUGH2 codes implemented with the AFC (Pruess et al., 1999, 1996; Finsterle, 1999) for this study. TOUGH2 is an integral finite difference code for simulating multidimensional coupled flow and transport of multiphase, multicomponent fluids in porous and fractured media. The iTOUGH2 is an inverse modeling code based on the TOUGH2 simulator.

Saturated Flow Simulations
Saturated flow simulations were performed to find the steady-state saturated flow rates and the equivalent permeability (ke) of the DFNM. The latter was used as the fracture permeability (kf) for the DCM. To calculate the ke of the DFNM block, we assumed that the DFNM consisted of a single homogeneous medium with constant injection of water at the top. All gridblocks were initially saturated with water, and flow rates were monitored at the bottom of the block. The ke was calculated with the saturated flow rate using Darcy's Law. We assumed flow through the matrix to be negligible because the matrix permeability is six orders of magnitude smaller than the fracture permeability, kf. Therefore, fracture flow governs the flow rate within the entire block.

Unsaturated Flow Simulations
Unsaturated flow simulations were performed (i) to compare the breakthrough curves of flow rates from the two conceptual models and (ii) to estimate the AFC parameter, {gamma}, by calibrating the flow rates calculated using the DCM against the synthetic data from the DFNM. The initial liquid saturations of the model blocks were set to 45% in the matrix and 2% in the fractures, which establishes approximate capillary equilibrium between the two media. Flow rates applied to the block models varied from 0.01 to 20% of the saturated flow rate. These flow rates were monitored at the top and bottom of the model blocks, until the difference between inflow and outflow rates became <1% of the injection rate (a balanced state). Water saturation distributions in the DFNM were also monitored regularly to observe possible flow focusing along fractures.

Unsaturated Transport Simulations
Transport simulations were only performed when the balanced state of flow was attained because at low matrix saturations, most of the tracer would be transferred into the matrix as the matrix imbibes water from fractures, and a tracer pulse breakthrough would not take place during a time frame reasonable for laboratory experiments.

When the balanced state was achieved, a finite amount of tracer was instantaneously introduced at the top of the model. The tracer mass exiting the block was monitored, and breakthrough curves were plotted. The cumulative tracer mass from the DCM was calibrated against the data from the DFNM for each flow rate to estimate the parameter {gamma}. Tracer distributions within the DFNM were also monitored.


   RESULTS AND DISCUSSION
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MODEL DEVELOPMENT
 NUMERICAL SIMULATIONS
 RESULTS AND DISCUSSION
CONCLUSIONS
 REFERENCES
 
Saturated Flow Simulations
Using saturated flow simulations, an effective permeability ke of 3.20 x 10–13 m2 was obtained for the DFNM based on a calculated balanced-state unit-gradient flow rate of 3.56 x 10–5 kg s–1. Neglecting the insignificant contribution from matrix flow, the effective permeability was used as the fracture permeability kf for the DCM. The saturated flow rate corresponds to an infiltration rate of 112 m yr–1, which is approximately four orders of magnitude greater than the average infiltration rate (5–20 mm yr–1) observed at Yucca Mountain (Montazer and Wilson, 1984; Flint et al., 1997). The flow rates used for the simulations ranged from 22.5 to 22 454 mm yr–1, corresponding to 0.02 to 20% of the saturated flow rate, in an attempt to simulate conditions comparable with natural conditions and, at the same time, more practically attainable in laboratory experiments. Henceforth, unless noted, the injection flow rates are represented by the percentage of the saturated flow rate.

Transient Unsaturated Flow Simulations
Numerical simulations were initiated with low initial liquid saturation (e.g., 45% of matrix saturation) without adding tracers, and water flow rates at the block bottom were monitored. Water saturation distributions in the entire model were also observed to find major flow paths. Flow simulation continued until the flow rates became balanced. The transient flow simulations provided breakthrough curves of water flow, as well as time requirements for laboratory experiments (Fig. 3 and 4) and transitional patterns of the liquid saturation distribution (Fig. 5 ).



View larger version (29K):
[in this window]
[in a new window]
  		Fig. 3. Water breakthrough curves of the DFNM with bottom capillary-barrier condition and various injection rates (top: 1–20%, bottom: 0.02–0.2% of saturated flow rate).

 


View larger version (30K):
[in this window]
[in a new window]
  		Fig. 4. Comparisons of water breakthrough curves in the DFNM with different bottom conditions and various injection rates (top: 1–2%, bottom: 10–20% of saturated flow rate).

 


View larger version (130K):
[in this window]
[in a new window]
  		Fig. 5. Distributions of water saturation in the DFNM with capillary-barrier bottom boundary condition when water reaches the bottom of model block (top left: 2% injection rate, t = 4.8 d; bottom left: 20% injection rate, t = 0.48 d; top right: 2% injection rate, t = 6.9 d; bottom right: 20% injection rate, t = 0.69 d). {gamma} = 0.

 
Time Requirements for Laboratory Experiments
Figure 3 shows normalized outflow rates increasing as a function of time for the DFNM with a bottom capillary barrier boundary condition. At injection rates of 0.02% (22.5 mm yr–1) or lower, the time to reach the balanced state exceeds 3 yr and is thus impractical for a laboratory experiment. As far as the time requirement for laboratory experiments is concerned, Fig. 4 illustrates that the different bottom flow boundary conditions and fracture–matrix interactions provide consistent time frames for flow breakthroughs. Given the result from the two-dimensional flow domain and assuming no significant discrepancies from three-dimensional flow systems not considered in the simulations, we concluded that flow rates >0.2% of the saturated flow rate should allow proposed laboratory experiments to be completed in a reasonable time frame.

Discrete Flow Behavior in the DFNM
Stepwise increases in the flow rates in Fig. 3 for lower injection rates (<2%) indicate discrete flow behavior. The first leap in the flow rates was caused by water flowing in the most preferential flow pathway, and subsequent leaps were induced by more flowing fractures as saturation in the fractures and matrix in the vicinity of flowing fractures increases. Water from the matrix would come last. At higher injection rates (>2%), water flowed rapidly through a few preferential flow pathways, and the flow rates quickly reached the injection rates, even though significant portions of the matrix were still relatively dry due to slow imbibition. Continuous observations are needed to capture discrete flow behavior at higher injection rates.

Another discrete flow behavior was found from liquid saturation distributions in the DFNM, which showed flow focusing (i.e., higher liquid saturation regions) even for the higher-injection-rate cases (>2%) (Fig. 5). Independent of bottom boundary conditions and fracture–matrix interaction, the advancing water front formed two preferential flow pathways. This discrete liquid saturation pattern resulted partially from variations in fracture aperture and capillary strength, as well as from nonuniform fracture connectivity.

The patterns of liquid saturation distribution are quite distinctive for different injection rates. At a lower injection rate (2%), most of the introduced water was imbibed into the matrix. When the water reached the block bottom, large portions of the matrix near the water flow pathways were nearly saturated (Fig. 5, upper left panel). Capillary forces impact the water flow and saturation distribution more strongly at the low injection rate than at higher injection rates. The higher injection rate (20%) resulted in higher saturation within a limited region immediately adjacent to fractures, with gravity-driven water flowing primarily through the fractures when water reached the bottom of the block (Fig. 5, lower left panel). Figure 5 (right panels) shows the liquid saturation distribution at later times when more flow through additional preferential flow pathways reached the bottom, corresponding to the second leap in the bottom flow rates shown in Fig. 3.

The simulations predict the occurrence of preferential flow pathways in a laboratory experiment with a 1-m3-scale block. This discrete flow behavior indicates that the DFNM block simulates multiple nonuniform preferential flow pathways through the fracture network and that the model as a surrogate of natural unsaturated fractured rock can be adopted for examination of the DCM implemented with the AFC.

Effect of Boundary Conditions on Unsaturated Flow
The free-drainage boundary condition at the bottom of the model block induced earlier breakthrough of injected water from the block, while the capillary-barrier boundary delayed the outflow until the weight of water accumulated in fractures exceeded the capillary pressure, as shown in Fig. 4. Even though less sharp increases appeared in the stepwise breakthrough curves, the free-drainage condition generated similar breakthrough curves and the same median breakthrough time. When the free-drainage condition was employed along with the higher matrix permeability at the fracture–matrix interface, the delayed breakthrough was due to higher imbibition of water into the matrix. Overall, the breakthrough curves became steeper as the matrix permeability increased with the free-drainage condition (Fig. 4). Figure 5 shows the water saturation distribution with the capillary-barrier boundary for two simulations with different flow rates after the same amount of water has been introduced to the DFNM. A broader region of the higher saturations (>90%) was observed for the lower flow rate, as would be expected because larger residence time allows imbibition to occur.

Note that capillary rise at the bottom of the model block with the capillary-barrier boundary condition (Fig. 5) does not occur within the block with the free-drainage condition (Fig. 6 ). Figure 6 compares the water saturation distributions with the bottom free-drainage condition (top panel) to the free-drainage condition with 100 times the matrix permeability (bottom panel) at the fracture–matrix connections. The free-drainage boundary condition resulted in the consistent pattern of two preferential flow pathways. The high imbibition with elevated matrix permeability causes a broader high-saturation region near the flow pathways (Fig. 6, bottom panel).



View larger version (82K):
[in this window]
[in a new window]
  		Fig. 6. Distributions of water saturation in the DFNM with different bottom boundary conditions when water reaches the bottom of the model block: free-drainage condition (top), free-drainage condition with 100 times higher matrix permeability at fracture–matrix connections (bottom). 2% of injection rates, t = 4.8 d, and {gamma} = 0.

 
Comparison of the breakthrough curves with saturation distribution patterns shows that the bottom flow boundary condition does not have much impact on the simulation results. This result would suggest that the laboratory experiments could be designed with a less complicated setup (i.e., a simple capillary-barrier condition will suffice).

Calibrations with the DCM
In Fig. 7 (top panel), the breakthrough curves for the DCM with {gamma} = 0 are compared with the curves from the DFNM. Discrepancies between the two models are apparent. The wetting front arrived at the bottom much faster in the DFNM at all flow rates, representing preferential flow through major flow pathways in the DFNM. Moreover, the DCM results did not show the stepwise increases in flow rates as observed in the DFNM, indicating that water seepage proceeded uniformly throughout the entire model width. The DCM was not expected to capture the step-wise flow rate changes, and the results suggest that the DCM does not accurately represent the discrete behavior of flow in unsaturated fractured rock with preferential flow pathways.



View larger version (39K):
[in this window]
[in a new window]
  		Fig. 7. Calibration of water flow rates to estimate {gamma} values for DCM with bottom capillary-barrier boundary condition and selected injection rates. {gamma} = 0 (top); {gamma} is estimated (bottom). Symbols are data from the DFNM and the symbols with line are from the DCM.

 
The {gamma} value was estimated from calibrating the flow rates calculated with the DCM against the flow rate data synthetically generated with the DFNM, while all other hydrologic parameters were fixed (Table 3). The {gamma} value was estimated by initially assigning a non-zero value. Estimation was constrained to the range between 0.01 and 1. Breakthrough curves showed better fits with the calibrated {gamma} values than with {gamma} = 0 (Fig. 7, bottom panel). With the calibrated {gamma}, the first appearance of water fronts in the DCM coincides more closely with the data from the DFNM, and the general trends of the DCM curves are more consistent with breakthrough curves from the DFNM, even though significant discrepancies are still found in the breakthrough curves. For different injection rates, the estimated {gamma} values vary from 0.39 to 1.0.


View this table:
[in this window]
[in a new window]
  		Table 3. Estimated AFC parameter,{gamma}, for simulations with different bottom flow boundary conditions and flow rates.

 
Because of the different characteristics of the two models, the DCM did not (and is not expected to) capture the discrete nature of liquid flow occurring in the DFNM (i.e., stepwise increase in flow rates). Capillary pressure in the DCM is by definition uniform throughout the entire bottom boundary in the DCM, but it can vary in the DFNM with different fracture properties. Consequently, DCM simulations represented by the average capillary pressure did not reproduce the DFNM's earlier water seepage through fractures having larger apertures and lower capillary pressures.

At low injection rates (<0.1%), capillary pressure predominantly controls liquid flow in fractures, and water exits the system only when the weight of accumulated water exceeds the capillary pressure. Reducing fracture–matrix interfacial area by increasing {gamma} was not sufficient to alter the water breakthrough curves because the capillary barrier at the block bottom would delay observation of water moving more rapidly through the block due to the reduction in interfacial area. Therefore, the presence of a capillary barrier at the outlet masks the observation most sensitive to changes in fracture–matrix interaction and thus prevents estimation of the related parameter {gamma}.

When the injection rate is large (>10%), a significant portion of injected water flows through a limited number of major flow pathways where the permeability is high. The imbibition rate of water into the matrix is too small compared with the high injection rates, and thus the outflow is relatively insensitive to matrix imbibition and even less sensitive to the fraction of matrix imbibition reduced by the interfacial area reduction factor.

In general, for injection rates greater than those constrained by the capillary-barrier effect, the {gamma} values were sensitive to flow rates as liquid saturation in fractures varied (Table 3). The estimated {gamma} values for the flow rates of 0.1 to 2% ranged from 0.39 to 0.59, and the flow rate estimation curves from the DCM with AFC implemented showed better agreement with the curves from the DFNM than those without the AFC (Fig. 7, bottom panel). Note that the AFC attempts to describe unsaturated flow and changes in fracture–matrix interface area as a function of flow rate using a single value for the {gamma} parameter. The fact that different {gamma} values were obtained for different flow rates suggests that the AFC continuum model is not capable of fully capturing the complexity of the flow processes in our discrete fracture network.

Figure 8 shows selected breakthrough curves of water flow under free-drainage bottom boundary conditions. Interfacial area reductions caused by the non-zero {gamma} yielded better fits between DCM and DFNM breakthrough curves (Fig. 8, top panel) compared with the zero {gamma} case. The estimated {gamma} values (0.38–0.49) are nearly the same as those determined using the capillary-barrier boundary condition (0.38–0.59); however, in general, the fits were not as good as those under the capillary-barrier boundary conditions, and the estimated {gamma} values vary with different flow rates.



View larger version (29K):
[in this window]
[in a new window]
  		Fig. 8. Calibration of breakthrough curve to estimate {gamma} values for the DCM with bottom free-drainage conditions and various injection rates (1 and 10% of saturated flow rates). Free-drainage condition (top), free-drainage condition with 100 times higher matrix permeability at fracture–matrix connections (bottom). Closed symbols are data from the DFNM, open symbols with dashed line are from the DCM with {gamma} = 0, and solid lines are the calibrations to estimated {gamma} values.

 
The model with the elevated matrix permeability provided sharp water breakthrough fronts without the AFC, and the water breakthroughs were highly sensitive to the AFC (Fig. 8, bottom panel). However, the elevated matrix permeability on fracture–matrix connection did not improve the fitting, indicating that at this model scale, larger interface area (i.e., more fractures) would not be a sufficient condition for the DCM to be able to predict flow behavior in the DFNM.

In summary, the fitting of the breakthrough curves is sensitive to {gamma} during a transient period, but the breakthrough curves from the two different models showed significant differences. The discrete flow behavior in the DFNM cannot be reasonably predicted by the DCM even with the AFC implemented, and the calibrated {gamma} values are not a unique property of the fracture network, but vary with flow rates.

Unsaturated Transport Simulations
A pulse of tracer was injected into the models once the outflow rates were balanced with the injection flow rates. Figures 9 and 10 show selected breakthrough curves (normalized total tracer mass exiting the system) and the distribution of tracer concentrations, respectively. The normalized tracer mass-breakthrough curves were plotted to compare the tracer transport in the different models with different bottom boundary conditions (Fig. 9). The DCM provided good matches to the DFNM in these balanced-state simulations, indicating that at the state, the DFNM acts more like a continuum model (such as the DCM). The differences in the early breakthrough curves were thought to result from mass transport through the limited fracture segments in the DFNM. Notice that there are some fracture segments that tracer did not reach in the earlier stage of transport (Fig. 10, top left panel).



View larger version (23K):
[in this window]
[in a new window]
  		Fig. 9. Calibration of tracer breakthrough curves to estimate {gamma} for the DCM with selected injection rates (2 and 10% of saturated flow rates). Capillary-barrier boundary condition (top); free-drainage condition (middle); free-drainage condition with 100 times higher matrix permeability at fracture-matrix connections (bottom). Closed symbols are data from the DFNM, open symbols with dashed line are from the DCM with {gamma} = 0, and lines are the calibrations to estimated {gamma} values.

 


View larger version (91K):
[in this window]
[in a new window]
  		Fig. 10. Distributions of tracer concentrations in the DFNM (free-drainage condition) 10% of injection rates at 0.03 d (top left); 10% of injection rates at 0.7 d (bottom left); increased matrix permeability with 10% of injection rates at 0.7 d (top right); 2% of injection rates at 1.7 d (bottom right).

 
The DFNM showed an early breakthrough with fast flows in major fractures and long tailing with significant imbibition of tracer into the matrix, and these phenomena appeared to result from the larger time allowed in the DFNM for the tracer transfer into the matrix or nonconductive fractures. The imbibed tracer became trapped in the stagnant flow region (Fig. 10, bottom left panel). The long tailing effect in the tracer breakthrough curves was more significant with longer residence times (associated with lower flow rates) (Fig. 10, bottom right panel) and the elevated fracture–matrix interaction (associated with the larger matrix permeability at the fracture–matrix connections) (Fig. 10, top right panel).

The breakthrough curves from the DCM were calibrated against the curves from the DFNM to estimate the {gamma} values. The normalized mass curves with DCM implemented with the AFC better capture the breakthrough of mass by limiting the fracture–matrix interface area. The calibrated {gamma} values range from 0.78 to 0.87 for 10 and 20% flow rates and nearly {gamma} = 1.0 for lower flow rates for both bottom flow-boundary conditions (Table 3). The calibrated {gamma} values also vary with flow rates and are not consistent with the values obtained from flow tests (Fig. 7 and 8).

For higher flow rates (>10%) with capillary-barrier and free-drainage boundary conditions, tracer transport is not very sensitive to {gamma} at the balanced state, and even without the AFC, the DCM captures the DFNM data quite well. This insensitivity arises from the short residence time (<1 d) of the tracer in the block, which has only a 1-m vertical travel distance. The sensitivity of transport to the {gamma} value is further lessened by slow imbibition of water when the matrix saturation is nearly in equilibrium with fracture saturation. For lower flow rates (<10%) with the capillary-barrier boundary condition and the free-drainage condition, and for all flow rates tested for the free-drainage boundary condition with elevated matrix permeability, the transport was sensitive to the {gamma} value, but the estimated {gamma} value reached nearly the upper limit (1.0), which means reduction in the fracture–matrix interface area reached the maximum at given saturation and the {gamma} value would not be sufficient to constrain the complicated tracer transport in the unsaturated fracture network. The early breakthrough can be captured by reducing the fracture–matrix interface area with the {gamma} value, but the later long tails in tracer mass likely caused by the slower transfer of tracer out of the matrix cannot be captured by the interface area reduction.


   CONCLUSIONS
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MODEL DEVELOPMENT
 NUMERICAL SIMULATIONS
 RESULTS AND DISCUSSION
 CONCLUSIONS
REFERENCES
 
Two models, a two-dimensional DFNM and a one-dimensional DCM, were developed on the basis of field observations and measurements and were compared to examine whether the DCM can simulate flow and transport in a 1-m3-sized fractured block represented by the DFNM. The AFC was incorporated into the DCM to examine whether the AFC improves the prediction capability for discrete flow behavior of the DFNM, which was used as a surrogate for unsaturated fractured rock. The simulations were also conducted to identify the most feasible design of laboratory experiments for demonstrating the AFC using a 1-m3-sized fractured block.

Of the two methods for demonstrating the AFC—varying the portion of the fracture network that contributes to water flow or varying the active interfacial area between fracture and matrix—we focused on the latter. This was accomplished by changing flow rates, causing different saturations and thus different interfacial areas, and by changing permeability at the fracture–matrix connections to increase interaction between the fracture and matrix.

Numerical simulations with the two conceptual models provided breakthrough curves of water and tracer with various plausible experimental conditions, including various injection rates and bottom flow boundary conditions for a proposed 1-m3-scale unsaturated fractured block experiment. The DFNM showed distinctive discrete flow patterns, such as stepwise increases in outflow rates and high water saturation along preferential flow pathways.

Simulated flow and transport using the DCM with varying flow rates and boundary conditions were compared with those with the DFNM and the AFC parameter {gamma} was estimated by inverse modeling. The results from the comparison can be summarized as follows:

  1. The DCM captures the flow pattern of the unsaturated fracture network better if the AFC is implemented. Earlier breakthrough of water at the bottom can be predicted using the AFC implementation. However, the estimated {gamma} value, which is expected to be a property of the fracture network only, varies with different injection rates.
  2. The free-drainage condition at the bottom boundary provides a better match to the DFNM data than the capillary-barrier boundary condition, but the calibrated {gamma} values were similar to {gamma} values for the capillary-barrier boundary condition when the injection rate is <10% of the saturated flow rate.
  3. Larger matrix permeability selected to mimic the larger interface area between fracture and matrix did not improve the prediction of water breakthrough with the DCM in the fracture network, indicating that the {gamma} was not sensitive to fracture–matrix interaction for expected laboratory conditions.
  4. The tracer transport in the fracture network was almost identically predicted by the DCM with or without the AFC when the injection rates were high (>10%). The AFC was not significant in the model where the fracture network acted as an averaged continuum at nearly saturated flow conditions, and the {gamma} values were not sensitive to the tracer transport at the balanced state of flow.
  5. For lower injection rates (<10%), the AFC helped capture the earlier breakthrough of tracer in the fracture network by the DCM, but the later tracer tailings resulting from slower diffusion of tracer coming out of the irregular matrix were not predicted by the AFC parameter calibration. Also, the estimated {gamma} values reached the upper limit, which suggested the AFC parameter was not a sufficient parameter to control the tracer transport in the fracture network.
  6. Numerical simulation provided a useful tool for the design of laboratory experiments that are planned for demonstration of the AFC. Based on the results from the numerical simulations, the AFC would be demonstrated by flow and transport tests using relatively low flow rates (<10%) when a 1-m3-sized block is subjected to laboratory experimentation.

Despite the difficulties of determining a {gamma} value from experimental data, performing a laboratory test with a large fractured block ({approx}1.3 m3) would provide a variety of data to add confidence in modeling flow and transport through unsaturated fractured rocks. Such a test, if performed correctly, could provide data allowing evaluation of the AFC in a three-dimensional fracture network. The use of modeling results would facilitate the design of a more focused and less resource intensive set of experiments.


   ACKNOWLEDGMENTS
 
The authors were helped greatly by Stefan Finsterle, Hui-Hai Liu, Grace Su, Boris Faybishenko, and Dan Hawkes in their careful reviews of this manuscript. This work was supported by the Director, Office of Civilian Radioactive Waste Management, U.S. Department of Energy, through Memorandum Purchase Order EA9013MC5X between Bechtel SAIC Company, LLC, and the Berkeley Lab. The support is provided to Berkeley Lab through the U.S. Department of Energy Contract No. DE-AC03-76SF00098.


   REFERENCES
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MODEL DEVELOPMENT
 NUMERICAL SIMULATIONS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 





This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Seol, Y.
Right arrow Articles by Ito, K.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Seol, Y.
Right arrow Articles by Ito, K.
GeoRef
Right arrow GeoRef Citation
Agricola
Right arrow Articles by Seol, Y.
Right arrow Articles by Ito, K.
Related Collections
Right arrow Flow
Right arrow Fractured Rock
Right arrow Vadose Zone Processes and Chemical Transport

HOME HELP FEEDBACK <