1. Introduction

Radioactive waste is very hazardous to humans and the biosphere. Therefore, disposal facilities for radioactive waste are designed to provide long-term isolation of the waste from the human environment by using the synergetic system of natural and engineered barriers [1-6]. Generally, to isolate the radioactive waste in the repository site, two fundamental methods must be considered. The first method is hydrologic isolation which the waste has to be kept out of contact with the water. The second consideration is what prevents the radionuclide migration out of the disposal facility by using the physical barriers with very low diffusion coefficients or high sorption coefficients such as concrete and clay. Based on these methods, the principle of radioactive waste disposal is to maintain the radiation doses and risk to as low as reasonably achievable from the disposal action. So, not only the individual doses but also the corresponding risks should be assessed in any optimization procedure. Hence, it is essential to quantitatively evaluate the performance of the waste disposal facilities by using safety assessment models. To give better confidence in the safety assessment process, various experiences have been studied in the development of computer codes for safety assessment models, the improvement of the database of radionuclide transfer parameters, and the development of scenario generation methodology [7-9].

Depending on the method or depth of considering various complex phenomena affecting nuclide transfer, safety assessment models can be classified into systemlevel and process-level. Traditionally, the system-level model that has been used in the safety assessment of the radioactive waste disposal facility simplifies complex phenomena by using several constants and focuses on the final dose calculation. Thus, this approach has many limitations in considering the evolutionary characteristics of the disposal system over time [10]. However, the process-level model carries out realistic and reliable performance evaluation by reflecting the various mechanisms and the change in site conditions. In the various process-level model, we chose COMSOL Multiphysics because it provides an efficient tool to handle the different modeling tasks, involving flow and transport. And it numerically enables implement the complex nuclide transfer mechanism owing to its flexible interface for preand post-processing [11].

Depending on the delivery regulation of radioactive wastes in Korea, fluidizable radioactive wastes such as concentrated wastes, spent resins, sludge, and dry active wastes have to be homogeneously solidified or immobilized, and then packaged into the waste container [12]. These treatments assist to limit the leaching of radionuclides and maintaining the structural integrity of the disposal site. Especially, the release rate of radionuclides from waste containers has to be low enough to secure the safety of a repository. So, the leachability index defined in the waste acceptance criteria indicates the degree of release of radionuclides and should be greater than 6 for the defined radionuclides such as Co, Sr, and Cs [13]. This index for each radionuclide can be obtained by performing the leaching test suggested in the ANS 16.1 [14]. In the case of the solidified waste, the leakage property of it is one of the key issues for the safety assessment model. In this regard, we introduce the mass transfer coefficient of the flux node in COMSOL Multiphysics to reflect the leakage property of the solidification.

In this study, we focus on the development of the safety assessment model for the underground silo of Wolseong Low and Immediate Level Waste (LILW) disposal facility in Korea. The simulations were done using COMSOL Multiphysics based on a finite element computer code. To obtain the mass transfer coefficient for the solidified waste form, we carry out leaching tests for solidified specimens using the method suggested in the ANS 16.1. Stable nuclides considered in leaching tests to simulate radionuclides are Co, Sr, and Cs, which are important radionuclides in the solidified operational wastes.

2. Methods

2.1 Model Definition

The hydrogeologic setting for the migration of radionuclides at the steady state shall be described by the conceptual model based on the underground silo #1 layout of the Wolseong LILW disposal facility in Korea, shown in Fig. 1. The model geometry consists of the waste container (WCon), the backfill (BFill), the concrete line (CLine), the excavation damaged zone (EDZ), and the host rock (HRock). The BFill is crushed stone which is a form of construction aggregate. Next, the HRock is predominantly composed of granite. The EDZ is the HRock zone around an underground excavation where irreversible changes have taken place. Such EDZ is inevitable because of the changed conditions resulting from the removal of rock during the construction process. Fresh groundwater enters the silo from the left boundary and moves right direction. The vertical flow of groundwater is not considered. The number of elements is 60,695 nodes.

Fig. 1

The conceptual model of the waste disposal facility.

2.2 Numerical Illustrations

A two-dimensional numerical model was developed by using COMSOL Multiphysics (ver.5.6) based on the proposed conceptual model in Fig. 1. The steady flow of groundwater was carried out in the Subsurface Flow Module of COMSOL Multiphysics. In the Subsurface Flow Module, the Darcy’s law is applicable for a saturated flow in porous media and expressed as equation (1):

where K is hydraulic conductivity (m·s⁻¹), u is the Darcy’s velocity or specific discharge vector (m·s⁻¹), and H is the hydraulic head (m).

The transport of radionuclides dissolved in groundwater was considered through various mechanisms such as advection, diffusion, dispersion, and adsorption [15]. The advection defines that the dissolved nuclides in groundwater are carried along with the flow direction of the groundwater. The diffusion determines that nuclides in groundwater will move from the area of greater concentration (WCon) toward the surrounding areas such as BFill, CLine, EDZ, and HRock which are less concentrated. Next, the groundwater containing the leaked and dissolved nuclide is not all traveling at the same velocity in porous media because the mixing event occurs along the various flow path. This mixing phenomenon is called dispersion. The dispersion is divided into the longitudinal and transverse dispersion. The former dispersion mixing occurs along the direction of the flow path. On the other hand, the transverse dispersion occurs when an advancing nuclide front spread in the direction normal to the direction of flow. Finally, as the nuclide species travel through a porous medium, they typically attach to (adsorb), and detach (desorb) from the surface of the solid phase by physicochemical processes, called adsorption. However, most of the adsorption processes and parameters are clearly unknown to the specific conditions of the waste form. COMSOL Multiphysics basically provides various adsorption models such as Langmuir, Freundlich, Toth, and BET. These default adsorption models typically require the related parameters. However, the value of the parameters for each model have not been provided under the condition of the disposal facilities. Therefore, the basic adsorption models pre-configured in COMSOL are difficult to apply as adsorption model for the disposal facilities. In addition to these default models, COMSOL provides user-defined modes. The adsorption model of the current safety evaluation model using the GoldSim code is implemented as a simple adsorption equation with an adsorption distribution coefficient. This model may be underestimated because it simply represents the removal effect due to adsorption. Therefore, no adsorption mechanism is considered in this study. Additionally, the radioactive decay of each nuclide was superimposed during these transport processes at the same time. Based on these migration mechanisms and decay of the radionuclide, the time-dependent species transport equation in porous media is defined at the Transport of Diluted Species Module in Porous Media and described by equation (2):

where ε_p denotes the porosity, c_i is the concentration of the nuclide (mol·m⁻³), D_i denotes the effective diffusion of the nuclide (m²·s⁻¹), u is the Darcy’s velocity or specific discharge vector (m·s⁻¹), λ_i is the half-life of the nuclides (s⁻¹).

The radionuclides’ mobility properties of the waste container are difficult to clearly define because they are generated from various generation streams [16]. Thus, we have simplified the nuclides released from the waste container to the surrounding domains through the boundary flux concept [17-19]. In the Transport of Diluted Species Interface of COMSOL Multiphysics, the flux node can be used to specify the species’ molar flux across a boundary. For example, the flux can occur due to chemical reactions or a phase change at the boundary and represent the transport to or from a surrounding environment. In this case, the prescribed mass flux rate corresponds to equation (3):

where J_o is the flux rate of external convection (mol·m⁻²·s⁻¹), k_c is a mass transfer coefficient (m·s⁻¹), c_b is the bulk concentration (mol·m⁻³), and c_i is the typical concentration far into the surrounding exterior domains.

The overall computational analysis process was summarized as follows: Firstly, the Darcy’s velocity is analyzed in the Subsurface Flow Module by introducing density, porosity, and hydraulic conductivity as the parameters. Next, the analyzed Darcy’s velocity field was set as the flow of the underground water for the disposal facility, and then the Transport of Diluted Species Interface was used to apply advection, dispersion, diffusion, and radioactive decay to calculate the concentration of radionuclides. In here, advection is calculated based on the Darcy’s velocity field. Diffusion is determined using the diffusion coefficient, while dispersion is calculated by longitudinal dispersivity and transverse dispersivity. In addition, the release characteristics of radionuclides from the waste container (WCon) were implemented using the mass flux rate with a mass transfer coefficient. The mass transfer coefficient for each radionuclide was determined from leaching tests.

2.3 Simulation Scenarios

As analysis scenarios, the migration of nuclides was evaluated for three scenarios derived from the combination of the intact engineering barrier (CLine) and the solidified waste form. Firstly, the small porosity and low hydraulic conductivity of CLine prevent contact between the WCon and groundwater by minimizing the percolation flux. Next, in the case of the solidified waste form, the nuclide transfer was delayed by the interaction between the nuclide species and the solidified medium. This delay mechanism was implemented with the mass transfer coefficient of the flux node at numerical modeling. Generally, the value of the mass transfer coefficient was determined depending on the type of nuclide and the solidified material. Among the major regulated nuclides in Korea, Co, Sr, and Cs, which are regulated in the acceptance criteria for the leachability index of the solidification, were preferentially carried out. Table 1 summarizes the three scenarios. Case A introduces the intactness of CLine and immediately release of radionuclide from WCon. Case B applies the degradation of CLine and immediately release of radionuclide in WCon. Case C includes the degradation of CLine and the solidified waste form at WCon. The solidified waste form was implemented by using the flux rate of radionuclide which simulated upon ⁵⁸Co, ⁶⁰Co, ⁹⁰Sr, and ¹³⁷Cs.

2.4 Determination of the Mass Transfer Coefficient

Considering the governing equations in COMSOL Multiphysics, the mass transfer coefficient is a diffusion rate constant that relates the mass transfer rate, mass transfer area, and concentration change as the driving force. To determinize the mass transfer coefficient of the solidified waste form, the cement specimen with stable nuclides was manufactured and the leaching test of that was performed based on the ANS 16.1 standard [14].

The solidified cement specimens were manufactured with a cylindrical shape (diameter = 50 mm, height = 100 mm) and water (W) to dried ingredient (D) ratio of 0.4 (W/D). The dried ingredient involves the commercial cement and stable nuclides source (Co, Sr, and Cs). The commercial cement used the Ordinary Portland Cement (OPC) produced by Ssangyong Cement Industrial Co. in Korea. Cobalt Nitrate (Co(NO₃)₂), Strontium Nitrate (Sr(NO₃)₂), and Cesium Nitrate (CsNO₃) were used as the stable isotopes of nuclides source. These materials were analytical reagent grade, respectively. Each specimen was prepared to contain the 10 mmol of one stable nuclide. Briefly, the certain amounts of source materials were dissolved in deionization water. The cement paste containing the stable nuclides was prepared by mixing OPC and each nuclide aqueous solution. The as-prepared cement paste was filled in the hand-made acryl mold with a cylindrical shape. The prepared cement specimens were cured at room temperature under the 90 ± 5% of relative humidity for 28 days.

The leaching tests of the obtained specimens were performed in compliance with the ANS 16.1 standard. In detail, the tests were implemented by immersing the specimens for 90 days by using the deionized water as the leachant. Each cement matrix and the leachant were placed in a polypropylene vessel. Because all surfaces of a specimen were contacted with the leachant, the specimens were perfectly submerged in the leachant. The leachant be entirely replaced and sampled at designated time intervals (2 h, 7 h, 1 d, 2 d, 3 d, 4 d, 5 d, 19 d, 47 d, and 90 d). The stable nuclides (Co, Sr, and Cs) concentration of leachate leached from the cement specimens were quantified by using an Inductively Coupled Plasma Mass Spectrometry (ICP-MS, Elan DRC II, Perkin Elmer, USA).

From the experimental data of the leaching test, the relationship between the mass transfer coefficient (k_c) and the leaching time is not linear, as shown in Fig. 2. The function that relates k_c and the time can be well fitted by the Belehradek equation (4),

Fig. 2

Fitting results for the mass transfer coefficient, using the experimental data from the reaching test of (a) Co, (b) Sr, and (c) Cs.

The fitted curve with red color for the mass transfer coefficient is plotted in Fig. 2. The yielded values of parameters (a, b, and c) were summarized in Table 2. We apply the above equation (4) as the governing equation of the mass transfer coefficient of Co, Sr, and Cr.

2.5 Model Assumptions and Boundary Conditions

For the migration model of nuclides, it is necessary to consider the dissolving and moving phenomenon in groundwater. So, the flow streamline of groundwater has to be considered. The following simplifying assumptions regarding fluid flow:

The migration of nuclides is based on the following assumptions:

The simulations are concerned with the regulated 14 nuclides such as ³H, ¹⁴C, ⁵⁵Fe, ⁵⁸Co, ⁶⁰Co, ⁵⁹Ni, ⁶³Ni, ⁹⁰Sr, ⁹⁴Nb, ⁹⁹Tc, ¹²⁹I, ¹³⁷Cs, ¹⁴⁴Ce, and ²³⁹Pu. These isotopes are contained in the radioactive waste container. However, the exact amount of each nuclide stock in the waste container is still unknown. To determine only the transfer tendency of each nuclide, the initial radioactivity of all isotopes collectively introduces 3.7×10¹⁰ Bq·m⁻³. Table 3 describes the hydrological parameters used in the saturated homogenous porous medium for the numerical model of flow and transport at the site of the Wolseong low and immediate level waste (LILW) disposal facility in Korea. These parameters were modified based on Low and Intermediate-Level Radioactive Disposal Facility Safety Analysis Report issued by Korea Radioactive Waste Agency.

3. Results and Discussions

3.1 Case A

Fig. 3(a) shows the flow streamlines of groundwater obtained through numerical simulation, which passes through the silo. The degree of Darcy’s velocity is displayed by the color legend with the red color indicating an intense flow velocity. The flow velocity of the EDZ domain represents higher than the other domains because it has a higher hydraulic conductivity than its surroundings such as HRock and CLine. Otherwise, Figs. 3(b) and 3(c) present that the inside of the silo was calculated with the lowest Darcy’s velocity due to the low permeability feature of the CLine. The value of Darcy’s velocity (streamline plot) and the specific radioactivity (Bq·m⁻³) of radionuclide (surface plot) simulated under a time-dependent study are represented in Fig. 4 over the 300-year compliance period. The color legend indicated the degree of the specific radioactivity, which set the red color as the highest specific radioactivity. The very low flow velocity in the BFill suggests that the migration of nuclides is expected to be mostly driven by diffusion. So, the nuclides spread throughout the BFill from the WCon to the CLine. Also, this trend of diffusion is no significant difference depending upon the nuclide because the diffusion coefficient of each nuclide was identically inputted in the same domain. Based on the result of Case A, the nuclides diffused from the WCon were kept inside the silo without the leakage of those while the integrity of the CLine is maintained.

Fig. 3

(a) The streamlines of groundwater flow, (b) the flow velocity of red line in the X direction, and (c) the flow velocity of green line in the Y direction for the waste disposal facility at Case A.

Fig. 4

The nuclide release of the waste repository facility at 300 year in the Case A.

3.2 Case B

In Case B, the scenario is examined where the CLine as the engineered barrier system has degraded. The concrete is an effective barrier to the penetration of groundwater due to its low porosity and high compactness. However, the degradation of concrete gradually occurred over time due to the various physical and chemical factors. The hydrological characteristic change of the concrete as a function of time is not clearly defined [20]. Therefore, the scenario of Case B assumes at the time at which the function of concrete is perfectly lost. In other words, the engineering barriers are assumed to immediately degrade such that the radioactive nuclides are released into the groundwater. After the degradation of concrete, it was supposed that the hydrological parameters of CLine were converted to them of BFill.

Fig. 5(a) represents that the flow streamlines of the groundwater pass through the disposal facility in the event of the engineering barrier degradation. The degree of velocity is represented by the color scale bar, in which the red color is attributed to the fastest Darcy’s velocity. Additionally, Figs. 5(b) and 5(c) show that the Darcy’s velocity magnification of the X direction flow in the horizontal line and the Y direction flow in the vertical line in Fig. 5(a). Comprehensively, the groundwater mainly passes through the bottom region of the WCon domain and the upper area of the BFill. The flux of groundwater inside the silo significantly grows up overall in comparison with Case A because the porosity and permeability of CLine are increased. From the increased velocity, the migration of nuclide is expected to be mostly driven by advection instead of other mechanisms.

Fig. 5

(a) The streamlines of groundwater flow, (b) the flow velocity of red line in the X direction, and (c) the flow velocity of green line in the Y direction for the waste disposal facility at Case B.

The specific radioactivity distribution (surface plot) of nuclides and the streamline of groundwater flow (streamline plot) simulated at the time, which is the highest flux rate toward outline (the boundary of the right end), are shown in Fig. 6, respectively. Fig. 6 was sorted in ascending order according to the half-life of nuclides from (a) to (n). Overall, the result of the simulation shows that the radionuclides transport in the same direction as the flow of groundwater. The plume of nuclides spreads from the boundary of the WCon to the HRock and flows slightly down. This result suggests that the transport of the nuclide was mainly dominated by advection. The plumes are diluted into deeper groundwater flow systems as they decayed away. The spread area of the plume was radically decreased as the half-life of a nuclide is shorter. In addition, the remained nuclides in medium region grow up with the increasing their half-life. Especially, in the short half-life of radionuclides such as ⁵⁸Co (Fig. 6(a)) and ¹⁴⁴Ce (Fig. 6(b)), the remained nuclides at HRock may be neglected compared to the other nuclides by undergoing radical radioactive decay. These results mean that the spread area and the specific radioactivity of nuclide be critically attributed to the half-life of radionuclide.

Fig. 6

The nuclides release of the waste repository facility at 2.4 year: (a) ⁵⁸Co, (b) ¹⁴⁴Ce, (c) ⁵⁵Fe, (d) ⁶⁰Co, (e) ³H, (f) ⁹⁰Sr, (g) ¹³⁷Cs, (h) ⁶³Ni, (i) ¹⁴C, (j) ⁹⁴Nb, (k) ²³⁹Pu, (l) ⁵⁹Ni, (m) ⁹⁹Tc, and (n) ¹²⁹I.

Fig. 7 shows the specific radioactivity flux rate of each radionuclide at outline as a function of time. In all nuclides, their flux rate increases steadily until they reach the maximum and then exponentially decrease. The decrease of the boundary flux rate broadens as the half-life of a nuclide is long. This slower decrease rate suggests that the release process takes place much longer. To directly illustrate the relationship between the half-life of nuclides and the flux rate, Fig. 8 was replotted as the maximum of the boundary flux rate versus the half-life of nuclides. The replotted graph depicts three distinguishable regions such as short half-life (half-life < 1 year; ⁵⁸Co and ¹⁴⁴Ce), medium halflife (1 year < half-life < 500 year; ³H, ⁵⁵Fe, ⁶⁰Co, ⁶³Ni, ⁹⁰Sr, and ¹³⁷Cs) and long half-life (half-life > 500 year; ¹⁴C, ⁹⁴Nb, ⁹⁹Tc, ¹²⁹I, and ²³⁹Pu). In the medium region, the flux rate maximum of the released nuclide has a positive linear relationship with half-life. Otherwise, this graph exhibits two plateau sections with the independence of half-life at the short and long half-life range. This result suggests that the relationship between flux and half-life decreases or disappears out of range from medium half-life.

Fig. 7

The evolution of nuclides flux rate as a function of time.

Fig. 8

The maximum flux for the half-live of nuclides at the outlet of model.

3.3 Case C

To consider the solidified waste form, the simulation of Case C introduced the obtained equation (4) as the governing equation of the mass transfer coefficient for ⁵⁸Co, ⁶⁰Co, ⁹⁰Sr, and ¹³⁷Cs. Fig. 9 shows the release fluxes of each radionuclide from the WCon as a function of time. The value of the release flux was much lower (above 1/1,000 times) in the comparison with Case B. Although ⁵⁸Co and ⁶⁰Co had same the mass transfer coefficient, the flux rate of ⁶⁰Co with a long half-life increase the maximum value and leakage time in comparison with ⁵⁸Co. Next, the half-life of ⁹⁰Sr and ¹³⁷Cs have similar values. However, ¹³⁷Cs, which has a high mass transfer coefficient, has a higher maximum value and longer leakage time between two nuclides. Comprehensively, the released amount of ¹³⁷Cs was significantly larger than that of others, indicating the synergistic effect of a high mass transfer coefficient and long half-life for the leakage of solidification. Thus, the flux rate of the solidified waste form has a proportional relationship with the half-life and the mass transfer coefficient.

Fig. 9

The evolution of nuclides flux rate from WCon as a function of time in Case C.

Fig. 10 shows the specific radioactivity distribution of ⁵⁸Co, ⁶⁰Co, ⁹⁰Sr, and ¹³⁷Cs under Case B and C simulated at the same time when the release flux is the maximum in Fig. 9. The specific radioactivity of Case C represents a very low value in comparison with one of Case B. This result suggests that the solidification with the cement delays and reduces the leakage of nuclides as long as it maintains its integrity.

Fig. 10

The specific radioactivity distribution of Case B ((a)–(d)) and C ((e)–(h)) simulated at the time of the maximum flux.

4. Conclusion

This study has developed new models to describe radionuclide transport for near-field at the underground silo of Wolseong LILW disposal facility in Korea. The migration of radionuclides has been simulated by using COMSOL Multiphysics as a two-dimensional finite element numerical model code. Simulation scenarios were designed by the three scenarios as maintaining the integrity of the concrete (Case A), the degradation of the concrete (Case B), and the solidification of waste (Case C). In Case C, the solidification of waste was implemented with the mass transfer coefficient of the flux node at numerical modeling. The mass transfer coefficient was induced by experimental data, which was carried out by using a leaching test based on ANS 16.1 upon Co, Sr, and Cs.

In the presence of the concrete as the engineered barrier, the inside of the silo was calculated with the lowest Darcy’ velocity due to the low permeability feature of the CLine. So, the migration of nuclides is expected to be mostly driven by diffusion. Also, the diffused nuclides were kept inside the silo while the integrity of the CLine is maintained. On the other hand, in the case of the degradation of concrete, simulated results show nuclides transport in the same direction as the groundwater flow at the same time. This result means that advection is the main mechanism of nuclide transportation in Case B. From the analysis of the specific radioactivity flux rate maximum, we found that the release of nuclide has a positive linear relationship with a half-life in the range of medium half-life. As the solidification, the interrupt effect of the nuclide release has a reverse relationship with the half-life and the mass transfer coefficient. Also, this effect of delaying and reducing the leakage of nuclides as long as it maintains their integrity.

In addition to the solidified properties, the waste containers still have the unknown migration mechanism of radionuclides. Also, various reasonable experimental data are necessary for the analysis of nuclide migration. However, there are realistic limitations in acquiring all the data to simulate the interaction between nuclides and engineered barriers in the environment of a disposal facility. Due to these limitations, in this study, adsorption was excluded from the various transport mechanisms of radionuclides for a conservative assessment. Therefore, additional studies on the accumulated experimental data and the detailed safety assessment model are warranted to perfectly evaluate the performance of the waste disposal facility. In this respect, although our model is a safety evaluation model with limitations, we believe that it will be used as a basic model for safety evaluation models applying various interaction mechanisms.