2. TASK PROGRESS

2.1 Active Tasks
2.2 Task Progress
2.2.1 Task 1.1.2: Data Updates
2.2.2 Task 1.2.1: 3D Hierarchical Fracture Model (HFM)
2.2.3 Task 1.2.2: Hierarchical Fracture Model Verification
2.2.4 Task 3.2.1: MS Windows 95 Analysis System
2.2.5 Task 3.2.2: Discrete Feature Analysis for TAGS Process
2.2.6 Task 4.1.2: Well Testing Data Acquisition
2.2.7 Task 4.2.1: Reservoir Model Implementation
2.2.8 Task 4.2.2: Reservoir Simulation
2.2.9 Task 5.1.2: WWW Site Updates
2.2.10 Task 5.2.1: Progress Reports
2.2.11 Task 5.2.2: Research Reports
2.2.12 Task 5.2.3: Presentations
2.2.13 Task 6: Management

2. TASK PROGRESS

2.1 Active Tasks

The following tasks were active during the quarter:

• Task 1.1.2 Data Updates
• Task 1.2.1 3D Hierarchical Fracture Model
• Task 1.2.2 3D Hierarchical Fracture Model Verification
• Task 3.2.1 MS Windows 95 Analysis System
• Task 3.2.2 Discrete Feature Analysis for TAGS Process
• Task 4.1.2 Well Testing Data Acquisition
• Task 4.2.1 Reservoir Model Implementation
• Task 4.2.2 Reservoir Simulation
• Task 5.1.2 WWW Site Updates
• Task 5.2.1 Progress Reports
• Task 5.2.3 Presentations
• Task 6 Management

2.2 Task Progress

This section describes progress during the quarter for each of the active tasks.

2.2.1 Task 1.1.2: Data Updates

During the quarter, Marathon Oil Company collected production and field testing data from the project study site and provided this data to Golder Associates. Data provided by Marathon included:

• Reservoir Connectivity (Tracer) Test Design and Results (Figures 2-1 through 2-8)
• Borehole fracture data (Figures 2-9 through 2-12)
• Core photography (Figure 2-13)

The data provided by Marathon during the quarter was assembled and is being incorporated into the project World-Wide-Web (WWW) site (Task 5.1.2).

2.2.2 Task 1.2.1: 3D Hierarchical Fracture Model (HFM)

During the quarter, MIT research under this task focused on completion and documentation of the 3D hierarchical fracture model (HFM) algorithms. Theoretical development during the quarter include:

• rotation and translation processes to orient and locate fractures with respect to spatial/geological information;
• algorithms to use StrataModel data in HFM realizations.

2.2.2.1 Rotation and Translation Algorithms

Algorithms for generation of fold-related fracture sets, implemented in the HFM were developed during the previous quarter (Dershowitz et al., 1997). The algorithms involve rotation of polygons according to specified relations to the local strike, dip, or normal vector of the fold surface at the polygon centers. The random rotation of polygons is used for cases when the local conditions modify the stress field so much that the most likely fracture orientation deviates from the most likely orientation defined by the general stress field. The 3D geometry of a fold is described by a cubic function. The coefficients of the cubic equation are derived through polynomial fit to the elevations of the contacts of folded strata in boreholes or on mapped outcrop exposures.

The orientations of fractures in folded strata are related to the varying strike and dip of the fold surface (Figure 2-14). The strike and dip of a fold surface at point P can be calculated as:

(2-1)

where n_f=(n_X, n_Y, n_Z) is the unit normal vector to the fold surface at point P (directional first derivatives of the function F(X,Y,Z)=0 that describes the fold), and __f and __f are the azimuth and latitude of n_f. To compare the orientation of a fracture-polygon to that of a fold, one can calculate the angle __between n_f and the unit normal vector of the fracture plane, n_p = (n_X,p, n_Y,p, n_Z,p). __can be calculated from the vector dot product:

(2-2)

To compare the strike of a fracture-polygon to the local strike of a fold, one can calculate the difference __strike between the azimuth of the polygon, __p, and the azimuth of the fold surface at the polygon center:

(2-3)

To compare the dip of a fracture-polygon to the local dip of a fold, one can calculate the difference __dip between the latitude of the polygon, __p, and the latitude of the fold surface at the polygon center:

(2-4)

Angles _, __strike, and __dip are checked against specified relationships of the fractures to the fold. For example, a small angle _ indicates that the fracture is subparallel to the fold surface. An angle __strike close to 90^o means that the fracture strikes approximately orthogonal to the local fold strike (i.e., parallel to the local slope of the fold surface). An angle __strike close to zero or 180^o means that the fracture strikes approximately parallel to the local fold strike (i.e., orthogonal to the local slope of the fold surface). A small angle __dip indicates that the fracture has approximately the same dip as the fold at that location.

If the orientation of a fracture-polygon does not conform to the specified relationships between fractures and fold, the polygon is rotated. For example, rotation by dip (so that the new fracture dip is subparallel to the local fold dip) is performed by assigning a new latitude angle __new to the polygon:

(2-5)

where ___dip is a small angle of allowed deviation of the fracture dip from the local dip of the fold. Rotation by strike is performed by assigning a new azimuth angle __new to the polygon:

(2-6)

where ___strike is a small angle of allowed deviation of the fracture strike from the local strike of the fold, or from the horizontal direction orthogonal to the local strike of the fold. The coefficient C_s is either C_s=0 or C_s=1.0, depending on whether the fracture is rotated to be concentric or radial to the fold, respectively. If necessary, a completely new orientation of the fracture can be generated (in a frame of reference where the local normal vector to the fold is assumed to be the mean polar direction). After calculating the new dip or/and strike, the fracture-polygon is rotated. Rotation of a polygon means that new 3D coordinates are calculated for every vertex.

2.2.2.2 HFM StrataModel Interface

An algorithm was developed during the quarter to read porosity values from the Yates Field StrataModel (Figure 2-15) and calculate the average porosity of the host rock that surrounds a numerically generated polygon. Then the polygon is retained as a fracture with probability P_f as a function of the average porosity:

(2-7)

where N is the number of geocells from the reservoir StrataModel intersected by the polygon, and n_i is the porosity of the i-th intersected cell.

In Tract 17, there is an additional influence of the rock lithology on the fracture intensity: "shales" are considered to be more ductile (hence less fractured) than the brittle dolomite between them. A numerical algorithm reads from the StrataModel not only porosity values, but also gamma ray (GR) data (indicator of shale content), and calculates the average porosity and the average GR of the host rock surrounding a given polygon. Then a polygon is retained as a fracture with probability P_f as a function of the average porosity and the average GR of the surrounding rock:

(2-8)

where N is the number of geocells from the reservoir StrataModel intersected by the polygon, n_i is the porosity and GR_i is the gamma ray response in the i-th intersected cell.

2.2.3 Task 1.2.2: Hierarchical Fracture Model Verification

Development of the HFM for Tract 49 was described in the previous quarterly report. During this quarter, MIT research extended the HFM verification simulation to Yates Field Tract 17.

Numerical modeling for Tracts 17 and 49 included the following:

• Inference of the model parameters from geology and field data, namely inference of number of fracture sets, fracture intensity P₃₂, fracture plane orientation PDF, fracture size variation, and relation of fractures to the field structure.
• Development of specific numerical algorithms.
• Results from numerical generations and comparison to field data.

The number of fracture sets, and their mean orientations, defined by the regional stresses, are the same in Tract 17 and Tract 49. Throughout the Yates field, field data consistently show that there are two major fracture sets. One set strikes to the northwest and the other one strikes to the northeast. The expected value for the strike of the northwest set is N50^oW, and for the northeast set is N40^oE (parallel to the two major regional depositional trends on the Central Basin Platform). Both fracture sets are composed of vertical or steeply dipping fractures. In the primary process of the numerical model, the mean pole direction (____) of fracture planes that belong to the first set is defined as (40^o, 90^o), i.e., orthogonal to a vertical plane striking N50^oW. The mean pole direction of fracture planes in the second set is defined as (-50^o, 90^o), i.e., orthogonal to a vertical plane striking N40^oE. A Fisher distribution with Fisher constant _=20 is used for the generation of fracture planes in each set. This spherical PDF (with the relatively high _=20) preferentially produces planes with poles at small angles with the specified mean pole orientation, as they occur in reality (mostly steep fractures in a narrow dip range).

At different locations in the Yates field, the predominant fracture strikes often deviate from the two regional directions since the majority of the fractures are also related to the curvature of the reservoir structure. The structure-defined variation of fracture orientations around the mean directions of the two major sets is relatively simple in Tract 49 and quite complex in Tract 17.

Fracture intensity P₃₂ in Tract 17 is determined from the mean spacing of significant fractures intersected by cores. Table 2-1 summarizes the spacing of significant fractures, identified in cores available from eight wells located inside or in the vicinity of Tract 17.

Table 2-1 shows the average spacing of intersections of significant fractures with the vertical wells in Tract 17 varies between 4.5 ft and 11.3 ft (1.4-3.4 m) where the rock porosity is n_ave<20%. The spacing increases where the rock porosity is n_ave>20% (see data for wells YU2416 and YU2509). When the smaller significant fractures are also considered, the average spacing is between 4.5 ft and 7.5 ft (1.4-2.3 m). There may be a tendency of increasing spacing with depth, suggested by the longest available continuous core (318 ft; 97 m) at well YU2433. Forty vertical fracture intersections in Tract 17 cores have heights from 1 ft to 44 ft (0.3-13 m), all but five having heights less than 10 ft (3 m).

Table 2-1 Tract 17 Fracture Spacing

Table 2-1 is based only on "significant fractures" identified in core. The values reported in Table 2-1 are:

• n_ave: average core porosity;
• N: number of large significant fractures identified at the given well;
• s_ave: average spacing of large significant fractures;
• SD: standard deviation of spacing between large significant fractures (given only where N>5).

Compared to Tract 49, there are more vertical or nearly vertical fractures in Tract 17. Predominantly vertical or very steep dips of the fractures in Tract 17 are also observed in the log analysis of the only horizontal well in the study area: well YU17D5. According to Tinker & Mruk (1995), "the quality of the data in this well and the confidence in the fracture picks are very high". The average spacing of 671 fractures intersections in the 1550 ft (472 m) long horizontal well is s_ave= 2.31 ft (0.7 m). Steep and nearly vertical fractures are more easily intersected by the horizontal well than by the vertical wells.

In Tract 17 the relationship of the fractures to the field anticlinal structure is more complex than that in Tract 49 (Dershowitz et al. 1997). Tract 17 is on the peri-anticline of the reservoir dome: a zone of transition from incompactible grainstones on the east side toward compacted mudstones to the west. In Tract 17, a northwest striking regional fracture set is sub-parallel to the major hinge of the peri-anticline and to a possible set of field-scale drape folds / faults. A northeast striking regional set is sub-parallel to the minor hinge of the peri-anticline (defining the cross-curvature of the fold) and to a second set of drape folds / faults.

The following parameters of the 3D fracture system model were derived for HFM simulations in Tract 17:

• Orientations of two fracture sets (inferred from regional geology): Vertical planes with Fisher distribution dispersion _=20, striking N50^oW and N40^oE. Fisher distribution dispersion _ was selected to produces close clustering of plane pole orientations around the mean pole.
• Modeling volumes:
  • Vertical logged wells: YU1711 and YU1755 inside Tract 17, and YU2511 in the vicinity of Tract 17. The modeling volume for a vertical well is defined by the following boundaries: a horizontal datum surface at 700 ft above sea level; four vertical planes (two striking south-north and two striking east-west), each at a horizontal distance of 300 ft from the vertical well; a quadratic upper surface, parallel to the San Andres top and located 100 ft above it. The coefficients of the quadratic function are calculated through a polynomial fit with MATLAB to the elevations of the San Andres top at the wells in the Tract 17 area.
  • Horizontal logged well: YU17D5 in Tract 17 (north-south trend). The modeling volume for the horizontal well is defined by the following boundaries: a horizontal datum surface at 800 ft above sea level (about 200 ft below the elevation of the well); two vertical planes, striking south-north, at a horizontal distance of 200 ft from the well; two vertical planes, striking east-west, at a horizontal distance of 800 ft from the midpoint of the well (i.e., 50 ft from the two ends of the well); a quadratic upper surface, parallel to the San Andres top and located 100 ft below it (i.e., at least 200 ft above the elevation of the well). The coefficients of the quadratic function are calculated through a polynomial fit with MATLAB to the elevations of the San Andres top.
• Fracture intensity (inferred from field core data): total P₃₂=0.8 [ft²/ft³], established through simulations to produce expected spacing in vertical boreholes E[s]=7 ft. The intensities of two fracture sets contribute to the total P₃₂: a major set with assumed intensity P_32,1=0.6 [ft^-1], and a minor set with intensity P_32,2=0.2 [ft^-1]. The N50^oW regional set is major in the vicinity of well YU1711 which is located on the crest of the major northwest striking anticlinal axis. The N40^oE regional set is major in the vicinity of the other three wells: YU17D5, YU1755, and YU2511, which are located in more peripheral areas of the peri-anticline.
• Fracture size distribution (approximate procedure, since field sampling of fracture sizes is not available): the mean of the equivalent radius is assumed =30 ft (related to the thickness of mechanical units in the San Andres formation). Marking by relative size eliminates the polygons with extremely large sizes, and produces fracture areas that have a standard deviation and median , where is the expected area of polygons marked as fractured. 65% of the produced fractures have areas smaller than the average area.
• Effect of porosity (approximate procedure since the exact effect of high porosity on fracture intensity has not been studied yet): fracture intensity in areas with porosity n>20% is assumed half of that where n<20%, i.e., P₃₂|n>20%=50% P₃₂=0.5(0.8)=0.4 [ft^-1]. A fracture is discarded with probability P_f=0.5 if the average porosity of the dolomite matrix around the fracture is n_ave>20%. This method results in discarding fractures with relatively small sizes which are the most likely to be located entirely in an area of high porosity. Thus in high porosity dolomite only relatively large fractures are retained which results in larger spacing of fractures intersected by the wells.
• Effect of shale (approximate procedure since the exact effect of high shale content on fracture intensity has not been studied yet): fracture intensity in areas with high shale content is assumed to be half of the fracture intensity in clean dolomite, i.e., P_{32, shale} =50%P₃₂=0.5(0.8)=0.4 [ft^-1]. Rock with high shale content (>10%) is defined by porosity n_ave<20%, and the gamma ray response GR_ave>40 API. Thus a fracture is discarded with probability P_f=0.5 if the average porosity of the dolomite matrix around the fracture is n_ave<20% and the average gamma ray response is GR_ave>40 API. This methods results in discarding fractures with relatively small sizes, since they are the most likely to be located entirely in the relatively thin shale layers. Thus the clean dolomite in Tract 17 includes numerous small and fewer large fractures, but only the relatively large ones cut across the thin shale layers.
• Relationship to the anticlinal structure (inferred from local geology): in Tract 17, where the slope of field anticline is very gentle, the relationship of fractures to the local slope of the structure is not as pronounced as it is in Tract 49. The shape of the field anticline is approximated by a cubic surface, fit to the shape of the Seven Rivers M horizon in Tract 17. The coefficients of the cubic surface are calculated through a polynomial fit with MATLAB to the elevations of the Seven Rivers M horizon in the wells in Tract 17. The relation of fracture orientation to the reservoir structure is more complex than it is in Tract 49: predominant fracture orientations in Tract 17 are related to the shape of the entire peri-anticline rather than only to the local slope of the structure. This relationship is accounted for by choosing which one of the regional sets (related to the shape of the peri-anticline) is dominant in the vicinity of different wells. Only a small percent of fractures in Tract 17 are related to the local strike of the fold structure at their centers. In the numerical generation, such fractures are rotated to be predominantly parallel to the slope of the cubic surface near wells YU1711, YU1755, and YU2511, and predominantly orthogonal to the slope (i.e., parallel to the structure itself) in the vicinity of the horizontal well YU17D5).

Table 2-2 summarizes results from numerical simulations of the fracture system around the three logged wells in Tract 17 (YU1711, YU1755, and YU17D5). These results reflect the last step of the simulations that established the correct fracture intensity P₃₂ of the system in Tract 17. The spacing of fracture intersections in the simulated horizontal well YU17D5 (shaded column in Table 2-2), which is essentially equal to the known actual spacing of 2.31 ft, confirms that the fracture intensity P₃₂ used in the simulations (based on fracture intersections with vertical wells) has been assumed correctly.

The horizontal sections in Figures 2-16a, 2-17a, and 2-18a illustrate the shape of the field anticlinal structure in the vicinity of the logged wells inside and near Tract 17. Fracture traces on hypothetical horizontal outcrops in Figures 2-16b, 2-17b, and 2-18b show how the numerically generated fracture system in Tract 17 relates to the field structure in the vicinity of the three logged wells in Tract 17 (YU1711, YU1755, and YU17D5). Figures 2-19, 2-20, and 2-20 show rosette diagrams of fracture strikes in simulated boreholes at wells YU1711, YU1755, and YU17D5. The strikes of numerically generated fractures are compared to the strikes of fractures, identified on the log profiles of the wells in Tract 17. Figure 2-21 illustrates the normalized dip distribution of numerically generated fractures compared to the normalized dip distribution of significant fractures identified in cores in Tract 17. The figure shows the fracture dip distribution in a simulation of borehole intersections with wells YU1711, YU1755, and YU17D5 (corresponding to the shaded row in Table 2-2).

Table 2-3 shows results from simulations in which the effect of high porosity and shale content on fracture intensity was also considered. Numerically generated fractures which lie entirely in regions with either average porosity n_20% or with average shale content of more than 10% are discarded with probability P_f=0.5. Thus, in the high-porosity dolomite the fracture system includes only relatively large fractures, whereas in dolomite with lower porosity both large and smaller fractures exist. Only the relatively large fractures cut across the thin layers with high shale content, whereas the intensity of small fractures in the shales is much lower than it is in the clean dolomite. The simulations in Table 2-3 illustrates the capability of the model to reproduce fracture intensity as a function of rock properties. Since the exact effect of the porosity and shale content of San Andres dolomite on the fracture intensity in the formation has not been studied yet, more precise simulations cannot be done at this stage of the Yates case study.

Table 2-2 HFM Numerical simulations of Tract 17 (near wells YU1711, YU1755, and YU17D5)

The values reported in Table 2-2 are:

• P₃₂: fracture intensity defined as cumulative fracture area per unit volume;
• : expected value of fracture equivalent radius;
• N: total number of fractures generated in a simulation;
• N_b: fractures intersected by the simulated borehole;
• s_ave and __s: average value and standard deviation of fracture spacing in the simulated borehole.

Dip distribution of the fractures in the shaded row is shown in Figure 2-22. The shaded column represents a test of P₃₂ (assumed on the basis of spacing of fracture intersections with vertical wells): compare to the actual spacing of 2.31 ft determined from log analysis of the horizontal well YU17D5.

Table 2-3 Tract 17 HFM Simulations with Intensity as a Function of Lithology

Table 2-3 reports the same parameters as Table 2-2. After generation of the fracture system, P₃₂ in dolomite of high porosity (n>20%) is reduced to half of P₃₂ in lower-porosity dolomite, and P₃₂ in "shales" (shale content>10%) is reduced to half of P₃₂ in clean dolomite.

2.2.4 Task 3.2.1: MS Windows 95 Analysis System

During the quarter, major effort focused on development of the spatial data analysis system. The name of the software under development was changed from "Fractal" to "Spatial" to better describe the type of spatial analyses being performed. The algorithms for Spatial are described in Dershowitz et al. (1997). The user interface developed for Spatial is described in Section 2.2.4.1.

2.2.4.1 Software Development Status

The status of active development of software is summarized in Table 2-4. All Windows 95 software are undergoing further improvements and should be considered Beta 1.0 versions.

Table 2-4 MS Windows 95 Analysis System

2.2.4.2 Spatial 1.0 User Interface

As a result of the increasing interest in Web-based software, the project team decided to implement spatial data analysis as a Java 1.1 application, which will ultimately be executable as a server application from the project web site. As a server application, users would not need to download and install Spatial in order to evaluate the software. Spatial 1.0 was coded using Java 1.1, and is currently compiled as a Windows 95 application.

Spatial 1.0 spatial data analysis software is designed to facilitate spatial analysis lineament map files, using the algorithms described in previous project reports (Figure 2-23). The user interface for Spatial 1.0 is illustrated in Figures 2-24 through 2-28. Spatial 1.0 provides the following four features for spatial analysis of lineament data:

• data gridding (Figure 2-24a)
• contour plotting of intensity, size, and termination mode to facilitate visual evaluation of spatial trends (Figure 2-24b)
• trend analysis of intensity, size, termination mode, and orientation to facilitate analysis of non-stationarity (Figure 2-25)
• dependency analysis for evaluation of trends with respect to location of primary sets (Figure 2-26)
• correlation analysis for comparison of spatial trends for different sets (Figure 2-27)

The flow chart for Spatial is illustrated in Figure 2-28, based on that for the previous software prototype, Fractal 1.0. The algorithm is described in Dershowitz et al. (1997).

2.2.4.3 Spatial User Interface

To execute Spatial, select the Spatial icon, and double click with the mouse. The general procedure for analysis is summarized in Table 2-5. The operation of the individual Spatial menu items are described in the next section. Navigation through Spatial is done using Microsoft Windows mouse conventions. In general, the left hand mouse button is used for making selections.

Table 2-5 Spatial Lineament Analysis Sequence

2.2.4.4 Command Summary

2.2.4.4.1 File Menu

New: Begin a new analysis, closing the current analysis

Open: Open a fracture lineament (.DAB) files.

Close: Close the current active analysis and all windows related to that analysis.

Exit: Exit Spatial.

2.2.4.4.2 Analysis Menu

Grid Size: Specify the size of the grid to be overlain on the lineament map for calculations.

Projection Angle: Specify the angle (clockwise from +Y to be used in calculation spatial trends.

Primary Set: Specify the primary set and dependent sets to be used in calculating spatial trends with distance from a primary set.

Calculate: Calculate values on grids. View options are not available until calculate has been selected.

2.2.4.4.3 View Menu

Contour Grid: Display color contour plots of intensity (P₂₁ (m/m²), length, and orientation. Based on this display, determine what direction should be used in looking for spatial trends (non-stationarity).

Spatial Trend: Display the variation in intensity (P₂₁ (m/m²), length, and orientation projected along a line at the specified projection angle. Interactively modify the projection angle using the compass provided.

Dependent Sets: Display the variation in intensity (P₂₁ (m/m²), length, and orientation with distance from fractures of the "primary set" along a line at the specified projection angle. Interactively modify the projection angle using the compass provided.

Set Correlations: Display the variation in intensity (P₂₁ (m/m²), length, and orientation projected along a line at the specified projection angle for fractures of primary and dependent sets. Interactively modify the projection angle using the compass provided.

2.2.4.4.4 Help Menu

Help Index: Not currently implemented.

Using Help: Not currently implemented.

About: Spatial QA, copyright, and license information.

2.2.5 Task 3.2.2: Discrete Feature Analysis for TAGS Process

The goal of this task is to develop the procedures for discrete feature modeling of the TAGS process appropriate for the Yates Field. During the quarter, work was initiated on development of the following procedures:

• procedures to derive appropriate DFN model parameters using the MS Windows Analysis System
• tracer transport modeling approaches for the study of reservoir connectivity as described by in situ tracer testing
• compartmentalization analysis using FraCluster
• procedures to apply heat particle tracking simulation to model TAGS
• procedures for derivation of parameters of ECLIPSE continuum simulation

These procedures are described below in Section 2.2.6 and 2.2.7, where the application of these procedures to Yates reservoir simulation is developed. Preliminary connectivity analyses were carried out for the in situ tracer testing carried out by Marathon to support the TAGS process, and a demonstration FraCluster calculation was carried out to understand Yates field compartmentalization behavior.

2.2.6 Task 4.1.2: Well Testing Data Acquisition

During the quarter, Marathon collected and processed well test and hydraulic response data, and provided the data for posting on the WWW server. Well responses provided by Marathon during the quarter are illustrated in Figures 2-29, 2-30, and 2-31.

2.2.7 Task 4.2.1: Reservoir Model Implementation

During this quarter, progress was made in preliminary development of the DFN model of the Yates field to be used in reservoir simulation. This model is based on the FracMan DFN approach, which has not yet incorporated the HFM model.

This section described the data analysis used to construct the preliminary Tract 17 DFN model for TAGS support simulations. This model is based on a quantitative analysis of fracture orientation and size distributions and intensity, and is complementary to the data analyses for the HFM as described in section 2.2.3.

Derived model parameters are summarized in Table 2-6. Visualizations of the preliminary DFN model is are provided in Figures 2-32, 2-33, and 2-34.

Table 2-6 Preliminary DFN Model for TAGS Support Simulations (Tract 17)

2.2.7.1 Stochastic and Deterministic Discrete Features

The preliminary DFN model for Yates simulations combines deterministic faults with stochastic discrete fractures. Large scale faults are located deterministically based on seismic surveys. Through seismic interpretation, the location of twelve major faults has been determined inside Tract 17. These are large scale features between 400 and 1500 m long, and have a large horizontal to vertical aspect ratio of between 2:1 and 10:1. These faults form two distinct subsets trending NW-SE and NE-SW and are assumed to be vertical. These features are represented deterministically in the fracture network model as their location and geometry is known.

The stochastic fractures combine fractures from all discrete features recognized in Tract 17 borehole logs. The geological and hydrogeological properties of these features is described below. These features are represented stochastically in the discrete fracture network model as the extent of individual features cannot be determined purely from borehole logs. This stochastic fracture set may also be further divided into subsets based on fracture orientation and this issue is discussed in the following section.

2.2.7.2 Orientation Distribution

The orientations of fractures logged in four wells in and near Tract 17 were analyzed to determine whether the data could be defined in terms of sets with specified orientation distributions, or whether a bootstrap simulation technique (Dershowitz et al, 1996) would be required to match the observed orientation patterns.

FMI orientation data from Wells YU1711, YU1755, YU17D5 and YU2511 were analyzed using the ISIS interactive set identification system. All of the wells are vertical with the exception of YU17D5, which is horizontal with a trend of approximately north-south. Stereoplots for each of these four wells are shown in Figure 2-35a-d, "Terzaghi" corrected for the bias due to the angle of the fracture relatively to the orientation of the wellbore. The three vertical wells show a similar pattern of fracture orientations, which are most clearly delineated in YU1711 (Figure 2-35a). There appear to be up to 6 sets whose prominence may been locally enhanced by the folding of the reservoir in the Tract 17 area. There are two sets that strike to the northwest, one dipping steeply to the northeast, the other dipping steeply to the southwest. There is another set of fractures striking to the northeast and dipping to the northwest, and a much less prominent set with the same strike, but dipping to the southeast. There may also be a significant set that strikes east-west and dips to the north, and another that strikes north-south and dips to the west. The horizontal well has intersected a prominent set that strikes east west and is nearly vertical, and a less significant set that strikes to the northeast and is also vertical.

One of the difficulties in interpreting these stereoplots is to decide whether to combine individual concentrations of poles into a single set, or to treat them separately. For example, there are three pole concentrations shown in the north-northeast portion of the stereoplot for well YU1711 (Figure 2-35a). At the highest concentration levels, there are three sets. At slightly lower concentrations, all three could be grouped together into a band that follows a stereoplot arc. The coalescence of distinct groups of poles at high concentrations into girdles at medium concentrations is a prominent feature of the stereoplots for fractures from the three vertical wells. Regional production data and flexure maps (Dershowitz et al, 1997) do not provide much evidence for an east-west striking set, which seems to form a significant concentration in all four wells, but there is good evidence from production data and flexure maps for northwesterly- and northeasterly striking sets.

The relatively greater importance of the east-west, northerly-dipping set in wells YU1755, YU1711 and YU17D5 may reflect the local folding of the San Andres reservoir. The structural contour map for the top of the San Andres in Tract 17 (Figure 2-36) shows that this set is approximately parallel to the structural contours and perpendicular to the contour surface at wells YU1711, YU1755 and YU17D5. At YU2511, where the contours strike northerly and dip to the west, the northwest, westerly dipping set has the highest concentrations of poles. Thus, there appear to be several joint sets whose relative importance may reflect local warping of the reservoir.

In order to evaluate whether concentrations should be lumped together or maintained as distinct sets, a series of detailed analyses were carried out on the data from each well. All data is taken from the San Andres unit. The results are summarized in Table 2-7:

Table 2-7 Orientation Analysis for Yates Field Unit

Table 2-7 shows that there is considerable variability in the mean pole trend and plunge for concentrations that appear to belong to the same set. For example, Set 5 in YU1711 and Set 2 in YU2511 differ in trend by almost 20 degrees. The east-west set varies over a similar range. The other sets, which are probably not enhanced by local folding also show variability over more than 20 degrees in trend. Based on Table 2-7, it seems impractical to fit a separate model for each set since the variability between wells would be on the order of 20 or more degrees. As a result, it was determined that numerical simulation of orientation distributions for the preliminary DFN model should be based on the bootstrap algorithm (Dershowitz et al., 1996).

The orientation bootstrap is illustrated in Figure 2-37. The bootstrap starts with the Terzaghi corrected data from each of the boreholes. Then, in each Monte Carlo simulation a value of orientation is selected from the Terzaghi corrected data, with a uniform distribution of probability. An orientation realization is then selected from a Fisher distribution with user specified dispersion _ with mean pole at the selected data point. The value of Fisher dispersion _ used determines how precisely the Monte Carlo population matches the Terzaghi corrected data. For high ________, the simulated orientations are almost identical to the data. For low _ (<10>, the simulated orientation distributions have a large degree of scatter about the measured values.

Figure 2-38a presents a stereoplot of Terzaghi corrected orientation data from the three vertical wells YU1711, YU1755 and YU2511. The data from the horizontal well has not been used since the northerly-striking sets are very difficult to intersect with such a well trajectory. Note that the fracture poles form four approximate girdles, with the most prominent ones being in the southwest, northwest and northeast quadrants. Figures 2-38b-d present bootstraps from this data set with _=5, 30, and 200. Based on visual comparison of stereoplots with a range of Fisher dispersion parameters, a value of Fisher dispersion _=30.0 was selected. This value will be used in preliminary simulations.

2.2.7.3 Fracture Size Distribution

A size distribution assumption was derived from the FMS and FMI data in the four wells YU1711, YU1755, YU17D5 and YU2511 using the method of LaPointe et al (1993). La Pointe et al. (1993) showed that the size distribution of a fracture population could be estimated from the relative proportion of fractures detected as a function of the number of pads on which the fracture was imaged (Figure 2-39). For a 4-pad FMS tool, this means that each fracture will be imaged on 1, 2, 3 or 4 pads. The larger the fractures, the greater the probability that they will be imaged on all four pads of the tool. For an 8-pad FMI tool, each fracture will be imaged on from 1 to 8 pads. An eight-pad tool was used for the vertical wells YU1711, YU1755 and YU2511, and a 4-pad tool was used for the horizontal well YU17D5.

In order to compare an 8-pad image log with a 4-pad log, it was necessary to combine adjacent pads on the 8-pad tool. This implies that the number of fractures imaged on 1 or 2 pads are combined, 3 and 4 are combined, 5 and 6 are combined and 7 and 8 are combined. This is not an exact correction, since the azimuthal coverage of a single pad from a 4-pad tool differs from the coverage of two adjacent pads on an 8-pad tool, but it is still useful.

Figure 2-40 shows the results for each of the four wells, for the mean and the median of the vertical wells taken together, and for all wells regardless of plunge. All of the intersection percentages as a function of the number of pads are very similar for all well and well groupings with the exception of the horizontal well, YU17D5, which has a higher percentage of 7 & 8 pad fracture images. The pad intersection percentages are summarized in Tables 2-8 and 2-9.

Table 2-8 Pad Percentages for FMI Log Data

Table 2-9 Pad Percentage Statistics for FMI Log Data

Tables 2-8 shows that wells YU1711 and YU2511 are very similar; YU1755 has a greater proportion of fractures that intersect all the pads of the tool. This could be due to:

(1) the westerly-striking, subvertical fractures that account for almost all of the fractures imaged in the horizontal well are larger than the other fracture sets;

(2) the azimuthal coverage of the 4-pad tool in the horizontal well differs from the azimuthal coverage of the 8-pad tools used in the three vertical wells;

(3) fractures in the vicinity of YU17D5 are larger than elsewhere.

Any or all of these factors may explain the difference in intersection percentages. However, considering the variability among the percentages for the three vertical wells, it probably due to spatial variability in the fracture sizes rather than to other causes. For purposes of estimating fracture sizes, the average for all wells was selected as the pad percentages for matching.

The fracture radius distribution was determined by comparing the measured pad percentage statistics from FMI logs against simulated pad percentages for assumed fracture size distributions (Figure 2-36). Several distributions were tested to see if they could match the observed vertical well pad intersection data. The best match was achieved with a truncated power law distribution with a fractal dimension of 1.81 and a minimum radius R_min of 1.0 ft. The distribution parameters and the pad intersection percentages for 25 realizations of the DFN model are shown in Table 2-10.

Table 2-10 Fitted Distribution of Pad Intersections (Power Law, D=1.81, R_min=1.0ft)

Note that this analysis was carried out for all fractures identifiable in FMI logs. It is likely that the fractures that play a significant role in interwell fracture connectivity represent the largest fractures only; small fractures probably play a very insignificant role at this scale. Thus the size should be adjusted to reflect the fact that the hydraulically significant fractures have a much larger minimum size cutoff. This issue was addressed below.

2.2.7.4 Fracture Intensity

The three dimensional stereological measure of fracture intensity P₃₂ is used in DFN modeling. P₃₂ is defined as the total fracture surface area per volume of rock containing the fractures. P₃₂. has the units of length^-1. For a fracture population with a defined orientation distribution and size distribution, the relation between P₃₂ and the number of fractures per unit length that would be intersected by a wellbore of a specific diameter and orientation is a multiplicative constant. This constant can be determined by creating a DFN model with any value of P₃₂, placing wells with the same orientations and diameters as the four Yates wells, and computing the number of fractures intersected per unit length. Table 2-11 and Figure 2-37 summarize the number of fractures per unit length (P₁₀) for the four wells as a function of the number of pads on which the fracture is seen.

Table 2-11 Fracture Intensity from FMI Data

Table 2-11 shows that the number of fractures per unit length (P₁₀) varies among the four wells and as a function of how many pads the fracture is imaged on. For all fractures, the number of fractures per foot varies from a low of 0.22 for YU1755 to a high of 2.2 for the horizontal well.

A range of interpretations of fracture intensity P₁₀ from Tract 10 data is summarized in Table 2-12. The mean wellbore length-weighted intensity for vertical wells is 0.70 fractures/ft. It is interesting to note that for the vertical wells, approximately half of the overall fracture intensity is made up of the biggest fractures, those that are imaged on 7 or 8 pads. If the conductive fractures are the biggest fractures, then this indicates that the conductive P₁₀ should be no greater than 0.35 fractures/ft. Evaluation of static spinner logs from Tracts 17, 49 and neighboring tracts show a spacing of flow anomalies on the order of from 40 ft to 200 ft, corresponding to intensity P₁₀ of 0.025 ft^-1 to .005 ft^-1.

Table 2-12 Fracture Intensity from FMI Logs and Spinner Surveys

The linear fracture intensity measure P₁₀ can be converted to the volumetric intensity measure P₃₂ by calculating the intensity P₁₀ (spacing) for a well in a DFN simulation with known P₃₂. Using a bootstrap fracture orientation distribution and a vertical well, the ratio P₃₂ / P₁₀ was determined to be 2.03. Thus an intensity P₁₀ of 0.70 ft^-1 corresponds to P₃₂ of 1.42 ft^-1 .

For reservoir scale simulation, the intensity required is the conductive intensity, rather than the geological intensity. One approach to relate geological and conductive intensities is by assuming that fractures below a specified size threshold are non-conductive, and using a truncated radius distribution with the corresponding intensity. The percentage of fractures below the radius truncation limit could then be used to calculate the ratio of conductive fracture intensity to geological fracture intensity.

Figure 2-41 illustrates the relationship between the radius distribution cutoff and the fracture intensity P_32c above that cutoff for a power-law radius distribution with dimension D=1.81 and minimum 1 ft. For a radius cutoff of 10m, P₃₂ would be reduced by 17.5%, and for a radius cutoff of 100m, P₃₂ would be reduced by 41.5%.

For simulation purposes, we have assumed that the intensity could adequately be represented by the 30 ft spacing features (P₁₀ = 0.033, P_32C = 0.07 ft^-1 ). A minimum radius threshold of 87 ft. was used to facilitate simulation by reducing the number of small fractures simulated at reservoir scale, while maintaining the fracture intensity P₃₂ . A more rigorous analysis of intensity and fracture radius threshold will be carried out in future quarters.

2.2.7.5 Transmissivity Distribution

The fracture transmissivity distribution can be derived from analysis of drill stem tests (DST's) or spinner surveys using Flare (Dershowitz et al, 1997). Alternatively, fracture transmissivity can be obtained by calibration of large scale flow and tracer testing and production histories. These analyses will be carried out during the next quarter. For the current simulations, a fracture transmissivity distribution was assumed based on the knowledge that the increase in head at the tracer injection well is small for the given injection rate. Preliminary scoping calculations then indicate that the fracture transmissivity should then be in the range 1 x10^-2 to 1 x 10^-4 m²/s for the fracture intensity assumed in the discrete fracture network model. Additionally, it is anticipated that the deterministic faults will on average be more hydraulically significant than the stochastic fractures. The base case model considers the deterministic faults to have a mean log transmissivity of -3.5 log m²/s and standard deviation of 1 log m²/s (corresponding to 4.5 x 10^-3 m²/s and 6.3 x 10^-2 m²/s on an arithmetic scale), while the stochastic fractures have a mean log transmissivity of -4 log m²/s and standard deviation of 1 log m²/s (corresponding to 1.4 x 10^-3 m²/s and 2.0 x 10^-2 m²/s on an arithmetic scale). There is considerable uncertainty in these values of transmissivity, hence a sensitivity study has been conducted to investigate this issue (Section 2.2.8).

2.2.7.6 Aperture

Fracture aperture is a key parameter for determining travel times within a discrete fracture network. Such calculations are made for the tracer test simulations described in Section 2.2.7.3 and aperture has been assigned by correlation with transmissivity. A relationship known as the `cubic law' is often used to describe this correlation, which has the general form of:

a = b * T ^0.33

where:

a = aperture (m)

b = constant appropriate to the flow medium and nature of fracture network

T = transmissivity (m²/s)

For the purposes of the preliminary tracer test simulations presented in Section 2.2.7.3, this equation takes the form of a = 0.011 *T ^0.33. This is a first order approximation of aperture and may be subject to calibration in future flow and transport modelling.

2.2.8 Task 4.2.2: Reservoir Simulation

In this task, reservoir simulations are carried out to demonstrate the practical application of the technologies developed during the project. During this quarter, work was initiated on development of the following model simulations:

• compartmentalization analysis using FraCluster
• tracer transport modeling approaches for the study of reservoir connectivity as described by in situ tracer testing

During the next quarter, reservoir modeling will be initiated for the Yates field:

• discrete feature network heat transport modeling for TAGS
• procedures for derivation of parameters of ECLIPSE continuum simulation

2.2.8.1 Compartmentalization Analysis

Compartmentalization directly affects oil production at the Yates Field by determining the amount of reservoir that a well can access as a function of oil column thickness. This section describes the oil production issues to be addressed through this analysis, and presents the preliminary analyses carried out during the quarter.

2.2.8.1.1 Background

A key component of the compartmentalization of the Yates Field is due to the effect of the oil column thickness. Changes to the oil column thickness and the location of the producing column relative to the well completions will result in changes in compartmentalization over time. The compartmentalization analysis developed in the following section addresses the relationship between oil column thickness and the producible oil compartment volume.

At Yates, FracCluster compartmentalization analysis will be used to study optimization of production rates and well spacings and to assess the benefits of converting from vertical to horizontal completions as the oil column thins. Figure 2-42 illustrates the dynamic compartmentalization process as the oil column thins.

Figure 2-42a represents a static oil column underlain by water and overlain by gas. All fluid contacts are approximately horizontal. Figure 2-42b shows that same system with active deep and shallow completions. Once production (flow) is initiated, the system behaves as dual porosity medium. The matrix which stores hydrocarbons also serves as a barrier to fracture network flow since matrix permeabilty is much lower than the fracture permeability. All fracture connections are fluid saturated (water, oil, or gas). During the initial period of production before there is any appreciable thinning of the oil column, contact lowering eliminates water connections while maintaining oil column thickness for cross-system oil connection. Gas connections continue to increase but these are not considered significant since there is no dynamic driving force to result in gas flow across the connections (only slow pressuring).

Continued co-production lowers fluid contacts and eventually leads to thinning of the oil column. Optimum contact lowering produces the maximum rate of oil movement from the matrix into the fracture system. If there is sufficient fracture connection within the oil column, the fracture network transports the oil across the network to the production well. As long as the oil column is thick, the oil production rate is a function of both column depletion and the oil transfer efficiency between the matrix and the fractures. Oil production at rates above the matrix oil mobilization rate can continue only as long as sufficient fracture network permeability within the oil column.

Eventually, the oil column will begin to thin with increasing oil production. Even though fractures might form geometric pathways that encompass both the oil and the gas columns, only those pathways entirely within the oil column can contribute oil to a producing well. The volume of matrix accessible through the fracture system will begin to decrease, as shown in Figure 2-42c. At some point, horizontal completions may be required to access sufficient volumes of oil stored in the matrix that would be difficult to capture through vertical completions.

Spill points provide an indication of the scale of reservoir compartmentalization. A spill point is a fluid contact location that is necessary to maintain horizontal movement of a fluid of interest. Figure 2-42c shows several oil spill points. The oil-gas contact would have to be higher than the two upper locations in order for the oil to move outwards through the fracture system towards wells located at the edges of the diagram. Alternatively, the oil-water contact would have to be lowered to the location of the lower oil spill point for oil to move outwards to the edges. Flow meter data can indicate something about the vertical spacing of spill points. At the Yates Field, observation wells in highly fractured areas over 600' from producing wells indicate effective spill points are about 25'-30' apart. In less perfectly connected areas, these spill points are greater than 50' apart. Within some vertical wellbores, connections may average 1-2 per 200' vertical section. The differences in the separation of effective spill points may be largely a function of equivalent radius; the network within the oil column elevation slice has a high effective radius.

Compartmentalization analysis can also be used to target areas for depletion processing. High-rate water disposal (withdrawal), which may double the rate of water-oil contact lowering, will rapidly eliminate water connections while oil connections grow, leading to a moderate oil production increase. Oil rate growth is planned while both fluid contacts are encouraged to "flatten" on a field scale. (Water shutoffs result in both increased water energy and connection.) Efforts to selectively withdraw oil from poorly connected areas will result in associated fluids (gas or water) production unless neither the gas nor water is well connected to the area of withdrawal. Low connection in all three phases indicates an area for depletion processing. Pressure depletion will result in hydrocarbon expansion for oil expulsion from the matrix into the fractures.

2.2.8.1.2 Preliminary Compartmentalization Analysis

A preliminary compartmentalization analysis was carried out using the FraCluster technologies developed during the last project year (Dershowitz et al, 1997). These simulations use compartmentalization analysis to investigate whether the thinning oil column at Yates will have a major affect on the drainage area of injection and production wells. This analysis and future extensions will support Marathon's Tract 17 well recompletion decisions, which will form the basis for Task 4.3, Technology Evaluation.

The compartmentalization analysis used the stochastic fractures from the DFN model derived in Section 2.2.7 above. Parameters are listed in Table 2-6. Visualizations are provided in Figures 2-32, 2-33, and 2-34. The compartmentalization analysis utilized a conductive fracture intensity P₃₂ of 0.05 ft^-1, corresponding to intensity P₁₀-ft^-1 of 0.025 and a mean conductive fracture spacing 40 ft.

Compartmentalization analysis was carried out for oil column thicknesses of 12.5, 25, 35, 50, 75 and 100 ft. These values correspond to current and expected future oil column thicknesses for Tract 17. Ten realizations were carried out for each oil column thickness. In each realization, FraCluster algorithms were used to calculate:

• the distribution of compartment volume (ft³);
• the distribution of compartment horizontal cross-sectional area (ft²).

The distribution of compartment volumes for a single realization of the DFN model with an oil column thickness of 35 ft is presented in Figure 2-43. The distribution of compartment volumes indicates that at this column thickness, well completions would have a high probability of intersecting sufficiently large oil compartments. However, even at this oil column thickness, a significant percentage of the oil column is contained in compartments of less than 10⁶ ft³ .

Figure 2-44 presents results of ten compartmentalization analysis realizations for an oil column thickness of 35 ft. Each realization is generated using the DFN model of Table 2-6, with a different random seed for the Monte Carlo simulation. As a result, this variation in compartment volume between realizations could be considered to represent the variation in compartment volume between different locations in Tract 17. The mean compartment volume generally varies between realizations by a factor of 2, but in one case varies by over an order of magnitude.

As expected, the mean compartment volume varies directly with the oil column thickness, as shown in Figure 2-45. The compartmental volume which could be accessed by a single vertical well decreases about 1.5 orders of magnitude as the oil column thickness decreases from 100 ft to 12.5 ft. This is expected to have a significant effect on completion design. The cross-sectional horizontal area of the cluster (Figure 4-46), which is a surrogate for drainage area, also decreases by about a factor of 3.

2.2.8.1.3 Discussion

One of the interesting implications of this preliminary analysis is that compartmentalization in the three phases (gas cap, oil column and water leg) can all be significantly different, even if they have identical fracture systems. This is because gas has to move through fractures in the gas column; oil only in the oil column; and water only in the water leg. Since the water leg might be quite thick, wells many connect easily over long distances. Directly above, however, the thin oil column might require a much more dense pattern. Recall that the entire co-production and TAGS process requires the following elements: inflation of the gas gap via gas injection wells; injection of steam into the gas cap to volatilize the oil-saturated matrix in the gas cap; withdrawal of oil in the oil column that has drained from the gas cap and been collected through gravity drainage in to the fracture system; and water withdrawal in the "sump" areas to assist lowering of the fluid contacts (gas injected builds both reservoir pressure and increases gas gap volume, both necessary for lowering gas-oil and oil-water contacts. This process also requires water withdrawal beneath the oil column that is being produced and re-injection elsewhere in the field).

This fracture compartmentalization has some other implications for TAGS. The drainage area or tributary drainage volume is a key issue for gas reinjection effectiveness and also for the volume of matrix that is heated by injected steam. Another important implication is that steam injector wells should be in different compartments from producers. If they are in the same compartment, then all that happens is that steam is cycled between injectors and producers, and only a minimal volume of matrix is heated. The ideal situation is for injectors to heat up matrix, which then drains oil into a separate fracture system that is produced by another well.

2.2.8.2 Preliminary Tracer Test Simulations

During 1996, Marathon carried out a series of tracer tests to improve understanding of reservoir connectivity in the Yates field. The first stage in DFN simulation for the TAGS process is to ensure that the DFN model developed can reproduce the patterns of connectivity seen in these tracer tests. Once this is achieved, these models will be used to demonstrate the heat particle tracking algorithm and the compartmentalization tools developed by the project.

The layout of the tracer testing is illustrated in Figure 2-1. The tracer experiment was carried out around a deepwater injection well located in the northern portion of the field. Both bromide and thiocyanate tracers were used. Two tracers were injected in different wells. This was done in order to evaluate multiwell communication in a cost-efficient manner, since monitoring wells could be sampled for both tracers at the same time. In addition, there is a background concentration of bromide in the reservoir, but only a negligible concentration of thiocyanate. This made it possible to evaluate the effect of background concentrations on the use of each tracer for future studies.

2.2.8.2.1 Field Measurements

Figures 2-2 through 2-7 present tracer breakthrough curves for injections in YFU 1711 and FYU 1755. Measurements indicate a minimum travel time of 2 days to YFU 1781 and 1782. Table 2-13 shows the tracer concentrations at monitoring wells during the tracer experiment. Major tracer breakthroughs were observed at four wells lying to the southeast of the injection well. Minor breakthroughs were observed in two wells lying to the east. These tests confirmed a dominant NW-SE directionality, and a much slower dispersion perpendicular to this direction.

2.2.8.2.2 DFN Model for Simulations

The initial discrete fracture network model was constructed using the geometrical and hydrogeological parameters summarized in Table 2-6. This model is referred to as the Basecase model, as it represents the initial conceptualization of the system which is then subjected to a series of sensitivity studies. The greatest uncertainty in the DFN model is in terms of the fracture transmissivity, and in particular the relative transmissivity of the deterministic faults and the stochastic fractures. This issue is addressed in the sensitivity studies below. The hydrogeological properties assigned to the fracture sets during the Basecase and sensitivity studies is given in Tables 2-6 and 2-13.

Table 2-13 Model Parameters for Sensitivity Studies

The DFN model used to for preliminary tracer test simulations was constructed within a cubical region centered on the tracer injection well 1711, with length and width of 1200 m and depth of 150 m. The top and bottom boundaries have been represented as no flow boundaries, while the lateral boundaries are represented as hydrostatic constant head boundaries. A model depth of 150 m was used as this is the saturated thickness of the reservoir, bounded to the top by the oil-water contact. The model length and width was selected to minimize the impact of the lateral boundary conditions on the simulations while ensuring the numerical efficiency of the model. The lateral boundaries have been rotated 45 degrees from North to maximize the coverage of production wells. Figure 2-47 illustrates the tracer test modeling region and identifies wells where major, minor and no tracer breakthroughs were observed.

The wells are represented in the model as internal boundaries, and are assigned appropriate injection and production rates. Nine wells have been included in the model, seven of which are active during the tracer experiment. It is assumed that these wells are open and penetrate the full thickness of the model. Production and injection rates have been calculated from averaged daily rates measured during the tracer test and are summarized in Table 2-14.

Table 2-14 Tracer Test Injection and Production Rates

The DFN Model used in these simulations is illustrated in Figure 2-48.

2.2.8.2.3 Simulations

The movement of bromide tracer in Tract 17 has been simulated in this study by using the particle tracking routine in the groundwater flow solver MAFIC. These simulations assume that Tract 17 has reached steady state flow conditions before the injection of bromide tracer at YFU 1711. Five hundred particles were released instantaneously into the steady state flow field at YFU 1711 and the movement of these particles under the steady state flow field was subsequently tracked for 25 days. It was also assumed that the bromide tracer is conservative (i.e., movement is not chemically retarded) and constant fluid density is assumed at all times.

Basecase Simulations

The movement of tracer under Basecase conditions is shown in Figure 2-49. This figure shows the location of particles 2 days after injection into the fracture network from Borehole 1711. In the Basecase model, breakthrough is simulated in Boreholes 1781 and 1782 as was observed during the tracer test in 1996, however in contrast to the actual test, breakthrough is also simulated in Borehole 1780. Figure 2-49 also shows that the simulated movement of tracer is dominated by the large deterministic faults. This is largely because these are more extensive, higher transmissivity features. The dominant trend in the stochastic fracture network is NW-SE, however due to the increased influence of the faults, this trend does not appreciably effect the simulated movement of the tracer.

In reality, the observed distribution of breakthrough of tracer shows a NW-SE trend, and this is not simulated in the Basecase model. Both the boundary conditions and the relative transmissivity of the deterministic and stochastic fracture sets may be important factors that influence the movement of tracer in the Basecase model, and these factors have been assessed in a series of sensitivity studies presented below.

A series of sensitivity studies have been conducted to assess the influence of model boundary conditions and relative transmissivity of deterministic faults and stochastic fractures. Four additional models were constructed to address these issues (Models 2a, 2b, 3a and 3b) and a summary of the parameters assigned in the sensitivity studies is given in Table 2-14.

Model 2a and b: Impact of Boundary Conditions

Model boundary conditions have been assumed in the Basecase model that reflect our geological and hydrogeological understanding of Tract 17. The influence of the top and bottom boundary conditions is likely to be significant given the proximity of the tracer injection location to these boundaries and this is assessed in Model 2a. Model 2b is a variant of Model 2a, and was constructed to assess the impact of a hydraulic gradient across the site on the simulated distribution of tracer.

Model 2a is identical to the Basecase model, however the top and bottom boundaries are assigned a constant head equivalent to the lateral boundaries. The tracer test was simulated using the Model 2a boundary conditions and the location of tracer 2 days after injection in illustrated in Figure 2-50a. This figure shows the tracer is mainly directed towards the top and bottom boundaries of the model, while the lateral movement within the model is limited. This is in strong contrast to the Basecase model illustrated in Figure 2-49, where tracer moves laterally away from the injection well. The results of this sensitivity study indicate that the Basecase model boundaries are more reasonable and better justified than Model 2a. The dramatic difference between the two models arises because the model thickness is small compared with its lateral extent. As the model is shallow the injected tracer is strongly effected by the proximity of the top and bottom boundaries and so the simulated distribution of tracer may be sensitive to the model thickness.

In Model 2b, again all the boundaries are set with a constant head condition, and a hydraulic gradient of 0.005 is applied across the model from NW to SE. Again the results are influenced by the top and bottom boundary conditions which limit the lateral movement of tracer. However, the impact of this relatively small gradient can be observed in the distribution of the tracer (Figure 2-50b), which when compared to Model 2a shows a NW-SE trend. The impact of a NW-SE gradient would be even more striking if the top and bottom boundaries were assigned no flow conditions, and this result indicates that a hydraulic gradient across the site could partly be responsible for the observed NW-SE distribution of tracer.

Model 3a and b: Impact of Deterministic Fault Transmissivity

As mentioned in previous sections there is considerable uncertainty in the transmissivity of the deterministic faults and stochastic fractures. This uncertainty may be reduced for the stochastic fracture set after the analysis of drill stem tests and spinner surveys is completed, however it is unlikely that the transmissivity of the deterministic faults can be derived from the data collected to date. In this case, the fault transmissivity becomes the subject of calibration to hydraulic events such as the tracer test injection event presented here. Calibration to such events reduces the uncertainty in our interpretation and verifies our