Engineering Geology 56 (2000) 185–195 www.elsevier.nl/locate/enggeo
Methodology development for modeling heterogeneous conductivity fields for a sandstone type uranium deposit, central Japan S. Mikake a, *, H. Yoshida a, K. Koide a, K. Yanagizawa a, N. Ogata b, K. Maekawa b a Japan Nuclear Cycle Development Institute (JNC), Tono Geoscience Center 959-31, Jorinji, Izumi-cho Toki, Gifu, Japan b Japan Nuclear Cycle Development Institute (JNC), Tokyo Office 1–9–13, Akasaka, Minato-ku, Tokyo, Japan
Abstract A methodology for modeling heterogeneous conductivity fields has been developed for the Tono sandstone type uranium deposit in Japan. The main purpose of this study is to test the hypothesis that the heterogeneity of uranium deposition is predominantly controlled by the heterogeneity of the fluid flow system. Three numerical techniques — Kriging, Multi-Kernel Modulation and Self-affine fractal simulation — have been studied. Krining and Multi-Kernel Modulation were applied in order to characterize the broad horizontal spatial structures. Self-affine fractal simulation was used to represent a high resolution vertical section through the ore body. Modeling of the flow and uranium concentration in the heterogeneous self-affine conductivity field was carried out. The results of this study show that: (1) fine resolution of the heterogeneous field can be applicable for the simulation by comparison with the calculated particle distributions and observed uranium concentration; (2) groundwater flow modelling is applicable to show that the location and magnitude of groundwater flow is constrained by heterogeneity of hydrological structure; (3) techniques such as the Self-affine fractal method can be used to simulate the influence of hydrological heterogeneity on the uranium concentration in sedimentary rock; and (4) the heterogeneity of the uranium concentration is most likely controlled by the spatial distribution of variation in the hydraulic conductivity field. © 2000 Elsevier Science B.V. All rights reserved. Keywords: Japan; Sandstone; Uranium
1. Introduction The Tono uranium deposit is distiruted in the Mizunami sedimentary basin in Tono area, central Japan (Fig. 1 and Table 1). The Tono area has been subject to numerous geological processes and * Corresponding author. E-mail address:
[email protected] (S. Mikake)
events such as sedimentation, erosion, uplift/ subsidence and marine transgression-regression within the last 20 Ma years. As a result, the sedimentary basin represents spatial heterogeneity both in terms of hydrological and mineralogical features ( Yusa et al., 1993). Groundwater flow in this system is primarily controlled by the hydrological structure which corresponds to the geometry of the geological
0013-7952/00/$ - see front matter © 2000 Elsevier Science B.V. All rights reserved. PII: S0 0 1 3 -7 9 5 2 ( 9 9 ) 0 0 14 2 - 8
186
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
Fig. 1. Geological setting of the Tono mine.
Table 1 Stratigraphic column for the Toki sedimentary basin
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
setting. In unconsolidated Tertiary sedimentary rocks hydraulic properties can be correlated with the heterogeneity of sedimentary sequences or the distribution of materials within various formations. The estimation of a heterogeneous conductivity field is of prime interest in order to understand groundwater flow in these sedimentary rocks. This study focuses on testing the validity of the assumption that spatial heterogeneity of the uranium concentration in the Tono area might be controlled by the heterogeneity of the hydrological flow system rather than by the heterogeneity of mineralogical differences. This paper describes a method for modeling the heterogeneous conductivity field by using several simulation techniques: kriging; multi-kernal modulation (MKM ) ( Williams, 1988); and self-affine fractal. The potential impact on the uranium distribution of a heterogeneous field simulated using the fractal method was investigated through particle tracking, and compared with the heterogeneous pattern of uranium deposition measured in the Tono area.
2. Geological setting 2.1. Geology The Tono area is located about 350 km west of Tokyo (Fig. 1). This area is geologically situated on the eastern part of the ‘‘Inner Zone of Southwest Japan’’. It lies within the Ryoke granite–rhyolite complex of Cretaceous age near its boundary with the Meso–Palaeozoic sedimentary rocks. A series of Neogene to Quaternary basins were formed on this pre-Tertiary basement. The Neogene to Quaternary sediments are divided into two groups, the Seto Group (Plio–Pleistocene conglomerates) and the Mizunami group (Miocene conglomerates, sandstones and tuffs). The Mizunami group is further separated into the Oidawara, Akeyo, Upper Toki Lignite Bearing and Lower Toki Lignite Bearing Formations ( Table 1). The granitic basement rock ( Toki Granite) is unconformably overlain by Mizunami Group. The Toki Granite rock and Mizunami
187
Group are displaced by a reverse fault (the Tsukiyoshi fault) with approximately 30 m of vertical displacement. The upper parts of the Toki Granite have been weathered and are the source for uranium deposits located in channel structures in the basal conglomerates of the Mizunami Group (the Lower Toki Lignite Bearing Formation). 2.2. Tono uranium deposit The Tono uranium deposit has been investigated by JNC as a natural analogue study of the uranium migration in a geological formation ( Yoshida et al., 1996). The mineralization age was inferred to be approximately 10 million years (Ochiai et al., 1989). The pattern of uranium concentration in the Tono deposit was hypothesized to be dependent upon the accessible groundwater flow paths rather than the compositional heterogeneity of rockforming minerals ( Yoshida et al., 1994). The characterization of the heterogeneous hydraulic conductivity field is therefore important for understanding of the evolution of the ore generation ( Yoshida, 1994). The uranium mineralisation is found just above the unconformity over an area several hundred meters wide, and several meters thick ( Yamakawa, 1991). This particular feature of the uranium deposit has been highlighted through a number of studies conducted since 1990 by JNC (e.g., Yamakawa, 1991; Yusa et al., 1993; Iwatsuki and Yoshida, 1999).
3. Approach for estimating the spatial distribution of hydrological properties The spatial distribution of hydrological structure was estimated by three methods: kriging; multi-kernel modulation (MKM ); and the selfaffine fractal method. Hydrological borehole data was subdivided into closed zones, and used as the basis for estimation. Kriging and MKM were used to estimate the spatial distribution of hydrological structure over a 500 m×500 m×200 m horizontal area. Nine boreholes have been drilled within this area for characterization studies. Rock mass property data were available from resistivity logging at
188
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
1 m intervals in these nine boreholes. Resistivity data has a strong correlation to hydrological data in the sedimentary rocks of Tono area (Ogata et al., 1992) and was used to develop a conductivity indicator. A self-affine fractal approach was used to simulate a 100 m×100 m vertical cross-section through the center of the Tono mine, coinciding with the uranium concentration deposit which intersected four of these boreholes. 3.1. Kriging The first approach involved estimation of the spatial distribution of hydrological structure by kriging. Classical kriging is a geostatistical approach developed in the mining industry to predict the amount and location of ore bodies. This approach is also recognized to be useful for designing the drilling interval (Asami et al., 1983). The estimation of properties at points where no observation exists is generated by the linear combination of a derived spatial correlation function applied at each of the known observation points. The correlation function used is itself an estimation of spatial covariance plot known as a variogram. 3.2. Multi-kernel modulation Following an evaluation of other relevant existing methods, we chose to use Multi-kernel modulation (MKM ). The MKM interpolates sparse and irregularly spaced data by modulating a polynomial trend surface with a linear summation of regular surfaces (kernels). This method is to modulate the least-squares trend surface, T(x, y), by adding to it an exact fit, R(x, y), of the residuals by the linear combination of functions which are asymptotic to the trend surface: Z(x, y)=T(x, y)+R(x, y)
(1)
Conceptually, this forms a synergy of two intuitively relevant factors with the mathematical basis. 3.3. Fractal approach Heterogeneous conductivity fields constructed using a fractal approach are defined by a power law variogram, i.e., fractal scaling law (Ohno,
1997). One method for simulating a self-affine fractal was used to design a heterogeneous conductivity field, which was then installed in laboratory tracer test experimental apparatus in order to validate models of flow and mass transport ( Hatanaka et al., 1996). The same approach was used to simulate uranium migration at Tono in order to build confidence in the applicability of this modeling approach. A fractal approach for the spatial variability of the borehole data was chosen for the geostatistical analysis. This approach has been applied to core logs and had been previously successfully applied to characterize the Tono permeability indicator data. The autocorrelation or spatial structure of heterogeneity of transmissivity T in geological media can be expressed by the spatial variogram c (h) defined as: c(h)=
1 2
{[y(x)−y(x+h)]2}
(2)
where x is a position, h is the distance between two measurements, h=|h|, y(x)=log T and {} denotes the spatial average over a large number of measurements. According to the field sites such as Tono mine and others, the spatial variability quantity c(h) defined in Eq. (2) may be represented by the following power law: c(h)=ah2P
(3)
for all displacements h less than a correlation length h°. For displacements greater than the correlation length where P and a are constants. This equation is usually referred to as a fractal scaling law and such a field y(x) is said to process statistical self affinity over a range of length scales. The fractal dimension is given by n+1−H in ndimensional space. 3.4. Fractal implementation A 100 m×100 m vertical cross-section through the center of the Tono area, which intersected by four of these boreholes, was chosen to be modeled (Clark et al., 1995). The location of these boreholes are shown in Fig. 2. This cross-section had a vertical extent from 100 m elevation to 200 m
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
Fig. 2. Location of the model region relative to the four boreholes.
189
190
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
elevation. This region was selected because it contains one of the most concentrated uranium deposits with a typical pattern of spatial heterogeneity. Also, this system can be represented using simple boundary conditions (i.e. no flux across the top and bottom boundaries) since it lies between two low permeability layers (Fig. 2). Hydraulic conductivity data were chosen from the available Tono data at 1m intervals along the length of each borehole. This produced a set of 299 data points for use in generating and conditioning conductivity fields. Taking the data at 1 m intervals meant that the data should not have any correlations induced by the measurement technique. These data were used within the macroaffinity code (Impey et al., 1994) to conditioning simulation of anisotropic fractal conductivity fields on a 50×100 cell rectangular grid covering the cross-section so that the flow field can be described with a 1 m (vertical ) and 2 m (horizontal ) resolution. A range of vertical and horizontal fractal dimensions were used to produce a total of 11 conductivity realizations, each of which fitted the measured
borehole data. The dimensions used were based upon previous analysis of the fractal dimensions and variance of the permeability data from individual boreholes in the Tono area (Humm and Impey, 1994). The vertical dimension (D ) was y varied from 1.6 to 1.9; the horizontal dimension (D ) was varied from 1.4 to 1.6. The generated x conductivity fields displayed differing degrees of anisotropy, but each had a general horizontal banded appearance, with these bands made up of stretched patches or lozenges of high and low permeabilities. The detailed heterogeneity characteristic of fractal fields was apparent. A head difference of 10 m which corresponds to a gradient similar to the present condition was applied across the region for the flow calculations. No-flow boundary conditions were imposed on the top and bottom boundaries as discussed above. The flow solver, part of macro-affinity, then calculated the pressure head field and Darcy flow field for each conductivity field. Three of the calculated flow fields were then used for the transport calculations. Transport calculations were made using a step
Fig. 3. Distribution of hydraulic conductivity using the kriging method.
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
input of unit mass applied at one of seven source injection points. The duration of the input step release was equal to the size of the time step used for the particle tracking, that is 100R years, where d R is a retardation factor less than several thoud sand. The simulations used had dimensions of (D , D )=(1.6, 1.6), (1.4, 1.9) and (1.6, 1.9), x y respectively. The first two of these fields illustrated the extremes of isotropy and anisotropy. The third field had the variance closest to the conditioning data variance. The validity of this fractal method was investigated using a particle tracking method.
4. Results and discussion 4.1. Estimation of heterogeneity of hydrological structure The horizontal spatial distribution of hydrological properties was estimated using kriging and MKM methods for a region that was 500 m× 500 m×200 m in the Tono area. The MKM method provides great continuity of hydrological
191
structure since estimate are related to observations that lie at relatively large separation distances ( Figs. 3 and 4). As a result, the MKM method is believed to better represent these sedimentary systems than the kriging method. First, there exists at least some statistical information in any sampling and that information can be used to bring the interpolation under control in the remotest areas. The trend surface can be viewed as an estimation of the global mean of the surface that is represented by the input data sample. Second, the influence of a control point on the interpolation for regions away from that control point should be a decreasing influence with distance. This requirement, leading to the introduction of asymptotic functions, guarantees that the estimated mean is the major determinant for remote area interpolation or extrapolation. The kriging and MKM methods are based on the statistical assumption that there are functions determined by spatial arrange of the random distribution of samples. These methods cannot be used to interpret a heterogeneous conductivity field
Fig. 4. Distribution of hydraulic conductivity using the MKM method.
192
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
Fig. 5. Distribution of hydraulic conductivity.
because the regular function is determined by the relationship between the variogram and the distance of the sample points. Fractal theory (statistical self-similarity) is, however, recognized to be
useful for estimating the heterogeneity of hydrological structures. The calculated hydraulic conductivity and Darcy flux field, where D =D =1.6, are shown y x
Fig. 6. Darcy flux magnitude field.
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
as an example in Figs. 5 and 6. Fig. 5 shows the hydraulic conductivity distributions. The Darcy flux magnitude field is shown in Fig. 6. The dark regions are low flux. The light regions are moderate flux. The highest flux is a dark stripe above lowest dot. The dots indicate particle sources. The high flow region calculated in the lower parts of the cross-section (a dark stripe above the lowest dot in Fig. 6) corresponds to a weathered granite zone through which uranium has been transported from its source regions. 4.2. Transport through the vertical section The goals of the transport calculations were to determine the applicability of the self-affine fractal method for the simulation of hydrological properties for sedimentary rocks and to test the extent to which the pattern of heterogeneity in the uranium deposition at Tono can be explained by controls on hydrological access or the location of flow paths. The upstream ends of seven individual channels were selected as source locations ( Fig. 6). It is clear that the largest contribution is from source No.2 (Fig. 6) which corresponds to the weathered granite and the adjacent conglomerate layer which
193
is also highly conductive. This is consistent with the geological inference described earlier. Another interesting feature is that the peak of source No. 2 comes one or two time steps earlier than that of the other sources. Again this looks reasonable if one takes into account that the channel, which runs through source No. 2 corresponds to the main flow path at the bottom of the sedimentary basin and the others are branches to the downstream part of the uranium deposits. Fig. 7 shows the corresponding particle distributions at 1000R years. The red tones have the highest d concentration and the white tones are the lowest concentration. Ellipses mark the observed uranium deposits along boreholes. The white areas are low uranium concentrations. The dark areas are intermediate values, which outline the light toned high concentration regions. The ellipses indicate the observed uranium deposits in the four boreholes. The high particle density zones show the location of the actual uranium deposits within the area. In addition, comparison with the flow field in Fig. 6 indicates that the high particle density zones coincide the high flow rate channels, in particular the channel which corresponds to the layer of weathered granite and adjacent conglomerate.
Fig. 7. Particle concentration distribution at 1000R years. d
194
S. Mikake et al. / Engineering Geology 56 (2000) 185–195
What is striking is that uranium is deposited at high concentrations in the region near the fastest channel but not actually in the channel (see the bottom of the two left-hand boreholes) both for the observed and the simulated distributions. One possible interpretation of this observation is that the peak flux through this channel occurred a couple of time steps earlier than the others and the plume of migrating uranium has already passed.
5. Conclusions Simulations of the flow and transport in the heterogeneous conductivity field were attempted to test the validity of the assumption that the spatially heterogeneous pattern of the uranium concentration is controlled by the hydrological access corresponding to the highly channeled flow structure. The fine resolution of the simulation, which is similar to the spatial scale of the uranium deposit enabled us to compare the calculated particle distribution with the heterogeneous and irregular pattern of the observed uranium deposition. Comparison of the simulated particle density estimates and observed uranium concentration is in qualitative agreement indicating that the selfaffine fractal simulation approach can be used to simulate the spatial heterogeneity of sedimentary systems in the Tono site. In addition, particle tracking simulations showed that groundwater flow occurred along paths of preferentially high hydraulic conductivity. Results indicate that in systems where the contrast in conductivity is sufficiently great, the location and magnitude of groundwater flow is constrained by heterogeneity of hydrological structure. Finally, correspondence of the spatial pattern of flow paths with the current uranium deposit concentrations supports the theory that uranium deposition was highly controlled by the arrangement of flow properties.
Acknowledgements The authors would like to express their thanks to JNC (Japan) for permission to publish the
results of this study. Thanks also to Dr E. Webb (Sandia, USA), Dr K. Tsubota (JNC ) and Dr. R. Metcalf for their constructive reviews and discussions.
References Asami, N., Kataya, K., Akiyama, S., Hirosaki, F., 1983. Applying geostatistics to design the drilling interval and pattern for Kutcho deposit in Canada.. Mining Geology 33 (2), 131–136. Clark, K.J., Takase, H., Impey, M.D., Humm, J.P., Maekawa, K., Ogata, N., Yanagizawa, K., 1995. Natural analogue study of uranium migration in the Tono mine. Sci. Basis Nucl. Waste Manag. XIX, 767–773. Hatanaka, K., Watari, S., Uchida, M., Takase, H., Impey, M.D., 1996. Experimental study on groundwater flow and mass transport in a heterogeneous porous medium. Sci. Basis Nucl. Waste Manag. XIX, 739–746. Humm, J.P., Impey, M.D., 1994. Fractal fitting of the Tono mine data.Intera report ID3249-7 Version 1. Impey, M.D., Williams, M.J., Humm, J.P., 1994. The MACROAFFINITY code (MACRO Version 1.2), theoretical background. Intera Report ID3618-1 Version 3. Iwatsuki, T., Yoshida, H., 1999. Groundwater chemistry and fracture mineralogy in the basement granitic rock in the Tono uranium mine area Gifu prefecture, Japan — Groundwater composition, Eh evolution analysis by fracture filling minerals. Geochemical Journal 32 in press. Ochiai, Y., Yamakawa, M., Takeda, S., Harashima, F., 1989. A natural analogue study on Tono uranium deposit in Japan. CEC Rep. Proc. 3rd NAWG., 126–138. Ogata, N., Ohsawa, H., Nakano, K., Yanagizawa, K., Nishigaki, M., 1992. Correlation between the electric resistivity and the hydraulic conductivity of sedimentary rocks and its application to the modelling of hydrogeological structures. Journal of the Japan Society of Engineering Geology 32 (6), 51–62, (in Japanese). Ohno, H., 1997. Application of fractal concept in geotechnology. Journal of the Japan Society of Engineering Geology 38 (3), 159–173. Williams, R.L., 1988. A technique for the geometric modeling of underground surfaces. SANDIA REPORT, SAND84-0307. Yamakawa, M., 1991. Geochemical behaviour of natural uranium series nuclides in geological formation. Proc. 3rd International Symposium on Advanced Nuclear Energy Research Global Environment and Nuclear Energy. pp. 150–158. Yoshida, H., 1994. Relation between U-series migration and microstructural properties of sedimentary rocks. Applied Geochemistry 9, 479–490. Yoshida, H., Yui, M., Shibutani, T., 1994. Flow-path structure in relation to nuclide migration in sedimentary rocks an
S. Mikake et al. / Engineering Geology 56 (2000) 185–195 approach with field investigations and experiments for uranium migration at Tono uranium deposit central Japan. Journal of Nuclear Science and Technology 31 (8), 803–812. Yoshida, H., Sakuma, H., Yusa, Y., 1996. The Tono natural analogue study Program. In: von Maravic, H., Smellie, J. ( Eds.), Sixth EC Natural Analogue Working Group Meeting. Proceedings of an International Workshop held in Santa
195
Fe New Mexico on 12–16 September 1994. European Commission Report EUR 16761 EN. European Commission. Yusa, Y., Ishimaru, K., Ota, K., Umeda, K., 1993. Geological and geochemical indicators of paleohydrogeology in Tono uranium deposits Japan. IN: Paleohydrogeological Methods and their Applications. Proc. NEA Workshop, Paris, 9–10 Nov. 1992. OECD, Paris, pp. 117–146.