A three-dimensional model of radionuclide migration from a canister of spent nuclear fuel

The development of a three-dimensional model of the near-field of a geological repository for spent nuclear fuel is described. The near-field comprises a canister where the spent fuel is enclosed, a layer of bentonite clay around the canister, a backfilled tunnel over the canister’s deposition hole and finally, the fractured rock adjacent to the bentonite and the tunnel. The main transport processes are diffusion, sorption, radioactive decay and groundwater flow. A mathematical model attempts to couple the water flow to the mass transport. This model has been simplified to reduce its computational complexity. Our results are compared with the results of a compartment model obtained from the literature. It is concluded from the agreement between our 3D model and the compartment model, that the resistance approach used in one-dimensional compartment models is robust enough for use in models for probabilistic risk analysis of long-term performance of a geological repository.


INTRODUCTION
Integrated assessments of repositories for high-level nuclear waste or spent nuclear fuel rely heavily on modelling of the long-term migration of nuclides throughout the engineered barriers of the repository (near-field), the geosphere (far-field) and the biosphere.
Assessing the long-term performance of the repository cannot be done without taking into account different types of uncertainties.Monte-Carlo (probabilistic) methods are very useful for tackling the analysis of uncertainty and they are often used nowadays in risk analyses of nuclear waste.However Monte-Carlo methods require the use of fast codes because they need to be run thousands of times to obtain significant statistics.To achieve enough speed, transport models are generally coded as one-dimensional models, i.e. compartment models.However, the process of dimensional reduction, can introduce different types of errors [1], so this aspect of model uncertainty should be addressed whenever the use of simplified models is necessary.The goal of this work is to access model uncertainty related to the dimensionality of two models of near-field migration (our 3D model versus the 1D model of Hedin [2]).Hedin developed his model to speed up Monte Carlo calculations for risk analyses 1 .The two models are conceptually close to one another.Apart from the dimensionality, the main difference between the models is the use of the resistance approach to mass transfer that Hedin introduced in his model.The resistance approach is based on earlier work by Neretnieks [3].In our work we do not use the resistance approach.
In this work we present the development of a 3D finite-element model for the near-field of a repository for spent nuclear fuel, and we compare our results with those obtained using the analytic solution of Hedin.

THE THREE-DIMENSIONAL NEAR-FIELD MODEL MODEL DESCRIPTION
The system to be modelled comprises a (perforated) canister containing the nuclear fuel, a protective barrier of bentonite clay in which the canister is embedded, and a plane parallel fracture in the rock adjacent to the bentonite (see Fig. 1).It is assumed that the canister pinhole has a constant cross-section over the first twenty thousands years.
Figure 1 The near-field model including the perforated canister, the bentonite clay, the rock and a fracture crossing the deposition hole.
The conceptual model for the fuel released from the canister follows the work of Hedin [2].Assume that at time t=0, a given species with decay constant λ has, an inventory M o .The canister corrodes with time and, after a certain delay time, t Delay , water enters in the canister through the corrosion hole, coming into contact with the fuel.At that time, the inventory of the waste is given by M(t Delay ) = M 0 e -λtDelay It is also assumed that a fraction α of the waste is instantaneously dissolved in the water and therefore, can readily migrate through the pinhole.The rest of the material embedded in the fuel matrix is released gradually by the matrix at a rate, P(t): where: M 0 -is the inventory at time t o (s) D F -is the degradation rate for the fuel matrix and is a constant (s -1 ) α -is the fraction of the waste that is immediately available for migration (-).
The main mass transport processes are diffusion, retardation and radionuclide decay in the canister, pinhole and bentonite clay and groundwater flow in the fractures.The release rate from the near-field is computed at the interface between the bentonite clay and the rock.
Hedin uses in his model a steady state flow corresponding to the Darcy velocity value q given in Table 1 (2.0 x 10 -3 m y -1 ).The steady state assumption in Eqs. 2 and 3 conforms to that of Hedin.
The advective velocity u of the groundwater in the fractures in the vicinity of the deposition hole can be computed by solving the three-dimensional steady state Navier-Stokes' equation assuming fluid incompressibility (Eq.2), together with the continuity equation (Eq.3): The Navier-Stokes' equation, above, is valid for a Newtonian, steady state fluid.The steady state flow is coupled to the mass transport that is described by the time-dependent advectiondiffusion equation (Eq.4), to which terms for sorption (retardation) and radionuclide decay have been added: where c i -is the nuclide concentration in the domain Ω i (mol m -3 ) R i -is the retardation coefficient (-) D i -is the diffusion coefficient (m 2 s -1 ) u i -is the advective velocity in the fracture (m s -1 ) λ i -is the radioactive decay constant (s -1 ) The domains Ω i correspond to the different regions of the model (canister gap, pinhole, etc.).
The model given by Equations 2-4 is a typical multiphysics model where the groundwater flow and the mass transport are coupled.We have, however, simplified this model, given that we have the necessary data from the water flow field.Because of this, we can disregard the rock adjacent to the canister, thereby avoiding the need to solve the Navier-Stokes' equation (Eq.2).Bentonite clay is almost impermeable to flowing water, its advective velocity being extremely low.As a result of this, the migration of the nuclides dissolved in the water of the canister gap, through the pinhole and the bentonite is controlled by diffusion.Disregarding the velocity u in Eq. ( 5), we obtain: (5) As stated above, the rock fractures are not modeled explicitly.The impact of the fractured rock is taken into account using a boundary condition, as explained in the fourth section.This approach reduces considerably the computational complexity of our 3D model.

MODEL GEOMETRY AND INPUT DATA
Figure 2 shows the geometry used to model the system.It is important to save memory and CPU-time and this can be done by exploring the symmetry of the system.The canister is 4.833 m high and has an inner diameter of 0.95 m.The thickness of the canister is 0.05 m.The cross section of the pinhole is 1.0x10 -6 m 2 .The bentonite clay is in the form of a cylinder with a outer radius of 0.875 m and is 6.833 m high.The volume of the canister gap where the nuclides are dissolved in, is 1.0 m 3 .The fracture aperture is 1.0x10 -4 m.The challenge in using the finite-element method to model our system lies in the huge variation of the occurring physical dimensions (four to five, orders of magnitude).
It is assumed that a pinhole will be fully developed after 300 years after the canister deposition.Water can at that time, penetrate into the canister and dissolve the nuclides allowing them to escape by means of diffusion through the pinhole.Those nuclides will continue to diffuse through the bentonite until they reach the water flowing in the rock fractures surrounding the borehole.It is assumed that the groundwater flow in the fracture corresponds to a mean Darcy flow of two millimetres per year (see Table 1).
The generic data used in the model is given in Table 1: The nuclide used in the model to calculate the breakthrough curves is nickel-63.Nuclide specific data includes the half-life, T 1/2 , the groundwater diffusivity, D water , the effective diffusivity of the bentonite, D Bent eff , the distribution coefficient in the bentonite, K Bent d , the initial inventory, M 0 , the fraction inventory that is immediately accessible to migration, α, the (constant) degradation rate of the fuel matrix, D F and the nuclide solubility, C sol .Those values are given in Table 2:  The different domains (canister gap, pinhole, bentonite clay blocks around the canister, bentonite clay blocks over the canister and bentonite clay blocks under the canister) are meshed separately.The total numbers of degrees of freedom are 132873.Figure 3 shows partial views of the final mesh used for the canister gap and the bentonite around the canister.
We have also used a less refined mesh in the vertical direction to check that our final results are mesh independent.

BOUNDARY CONDITIONS
To mimic the impact of the fractured rock on the bentonite clay we balance the diffuse flux out of the bentonite buffer with the advective flux in the rock.If the average concentration outside the buffer is C m , the the advective flux can be approximated by where L is the vertical height of bentonite buffer, R is the outer radius of bentonite buffer and q is the Darcy velocity in the host rock.It is assumed in eq.6 that the relevant flow through the rock is over an area equal to the cross-sectional area of the bentonite buffer.Taking again an average, the diffusive flux from the bentonite buffer is Besides the boundary condition given by eq. 9 above, we need one equation expressing insulation/symmetry for the relevant surfaces.This condition is given by -n .(-Dٌc) = 0 (10)

RESULTS AND DISCUSSION
The estimation of the release of the Ni-63 radionuclide is performed for a fracture density corresponding to one fracture per canister.The release rate obtained is shown in Fig. 4. Our model is three-dimensional and the model of Hedin is one-dimensional.A discrepancy between the two curves is therefore expected.It reflects a conceptual-or model uncertainty.Different types of uncertainties are common (and usually large) in the field of radionuclide transport.Therefore our result agrees well with the breakthrough curve of Hedin [2].The distribution of the nuclides concentration in the bentonite is shown in Fig. 5 for three different times: 369 years, 439 years and 856 years.The fact that the pinhole has a so small cross section results in a plume dispersion that is spheric in form, following the assumption of Neretnieks [3].It is therefore to expect that the resistance formulas used in the 1D analytical models should be in good agreement with the 3D numerical results.Fig. 6 shows a horizontal cut of the radionuclide concentration profile at the level of the pinhole.Here we can also observe the spherical form of the plume.Fig. 7 shows the corresponding vertical slice through the bentonite.The diffusive flux streamlines that can be observed from a top-down perspective in the bentonite surrounding the casing, are shown in Fig. 8.

SUMMARY AND CONCLUSIONS
In recent years, compartment models have been the workhorses of radionuclide transport studies in certain programs of assessments of the long-term performance of geological repositories for spent nuclear fuel.These models are a way to resolve the need for fast computer codes to be used in probabilistic risk analysis.However, the indiscriminate use of compartment models to predict the evolution of complex environmental systems raises a number of questions, from their numerical accuracy to transparency and QA, subjects that are outside the scope of this paper.Another aspect of modelling complex systems is the issue of model (or conceptual) uncertainty.This aspect is indirectly addressed in this paper in the sense that the use of alternative models to describe the behaviour of a system can give different results.However, this does not imply that any of the models are wrong, but rather, that there is sufficient leeway for different assumptions that are reasonable and that lead to different interpretations of the input data and consequently to different results.This situation is not unusual in the field of nuclear waste assessments and is tackled by probabilistic methods.
In this paper we focus on the methodology behind the set-up of (more or less) black boxes representing the different compartments of a 1D compartment model and using equivalent resistances to couple the mass transfer between the boxes.Using a case study of radionuclide transport in the near-field of a repository, we compare the results obtained by our 3D finiteelement model with the ones obtained by a 1D compartment model.We have been very careful to use assumptions that are the same or very close to those used in the 1D model.It is shown that, for a pinhole of small cross-sectional dimensions representing a perforation in the canister containing the spent fuel, the results of the two models are in good agreement.Our results are valid up to 20 000 years only because the present version of our model has some limitations.To improve the three-dimensional model, we need to a) allow the pinhole cross-section to change as a function of time because it is this feature that does not allow us to compare our results with those of the 1D compartment model after 20 000 years, b) to introduce the Navier-Stokes' equation for groundwater flow and to couple it to mass transport (making the model a true multiphysics model) to check the accuracy of the use of a mean water flow field to describe the impact of the fractured rock that we actually mimic using a boundary condition, c) to introduce a tunnel over the deposition hole of the canister, d) to introduce the disturbed zones at the floor of the tunnel and at the interface between the bentonite and the rock and e) to introduce matrix diffusion in the nearby rock.These features will be included in the next version of the model.
It has also been shown that it is feasible nowadays to use general commercial software to model large systems using finite elements in 3D.Quality assurance of numerical algorithms and the model transparency that can be obtained by using the real geometry of the systems to be simulated will be strongly enhanced.The availability of such FEM packages capable of true multiphysics modelling represents an important advance in the numerical modelling of complex systems.

Figure 2
Figure 2 The geometric model showing the canister, the pinhole and the bentonite.Domain symmetry is used to simplify the model.

Figure 3 -
Figure 3-Partial views of the final mesh used for the canister gap (a) and the bentonite around the canister (c).
the bentonite buffer diffusivity.Balancing these fluxes locally at all points on the surface gives

Figure 4
Figure 4 The breakthrough curve of Ni-63.

Figure 5 -
Figure 5-The evolution in time of the concentration in the bentonite of the Ni-63 nuclide escaping from a pinhole with a cross-sectional area of 1.0E-6 m 2 .

Figure 7
Figure 7 The concentration profile of Ni-63 on the vertical plane through the pinhole.

Figure 8
Figure 8 Top-down view of streamlines showing the diffusive flux leaving the pinhole.

Table 1 :
Generic data