Uranium is a natural radioactive element which is available in three isotope forms of U²³⁸, U²³⁵ and U²³⁴ in the environment. After World War II, uranium was mined extensively for production of nuclear weapons and later for nuclear power plants to produce energy [1]. In India, 2.7% of the total electrical energy is produced from nuclear power plants [2]. Uranium element is extracted from uranium ore and the uranium ore residues are deposited in ponds. These ponds are referred as uranium tailing ponds, which contains the radioactive decay products of uranium and heavy metals [3, 4]. The major radiation risk for living organism from uranium tailings are exposure to gamma radiation, windblown radioactive dust dispersion and radon gas and its progenies [2]. The major environmental impacts of uranium tailings are failure of tailing dams, exposure to radiation, surface and groundwater contamination due to leakage of radioactive and other toxic elements [3–5]. There are many studies executed earlier to treat various contaminants from water at various levels. Mercury sorption experiment has been performed by silica/carbon nanotubes and silica/activated carbon [6]. Various types of water treatment and recycling techniques have been discussed in terms of their basic principles, applications, costs, maintenance and suitability [7]. Saleh et al. [8] carried out a study on preparation and investigation of a multifunctional material that can be used in wastewater treatment for removal of arsenic. Further, Saleh [9] used Nanocomposite of carbon nanotubes/silica nanoparticles for adsorption of Pb(II). The activated carbon developed from waste rubber tires is identified as an efficient material for removal of Chromium from wastewater [10]. The experimental results demonstrate that the combining of silica and nanotubes is a promising alternative material, which can be used to remove the mercury from wastewaters [11]. Later, Zare et al. [12] conducted kinetic and thermodynamic studies for the efficient removal of radioactive uranium from solvent phase using nanoparticles, which exhibited excellent thermal stability and large surface area. Tawfik et al. [13] used poly-ethyleneimine modified activated carbon/Fe as an effective magnetic adsorbent to remove uranium ions from aqueous solution as a function of batch adsorption parameters; and it was concluded that the PAF magnetic sorbent could be considered as a promising and effective adsorbent for the purification of wastewaters from uranium ions. However, understanding the processes involved on the migration of radionuclide from uranium tailing pond to sub-surface is a necessary prerequisite for forecasting the long-term performance and prediction of groundwater quality.

Modeling the fate and transport of radioactive wastes has been known to be a valuable approach in accessing the radiological impact assessment of uranium tailing ponds in sub-soil system. The long-term impacts on release of uranium from tailing pond have been studied by many researchers using modeling approach [14–18]. Groundwater contamination due to leaching of uranium and its progenies from tailing ponds is a major concern. This may lead to unfit the groundwater for human consumption. The naturally occurring U, Th, Ra and Rn nuclides transport in groundwater are mostly controlled by physicochemical processes of weathering, decay and sorption at the water-rock interface [19, 20]. Rotter et al. [21] developed a biogeochemical reactive transport model for uranium immobilization in both single and dual-porosity porous media to study the bioremediation efficiency on uranium. Recently, a decay-chain transport of uranium in a heterogeneous anisotropic saturated porous medium is performed [22]. This study identified the necessity of monitoring the short-lived progeny radionuclides apart from their long-lived parents in the groundwater in the vicinity of uranium tailing ponds.

It is observed that many modeling studies are performed on transport of uranium in the saturated porous medium. The modeling studies to investigate the transport of radionuclide in unsaturated condition are very limited. Qian et al. [23] simulated the vertical transport of Sr⁹⁰ in the unsaturated Chinese loess under artificial rain conditions. Merk [24] studied the radionuclide transport from the concrete rubble when the radionuclide dissolves in the infiltrated water using HYDRUS-1D. Sanchez and Thorne [25] presented a mathematical model to study the behavior of U²³⁸ series radionuclide entering soils in solution and their uptake by plants along with the decay chain of uranium in soil-plant system under different hydrological regimes. Nair et al. [26] developed a multi-compartmental source term model to assess the radionuclide leaching in a saturated or unsaturated porous medium such as a uranium tailing ponds.

Vadose zone thickness can extend several meters and the soil hydraulic parameters are nonlinear [27–31]. Few studies indicated the vadose zone thickness can extend up to 50 m [32]. The flow and transport of radionuclide is caused by unsaturated matrix permeability, porosity and sorption coefficient [33]. Further, experimental and modeling study confirms that the radionuclide migration in unsaturated zone is highly dependent on unsaturated van Genuchten parameter [34]. The modeling as well as experimental studies in field scale also discussed earlier on the migration of radionuclide in large vadose zone thickness [35–37]. Moreover, the transport behavior of uranium from the tailing pond has extremely affected under various soil types such as sand, silt, clay, sediments and sandstone [19, 20, 23, 24].

The mill tailings are expected to include all the radionuclide in the ²³⁸U decay series after the uranium extraction from ore [22]. The ²³⁸U decay chain contains four long-lived radionuclide (²³⁸U, ²³⁴U, ²³⁰Th, ²²⁶Ra); five short-lived nuclides (²³⁴Th, ²²²Rn, ²¹⁰Pb, ²¹⁰Bi, ²¹⁰Po); and five very short-lived nuclides (^234mPa, ²¹⁸Po, ²¹⁴Pb, ²¹⁴Bi, ²¹⁴Po). On the other hand, ²²²Rn can escape from the ore and other short-lived progenies may decay to innocuous levels during the ore extraction processes. Thus, primarily uranium tailing ponds consists of long-lived radionuclide of ²³⁸U, ²³⁴U, ²³⁰Th, ²²⁶Ra. Subsequently, ingrowths of progenies will happen in the tailing ponds as well as in the percolated water.

Many modeling studies have been investigated on the migration of various radionuclides in the unsaturated porous medium. However, all long-lived, short-lived and very short-lived ²³⁸U nuclide transport in unsaturated porous medium has not been analyzed extensively. Since, the uranium tailing ponds are located in the soil surface a suitable numerical model has to be developed for vertical transport of multi-species radionuclide incorporating the nonlinear soil hydraulic parameters under unsaturated condition. Hence, the motivation of the present work is to develop a one-dimensional decay chain transport model to predict the concentration of radionuclide from uranium tailing pond due to vertical leaching under unsaturated condition. Further, the developed model is extended to evaluate the transport of radionuclide in different soils which are generally available beneath the uranium tailing pond.

Uranium tailings ponds are generally located in unsaturated zones [26]. The vertical transport by soil water, of radionuclide from uranium tailings ponds, is of major concern in assessing the vulnerability of groundwater contamination. Fig. 1 illustrates the conceptual understanding of uranium transport in an unsaturated zone.

The equation of mass conservation can be coupled with the Darcy’s law to produce the Richard’s equation which governs the water flow in an unsaturated porous medium [38–40] that is used for this study is given in Eq. (1):

Where: $C(h)=\frac{\partial \theta w}{\partial h}$ is specific moisture capacity (1/cm); θ_w is the volumetric water content (cm³/cm³); h is the pressure head (cm); K is the unsaturated hydraulic conductivity (cm/s); t is the time (s); z is the vertical coordinate (cm) positive downward.

In unsaturated porous medium, the part of pore space volume is filled by air. Hence, the moisture content θ_w is generally smaller than the porosity. The water pathway modeling in unsaturated porous system is considerably complicated as compared with saturated porous system. The water content (θ_w), pressure head (h) and hydraulic conductivity (K) are highly correlated in the unsaturated case. The hydraulic conductivity in unsaturated porous medium is not only depends on the matrix material but also on the local pressure head [24]. Widely acknowledged empirical relationships were arrived by van Genuchten [41] and described in Eq. (2) and (3).

where η is an empirical constant affecting the shape of the retention curve, θ_s is the saturated water content, θ_r is the and residual water content, α and β are both empirical constants.

The decay chain transport model for the parent radionuclide for a one-dimensional and homogeneous unsaturated porous medium can be represented as [22] expressed in Eq. (4):

where N₁ is the concentration of the parent radionuclide in the aqueous phase (atoms L⁻³ or mol); D is the hydrodynamic dispersion coefficient (L²T⁻¹) (D = α_L×v); α_L is the longitudinal dispersivity (L) can be produced by different pore sizes along the flow paths and/or by tortuosity; v is the pore water velocity (L/T) can be obtained by dividing the Darcy velocity by the effective porosity (which is expressed as v = q/θ_w); q is obtained from Darcy’s law $q=-K\left[\frac{\partial (h-z)}{\partial z}\right]$; λ₁ is the radioactive decay constant (T⁻¹).

Similarly, the decay chain transport equation for the i^th member of the decay chain can be written [22] as given in Eq. (5):

where i > 1 to M; M is the number of total nuclides involved in the decay chain. The ingrowth of the progenies from the previous parent radionuclide is included in the last term of Eq. (5). The sorption of radionuclides is described by linear and reversible. Hence, the retardation factor (R_i) of the radionuclide ‘i’, which can be expressed in Eq. (6):

where k_di is the distribution coefficient of radionuclide ‘i’ (L³M⁻¹); ρ_b is the bulk density of the porous material (M L⁻³). Fig. 2 illustrates the schematic representation of solving water flow and transport of radionuclides in an unsaturated porous medium.

The developed model is applied to a hypothetical uranium tailing pond. The soil properties used for this study are shown in Table 1. The radionuclides used in the decay chain transport model and their properties are given in Table 2. Further, the developed model is extended to various soil types and the corresponding soil hydraulic properties are presented in Table 3.

Fig. 5 provides the spatial variation of water content and the corresponding concentration profile of ²³⁸U after 1, 10 and 50 d in the unsaturated subsurface system in the presence and absence of sorption process. It is observed that the volumetric water content profile abruptly reduced from 0.207 (cm³/cm³) to nearly 0.1 (cm³/cm³) within the shallow depth (up to 50 cm) in one day. However, the water content profile remains constant up to 350 cm and 1,450 cm for 10 d and 50 d, respectively, i.e., the unsaturated system up to these depths continuously receives water and further reduces to the low water content (0.1 cm³/cm³). Further, it is observed from Fig. 5 that the movement of ²³⁸U seems similar to the profile of water content after one day in the absence of sorption process. The concentration of ²³⁸U starts at 0.129 Bq/L from the soil surface (uranium tailing pond) and is reduced to zero at 40 cm depth after 1 d. Conversely, significant delay is observed in ²³⁸U concentration profile (250 cm depth after 10 d and 900 cm depth after 50 d) as compared with the corresponding water content profile (350 cm depth after 10 d and 1,450 cm depth after 50 d). Additionally, the results of simulation indicate that there is strong delay in movement of ²³⁸U in the presence of sorption. This is because of the distribution capacity of ²³⁸U on soil particles (5.00 × 10⁻⁴ L/mg) which hold huge amount of ²³⁸U concentration as solid phase, which retard the movement ²³⁸U in the aqueous phase. Due to the sorption process the concentration of ²³⁸U becomes zero form 0.129 Bq/L from its source at 1.5, 2.5 and 5 cm depth after 1, 10 and 50 d, respectively. A closer observation between the presence and absence of sorption process implies that there is approximately 99% reduction in depth of penetration on ²³⁸U concentration with the sorption process. It can be concluded from Fig. 5 that the water content variation in the unsaturated soil is directly affecting the movement of ²³⁸U concentration. The presence of sorption leads to retard the movement of ²³⁸U concentration, which actually aids to minimize the ²³⁸U concentration in the deeper zone.

Fig. 6 shows the evolution of the water content for different infiltration velocity and its associated ²³⁸U transport in unsaturated soil in the presence and absence of sorption process for 10 d. The water content starts at 0.17 m³/m³ from the soil surface and gradually reduced to 0.1 m³/m³ at the depth of approximately 100 cm for the infiltration rate of 1 m/y with the saturated hydraulic conductivity of 2,907 m/y. Moreover, it can be observed from Fig. 6 that the water contents at the soil surface are 0.17, 0.207 and 0.22 m³/m³ when the infiltration rates are 1, 10 and 20 m/y, respectively. The higher infiltration rate eventually attempts the water content to move deeper in the unsaturated porous media. For example, as the infiltration velocity increases (10 m/y and 20 m/y), a relatively large amount of water content is retained in the soil surface (0.207 m³/m³), and subsequently the depth of movement of water content is also increased (maximum water content is observed up to the depth of 250 cm and 500 cm for the infiltration rate of 10 m/y and 20 m/y, respectively). Further, the results show that the concentration of ²³⁸U is reduced from 0.129 Bq/L at the soil surface to zero at the depth of 75 cm, 250 cm and 400 cm for the infiltration rate of 1, 10 and 20 m/y in the absence of sorption process. This shows the typical water content variation in the unsaturated soil is influencing the transport of ²³⁸U and the variation in infiltration is also impacting the depth of penetration. The similar trent is also observed for ²³⁸U transport in the presence of sorption process with considerable delay. This delay can be explained by the retradation induced by the solid phase of ²³⁸U during the sorption phenomina. It can be concluded from Fig. 6 that, the behavior of variation in the water content distributions during various infiltration rates is significantly affecting the movement of ²³⁸U in the unsaturated sub-surface system.

Fig. 7 provides the spatial distribution of ²³⁸U and its progenies in the absence of sorption process in the unsaturated sandy soil with the assumed infiltration rate of 10 m/y during 1, 5 and 10 y. It is observed from Fig. 7 that there is a maximum concentration (0.129 Bq/L) observed for ²³⁸U and ²³⁴U at the soil surface (uranium tailing pond) and it is further reduced to zero concentration as the depth increases. Conversely, all other progenies initially start from very low concentration (nearly 0 Bq/L) at the soil surface and gradually increase to their corresponding peak concentration with depth and finally become zero concentration followed by its asymptoticity. The order of peak concentrations is observed such as ²³⁸U = ²³⁴U > ²³⁰Th > ²²⁶Ra > ²³⁴Th > ²¹⁰Pb > ²¹⁰Bi > ²²²Rn > ²¹⁰Po. As soon as the concentration of ²³⁸U is available, the decay starts, gets converted into ²³⁴Th, and subsequently into ²³⁴U concentration. As ²³⁴U is abundantly available from the source at the soil surface; its concentration profile is similar to ²³⁸U. Further, the decay of ²³⁴U produces the remaining progenies in the decay chain reaction. The lowest ²¹⁰Po concentration is due to the last product in the decay chain and relatively its own high decay rate (1.83 y⁻¹). It is keenly observed from Fig. 7 that, the second lowest peak concentration of ²²²Rn due to the highest decay rate (6.6 × 10¹ y⁻¹) among all the uranium progenies in the uranium decay chain reaction. All the results present in Fig. 5 are independent of sorption process and the variations in the concentrations are observed only based on the respective decay rate. It is further observed that all the ²³⁸U and its progenies are available up to the depth of 800 m after 10 y and dangerous for the deeper zone in the unsaturated sandy soil. Thus it can be concluded from Fig. 7 that the concentration of ²³⁸U and ²³⁴U are continuously arrived the deeper depth from the soil surface (uranium tailing pond) for longer time interval and an enhanced peak concentrations are observed at a various depths (50, 250 and 550 m) at various time interval (1, 5 and 10 y) for all other ²³⁸U progenies except ²³⁴Th concentration.

The spatial distribution of ²³⁸U and its progenies during various intervals in the presence of sorption process in unsaturated sandy soil has been illustrated in Fig. 8. Since the sorption process of all the ²³⁸U and its progenies is considered in this simulation, a significant delay in movement and the reduction in peak are observed in all ²³⁸U and its progenies concentration. It is observed that there are approximately one to two order of magnitude reductions in peaks for all species except ²³⁸U and ²³⁴U from Fig. 7 and 8 at 10 y. This peak reduction is due to the consideration of sorption process. Moreover, the maximum concentration of ²³⁸U and ²³⁴U in the soil surface is also abruptly reduced when depth increases, but which is moving long depth in the case of absence of sorption process. It is further observed that the peak concentration of ²²²Rn is also reduced to one order even though the distribution coefficient is zero as compared to all other species. This reduction happens due to the distribution coefficients of the previous species which reduce its corresponding concentration in the aqueous phase of the decay chain and in turn, lead to the considerable reduction in peak ²²²Rn concentration. Essentially, the sorption process hinders the movement of all the species in the aqueous phase due to its corresponding distribution coefficient which helps to attract the species from the aqueous phase to the solid phase. The order of peak concentrations in the sorption process is similar as the absence of sorption process such as ²³⁸U = ²³⁴U > ²³⁰Th > ²²⁶Ra > ²³⁴Th > ²¹⁰Pb > ²¹⁰Bi > ²²²Rn > ²¹⁰Po except the order of reduction in peak concentration for ²³⁸U and its progenies. Further, it is observed from Fig. 8 that the concentration is virtually negligible after 1 y and 5 y as compared with 10 y. Moreover, the depth of penetration of ²³⁸U and its progenies is reduced to less than 10 m and the peak is observed at nearly 2 m depth due to the presence of sorption process, whereas the ²³⁸U and its progenies move up to 800 m and the peak is observed at nearly 550 m depth when the sorption process is not considered after 10 y. Since, the half-lives of the parent radionuclides are very long [22], it is great possibility that the concentration of ²³⁸U and its progenies can move in deeper depth even though the sorption process is considered.

Fig. 9 provides the spatial distribution of ²³⁸U and its progenies in the absence of sorption process in the unsaturated silty soil with the assumed infiltration rate of 3 m/y during 1, 5 and 10 y. The water content is varied from 0.438 m³/m³ at the soil surface and gradually reduced to 0.099 m³/m³ at the depth of approximately 12 m, 50 m and 100 m (not shown in the Fig. 9) for the infiltration rate of 3 m/y with the saturated hydraulic conductivity of 22 m/y after 1, 5 and 10 y, respectively. It can be noted that the concentration of ²³⁸U and ²³⁴U is following similar trend with maximum inflow concentration of 0.129 Bq/L at the soil surface and reducing to zero at 10, 45 and 90 m depth after 1, 5 and 10 y, respectively. As compared with the sandy soil (Fig. 7) the concentration movement with respect to the depth is significantly reduced in the silty soil, i.e., the concentration of ²³⁸U and its progenies are become zero at the depth of 100, 450 and 800 m in the sandy soil whereas it is reduced to 10, 45 and 90 m in the silty soil (approximately 90% reduction in movement). Similarly, the peak concentrations are also observed in the shallow depths in the silty soil as compared with sandy soil (For example, the peak concentration of ²³⁰Th is 50 m in the silty soil and 550 m in the sandy soil for 10 y). This reduction in movement in silty soil is primarily due to the decrease of two orders in hydraulic conductivity and other changes in soil hydraulic parameters such as water content, pressure head and van Genuchten parameters (Table 3) as compared with the sandy soil. However, the order of peak concentrations is relatively similar which is observed in the sandy soil such as ²³⁸U = ²³⁴U > ²³⁰Th > ²²⁶Ra > ²³⁴Th > ²¹⁰Pb > ²¹⁰Bi > ²²²Rn > ²¹⁰Po. Further, it is keenly observed that the maximum concentration of ²³⁸U and ²³⁴U are reducing with respect to depth in the silty soil, while the same (maximum) concentrations are observed for larger depths in sandy the soil (Fig. 7). For instance, the concentration of ²³⁸U and ²³⁴U are constant up to the depth of 250 and 550 m and slowly becoming zero in the sandy soil (Fig. 7), whereas the maximum concentration starts declining immediately just below the soil surface (uranium tailing pond as source) in the silty soil (Fig. 9) due to its corresponding soil hydraulic properties. Thus, it can be concluded from Fig. 9 that the ²³⁸U and its progenies in the absence of sorption process in the unsaturated silty soil is significantly reduced its depth of movement predominantly due to the reduced unsaturated hydraulic conductivity and other soil hydraulic parameters such as water content, pressure head and van Genuchten parameters, which in turn, the enhance availability of ²³⁸U and its progenies are in the shallow depths. However, the concentration of ²³⁸U and its progenies are speared up to 90 m depth after 10 y, which shows that the groundwater aquifers may get affected if the water table is shallow.

Fig. 10 describes the spatial distribution of ²³⁸U and its progenies in the presence of sorption process in the unsaturated silty soil. The ²³⁸U and ²³⁴U concentrations are released from the uranium tailing pond at the soil surface which have higher concentration (0.129 Bq/L) and become zero at 1 m after 1 y but very marginal downward movement is observed after 5 and 10 years. Similarly, other progenies of ²³⁸U are also shown a low depth of transport in the unsaturated silty soil. This significant retardation in movement on ²³⁸U and its progenies is predominantly happened by the sorption process of all the uranium species and the low hydraulic conductivity of soil. Approximately 90% reduction in the depth of penetration is observed in the sorption process when compared Fig. 10 with 9; (zero concentration is observed nearly at 1, 3 and 5 m in the presence of sorption process and 10, 45 and 90 m in the absence of sorption process after 1, 5 and 10 y); on the other hand almost 99% reduction in depth of penetration is observed in silty soil (Fig. 10) as compared with the sandy soil (Fig. 8) in the presence of sorption process (i.e., zero concentration is observed nearly at 1, 3 and 5 m in the silty soil and 100, 450 and 800 m in the in sandy soil with the sorption process after 1, 5 and 10 y). Moreover, it is observed from Fig. 9 and 10 that the peak concentration of all the uranium species except ²³⁸U and ²³⁴U are reduced nearly one order in the presence of sorption process as compared with the absence of sorption phenomena. It can be concluded from Fig. 10 that the concentration of ²³⁸U and its progenies vanish at very shallow depths (approximately before 10 m) and offer a positive impact to stop the groundwater contamination.

In the present study, a numerical model is developed to describe the transport of ²³⁸U and its progenies in the vadose zone from uranium tiling pond situated on the soil surface. The one-dimensional coupled water flow and ²³⁸U and its progenies decay chain reaction is developed using implicit finite difference scheme. In addition, this study is also extended to perform the movement of ²³⁸U and its progenies concentration in different soils such as sand and silt. From the developed model, the following conclusions are arrived:

In over all, the presence of sorption process in the silty soil significantly helps to avert the groundwater contamination, especially at deeper depths and for larger span. Future studies can be intended to integrate other mechanisms such as evaporation, root uptake, surface runoff, degradation in uranium pond and other climatic conditions to achieve the combined effect of all the coupled processes on the equivalent system to develop a comprehensive model.