Long-Term Evolution of Uranium Mobility within Sulfated Mill Tailings in Arid Regions: A Reactive Transport Study

: Management of mill tailings is an important part of mining operations that aims at preventing environmental dispersion of contaminants of concern. To this end, geochemical models and reactive transport modeling provide a quantitative assessment of the mobility of the main contaminants. In arid regions with limited rainfall and intense evaporation, solutes transport may signiﬁcantly differ from the usual gravity-driven vertical ﬂow. In the uranium tailings of the Cominak mine (Niger), these evaporative processes resulted in the crystallization of gypsum, and to a lesser extent jarosite, and in the formation of surface levels of sulfated gypcrete, locally enriched in uranium. We present a fully coupled reactive transport modeling approach using HYTEC, encompassing evaporation, to quantitatively reproduce the complex sequence of observed coupled hydrogeochemical processes. The sulfated gypcrete formation, porosity evolution and solid uranium content were successfully reproduced at the surface and paleosurfaces of the tailing deposit. Simulations conﬁrm that high solubility uranyl-sulfate phase may form at the atmospheric boundary where evaporation takes place, which would then be transformed into uranyl-phosphate phases after being watered or buried under fresh tailings. As these phases usually exhibit a lower solubility, this transition is beneﬁcial for mine operators and tailings management.


Introduction
Management of mill tailings and mine waste constitutes a sensitive operation with potential serious environmental hazards [1][2][3]. This requires mine operators to perform careful monitoring [4] and risk assessment regarding the storage and management of the waste. In particular, mill tailings originating from uranium mines pose an additional radiological hazard, in part due to potential radon exhalation [5] and subsurface contamination [6][7][8][9].
The physical and chemical nature of these materials depend not only on the nature of the ore but also on the chemical extractive processes as well as the storage conditions. While waste rock typically consists of rock with larger grain sizes depleted in the elements of interest, mill tailings result from the fine crushing and chemical processing of the ore. Sulfuric acid leaching is widely used for uranium due to its high solubility in acidic conditions, which is further increased by the formation of sulfated uranyl aqueous complexes [10][11][12][13]. The interactions between uranyl and sulfate ions in solution therefore act as one of the major control mechanisms for uranium mobility. The resulting material consists of a complex structure combining ore minerals, partially dissolved by chemical treatment, and secondary minerals that depend on the processing method (see [14][15][16][17][18] for instance). Mill tailings thus consist of complex and reactive objects whose characterization is challenging due to the large variety of secondary mineral phases, which are often amorphous or below typical detection limit. In addition, uranium migration in subsurface levels depends on hydrobiogeochemical [19][20][21] and climatic factors. In arid regions with limited rainfall and significant evaporation, several hydrogeochemical processes interact to result in complex migration pathways.
Déjeant et al. [14] studied the evolution of uranium tailings from the Cominak mine (Niger). Cominak was an underground mine whose operation started in 1978, with a uranium grade of 4300 ppm and a tonnage of 44,000 tons, with uraninite as the major bearing phase of a sandstone deposit. The ore processing mainly consisted in the crushing and chemical treatment with sulfuric acid, with a recovery rate between 90% and 96% [15]. The generated mill tailings were deposited as a pile, which was initially watered. After 10 years, Déjeant et al. [14] highlighted uranium reconcentrations associated with evaporation processes near the atmospheric boundary. Surface formation of salt-enriched crusts (gypcrete, calcrete) have been observed under other evaporitic circumstances [22][23][24][25] and are sometimes also associated with uranium [26].
Reactive transport modeling has been successfully used in several mining applications, e.g., the prediction of acid mine drainage [27][28][29][30][31][32][33], the study of uranium migration in the context of mine remediation [17,[34][35][36][37][38][39][40][41] or for uranium in situ recovery [42,43]. However, modeling evaporitic systems has been challenging for reactive transport simulators, which usually do not include a water balance. Evaporation leads to a loss of water, resulting in a local increase in the solute concentration, potentially reaching evaporite's solubility. This cannot be reproduced by the typical Darcy approach used by most simulators [44] and requires a water-vapor phase change. Gran et al. and Nakawaga et al. [45] have provided a non-isothermal approach to model evaporation processes by modifying the traditional soil water characteristics. After 12 days of evaporation, their results yielded a narrow zone (≈5 cm) of enriched fluid chemical compositions. Lahrouch et al. [46] also provided an experimental and numerical study to model evaporitic systems and mostly focused on the aqueous transport of cations. The method developed in [47] is suited to describe the three phase equilibrium of water and can therefore be used to simulate evaporitic processes.
In this research paper, we use HYTEC [48][49][50] to simulate the evolution of uranium mill tailings studied by Déjeant et al. [14,51]. Reactive transport simulations will provide a quantitative and thermodynamic assessment of the long-term mobility of uranium within sulfated mill tailings. Simulations will be calibrated on observed and published data of the tailings from the Cominak mine. Section 2 presents a summary of the observed features within the mill tailings, as well as the numerical and modeling approach. Section 3 presents simulation results, while discussion and concluding remarks are given in Sections 4 and 5.

Geographic Zone of Interest
The tailings that were studied originate from the uranium mine of Cominak (CK), Niger (18°45 56 N, 7°20 23 E). The uranium deposits belong to the Tim Mersoi sedimentary basin and are located in the Arlit region of the Sahel desert. This zone is extremely arid with a mean annual rainfall of 100 mm/y. The potential evaporation rate, e.g., during the watering of mill tailings, is very high due an annual mean temperature above 30°C and the dryness of the air. Evaporation rates in the range of 3000-5000 mm/y have been observed under such conditions [52][53][54].

Materials
The tailings from the CK mine result from the crushing and acidic leaching of the ore. The mineralogical evolution from the ore to the stored tailings was already assessed in previous papers [14,15,51,55]. The next section provides a summary of the observed mineralogy, which constitutes the starting point towards the construction of our model.

Mineralogy
Déjeant [51] identified the major mineral phases inherited from the ore: quartz, kaolinite and microcline. Furthermore, minor refractory phases were identified, such as Fe-oxides, apatite and sulfide minerals. Uranium is mostly found within uraninite and coffinite relicts. Part of the ore minerals dissolved during the chemical processing step. Despite the acidic and oxidizing conditions, U(IV) oxides and sulfide minerals persisted, indicating their stability and low reactivity. After being stored on site, secondary hydrated sulfate phases were observed in the tailings, mainly gypsum and, to a lesser extent, jarosite. The high amount of sulfuric acid used in the chemical treatment allowed these sulfate phases to precipitate.

Tailings Evolution
Déjeant et al. [14] described the vertical distribution of uranium in the tailings after 10 years of disposal. Figure 1a shows the evolution of the vertical distribution of the uranium content [15], with visible significant enrichment at the top level, close to atmospheric conditions, as well as in the paleocrust at a 15 m depth. These enrichment levels correspond to low porosity levels ( Figure 1b). Evaporation-promoted local mineral reprecipitation near the atmospheric boundary resulted in a sulfate-cemented crust, with a porosity reduction from 45 to 35% (Figure 1b), mostly consisting of gypsum that precipitated in the pore space ( Figure 1d). Several other phases were observed in smaller quantities, such as a magnesium-sulfate phase [51]. These levels also corresponded to the levels of enriched elemental composition, particularly in S and P (Figure 1c). U speciation has been characterized at these locations using SEM imaging [14,55]. In the sulfated-cemented crust, uranium is mostly found as a uranyl-sulfate phase. while it mostly consists of uranyl-phosphate phases in the paleocrust (Figure 1e,f). Uranyl phosphate minerals are commonly encountered as alteration products of uraninite and coffinite in oxidizing conditions [56][57][58][59][60].

Methods
In this section, we detail our reactive transport approach from the mathematical and numerical description to the chosen reactive system.

Mathematical Description
Most reactive transport algorithms do not include water conservation. Hence, removing water is often accompanied with removing a solute. In an evaporitic system, the removal of water leads to a concentration of dissolved metals that may reach the solubility of evaporites.
In general, the mass balance equation of water can be written as: where n v and n l denote the mass (in mol·dm −3 bulk ) of water in the gas and aqueous phases, Q r denotes the chemical reaction term (e.g., consumption or release of water molecules by minerals), ω is the porosity, and ρ l,n denotes the molar density of water (mol·dm −3 liquid ). This equation is based on the hypothesis that water vapor is subject to diffusion (D g ), while liquid water only moves by advection following the Darcy velocity − → u l ·g 0 denotes the vapour concentration within the gas phase (mol·dm −3 gas ). Darcy velocity is written as: with k r,l the relative permeability (using van Genuchten model), K the hydraulic conductivity (m ·s −1 ) and h the hydraulic head (m). This equation can be divided into two parts: where the first equation corresponds to the traditional Richards equation and the positive term Q evap stands for water evaporation. Most reactive transport algorithms decouple the resolution of the Richards equation from the reactive transport problem. For water and solutes, the reactive transport problem is written as [47]:  The top one is in contact with the atmospheric condition, which is subject to evaporation. The lower one corresponds to the old atmospheric boundary, which had been subject to evaporation for ≈10 years. Redrawn from [14]. (b) Porosity profiles showing two levels of reduced porosity corresponding to the levels of uranium enrichment. Redrawn from [15]. (c) The elemental profile of S and P within the tailings.
In the latter, T i,l describes the total concentration of component i within the liquid phase (mol·dm 3 liquid ), S l , S g denote the fluid phase saturations, and D = D l + α|u| combines the effective diffusion (D l ) and dispersion (α) coefficients] in the liquid phase. Q ri is the chemical source/sink term constrained by the geochemical system and describes the mineral dissolution/precipitation rates, as well as desorption/sorption processes.
The approach from Seigneur et al. [47] is suited to the description of such systems, as water is part of the geochemical system. One additional level of coupling is required to link the water vapor concentration, i.e., the relative humidity RH, to the liquid capillary pressure [61]: where R is the ideal gas constant, T is the temperature, and p cap is the capillary pressure. The latter is computed using empirical relations such as the common Van Genuchten law [62]: where p (VG) e is the air entry pressure. Hence, the unsaturated flow problem is solved traditionally. Then, the two-phase reactive transport problem is solved using the approach given in [47], with the water vapor concentration being dependent on the liquid saturation. Then, source terms for water are computed to solve the unsaturated flow problem at the following timestep. Given the fast dynamics of vapor diffusion and the explicit treatment of the source terms, small timesteps are required.

Evaporation and Flow
To model the evaporation, a constant dry atmospheric boundary is modeled at the top. As no liquid water is allowed to be exchanged with the atmospheric boundary, only evaporation and diffusion of water vapor can occur near the atmospheric boundary. Within the pile, water movement is driven by capillarity in the liquid phase and diffusion in the gas phase. In this system, water evaporation near the top decreases the water content, which triggers an ascending flux towards the atmospheric boundary by capillarity. Solutes are then concentrated as water evaporates, and minerals reach their solubility saturation limit and, consequently, precipitate, which leads to the observed decrease in porosity.
While the simulation of the capillary ascension of water is straightforward, its full description with the geochemistry is numerically challenging over the whole column for multiple reasons. First, hydraulic parameters (permeability, unsaturated parameters) are poorly constrained. Second, it is difficult to assess to which levels the volumetric water contents decrease. For very low water content, a decrease in relative permeability prevents the capillary flux from being maintained over a long period of time. Hence, volumetric water contents may decrease to near-zero values, which constitutes a numerical challenge for any geochemical code. For these reasons, we focus on a smaller (3 m) column, for which the diffusive properties of water are set to trigger a 2 m/y ascending flux with a constant composition at the low boundary. Furthermore, within 10 years, such a flux corresponds to the total amount of water initially available if the whole pile was initially saturated in water.
Simulation parameters are given in Table 1. The high hydraulic conductivity K was chosen mainly to compensate the very low relative permeability induced at very low water content (saturation below 10%). Due to limited data to constrain the Van Genuchten parameters, they were chosen to be comparable to coarse sandstones. At the top of the domain, a constant-relative humidity of 2% was chosen in an atmospheric cell with zero permeability, preventing liquid water from being exchanged with the atmosphere and only allowing evaporation at the top. Gas diffusion was increased by a factor 10 at the top to induce local evaporation near the atmospheric boundary.

Parameter
Value Unit

. Geochemical System and Thermodynamic Database
In this study, we used the Prodata thermodynamic database (v1.4) dedicated to environmental geochemistry in mining contexts [13] in combination with the truncated Davies model of activity corrections of dissolved ions. This relation is indeed applicable up to moderate ionic strengths I of about 0.5 M, which was never reached in the model, except in the upper desaturated zone close to the atmosphere (I ≈ 0.75 M). Higher ionic strengths were prevented due to the control of water chemistry by the precipitation of solid phases. Ca was controlled by the formation of gypsum. Similarly, the mobility of Mg is also controlled by the observed precipitation of sulfate-magnesium phases. For the sake of simulation simplicity, gypsum will be included as a surrogate for these Ca-Mg phases. Regarding alkalis (represented by K in our model), limited concentrations are modeled by their incorporation in the uranyl-sulfate phases and jarosite. The sulfate concentration is controlled by the neoformation of all these minerals.
To represent the tailing mineralogy, we have chosen the minerals described in Table 2, with composition and thermodynamic data obtained from Prodata [13]. Other phases were neglected, either due to their apparent inert character (e.g., pyrite encapsulated in quartz grains) or due to their limited influence on simulation results (e.g., kaolinite).
With respect to the uranium minerals, an analysis of the observed uranyl-phosphate phases in the Cominak tailings revealed an atomic U:P ratio of 1 [55]. Hence, this ratio eliminates autunite. Several UP minerals with U:P = 1 were investigated in the database, and the solid UO 2 HPO 4 :3H 2 O(cr) was selected. Regarding uranyl-sulfate minerals, hydrated uranyl-sulfate minerals (UO 2 SO 4 :nH 2 O) can be found in the Prodata database [13]. However, these minerals are very soluble and unlikely to form. Zippeite phases have been found as secondary mineral phases in mill tailings and wastes [55,63]. Zippeite are hydrated uranyl sulfates M(UO 2 ) 4 (SO 4 ) 2 O 2 (OH) 2 :(H 2 O) 4 containing monovalent (K + , Na + ) or divalent (Mg 2+ , Ca 2+ ) cations that generally form under low pH conditions. Sharifironizi et al. [63] provided thermodynamic data for the K-zippeite (see Table 2). This mineral was considered as a proxy for uranyl-sulfate formation inside the gypcrete.
The chemical composition of the tailing water was poorly known. Values of pH around 4 were nevertheless measured, which is consistent with the fact that the tailings were not neutralized with lime after the sulfuric acidic U extraction process. The model considered several reasonable hypotheses to set the total concentration of the main cations and anions. The initial concentrations in calcium and sulfates were fixed by the equilibrium with gypsum and jarosite-the two main sulfate minerals that were clearly identified in the tailings. The initial concentrations in the minor species silica and phosphates were fixed by the equilibrium with hydroxyapatite and quartz, respectively. For uranium, initial total dissolved concentrations were set as 50 mg/L, within the range of observed concentrations between 10 and 500 mg/L.
Atmospheric O 2 gas was assumed to be uniformly present in the tailings columns due to the strong water desaturation of the system, allowing for air diffusion, as well as the lack of any vegetation cover at the tailing surface and of any organic matter inside the tailings. Therefore, the redox potential was constantly oxidizing, which favors the highest redox states of uranium and iron (U(VI) and Fe(III), respectively) of the mineral phases.
Sensitivity analyses were performed on the initial mineralogy and reactivity description. The dissolution of apatite, jarosite and quartz was modeled under kinetic constraints during the reactive flow modeling. A kinetic description of uraninite was also motivated by observations of the remaining crystals in the tailings, which could not be obtained by describing it as reacting at equilibrium since the redox conditions were oxidizing. We used the following simple kinetic law for the dissolution of a mineral M (concentration (M)) where k is the intrinsic rate constant (mol·m −2 ·s −1 ); A v the volumetric surface (m 2 ·L −1 ) that can be calculated from the specific surface A s (m 2 ·g −1 ); Q the species activity product; and K the solubility product. Table 3 reports the kinetic parameters used in the simulations. The specific surfaces of quartz, apatite and uraninite are compatible with the particle size produced during the ore processing, whereas jarosite is a secondary phase that formed during the tailing alteration. A first estimation of the kinetic rate constants was performed from the wellknown compilation of Palandri and Kharaka [64]. However, the rate constants of apatite and uraninite had to be decreased by two orders of magnitude in order to best simulate the measured bulk chemistry of the P and U profiles. It is a common observation that in situ kinetic rates are generally lower than the rates obtained for fresh materials under lab conditions. Among other reasons, passivation by protecting diffusive layers can explain such a decrease, as discussed in Section 3.3 for the apatite grains.

. Simulation Descriptions
In this study, four reactive transport simulations with increasing complexities are performed to understand the hydraulic and geochemical interactions at play within the Cominak tailings ( Figure 2). First, we look at the evolution of a representative batch geochemical systems under evaporating conditions, with an emphasis on pH, porewater composition and secondary mineral saturation indices (Section 3.1). Second, we focus on the formation of the gypcrete at the atmospheric boundary through modeling the capillary rise and solute concentration due to near-surface evaporation (Section 3.2). Third, we look at the apatite evolution within our system. While apatite is fairly soluble in acidic conditions, modeled kinetic rates are relatively low. We model the progressive dissolution of apatite grains in a geochemical environment comparable to the gypcrete (Section 3.3). To this end, the reactive transport simulation results from the gypcrete formation are used to constrain the pore water chemistry, which is in contact with the apatite grains. This allows assessing whether the precipitation of uranyl phosphate plays a passivating role, providing a physicochemical explanation for the rather limited apatite dissolution rates. Finally, we study the evolution of the gypcrete after being buried under a layer of fresh tailings (Section 3.4). Final conditions of the gypcrete formation simulations are used as initial conditions, with fresh tailings put on top, as schematically shown in Section 2. This has the effect of suppressing evaporation at the former gypcrete location. This new simulation follows the same methodology and parameters as simulation 2 and is used to predict and understand the long-term mobility of uranium.

Geochemical Evolution during Evaporation
Batch (0D) simulations were performed to investigate the impact of evaporation on the geochemical system. The results are displayed in Figure 3 and reveal that only considering evaporation may lead to significant acidification by increasing the proton concentration, lowering the pH by 1-2 pH units. However, since the coupling between jarosite dissolution and ferrihydrite precipitation may constitute a powerful buffer couple, and it was included in the geochemical model. Furthermore, several secondary uranium-sulfate phases were considered from the Prodata database (v1.4) that never reached saturation, even after a reduction in the water content of a factor 100. This is due to the increasing importance of the aqueous complex UO 2 SO 4 , which prevents reaching the saturation of these phases. Conversely, zippeite solubility could be reasonably reached, partially due to the presence of jarosite-K.

Gypcrete Formation
Simulation results of the coupled evaporation-transport and reactions within the tailings are given in Figures 4 and 5. The mineral system and reactivity are described in Tables 2 and 3.   Figure 4a shows the spatial profile of the upward Darcy velocity induced by evaporation. The velocity vanishes at the top due to evaporation, and the in-depth profile is uniform, which is consistent with observations and models obtained by Nakagawa et al. [46]. Figure 4b,c shows the hydraulic evolution of the system: water content progressively decreases in the system, while porosity progressively decreases near the top atmospheric boundary. The decrease in porosity is associated with the reprecipitation of gypsum, as represented in Figure 4e. Hence, evaporation acts in favor of precipitation in two ways, combining a steady arrival of solutes by the induced ascending flow rate and a local concentration. The pH is well buffered around 4-4.5 over the whole duration (Figure 4d). Figure 5 provides a temporal analysis of the geochemical evolution in the node just below the atmospheric boundary. While the temporal evolution may not be representative of reality, as the modeled processes may have occurred on a shorter timespan, the main observed processes are well reproduced by the model. The simulations show a progressive decrease in porosity and a significant increase in U content in the solid phase to reach values over 4000 ppm, which is comparable to what was measured experimentally. The relative importance of zippeite and uranyl-phosphate (Table 2) minerals is strongly dependent on the amount of available phosphates, which is controlled by the kinetic description of apatite. Figure 5f shows a steady increase in U content in the solid phase with time. It must be noted that this evolution will eventually be limited by the total amount of water that can be evaporated. Considering a total height of saturated tailings, such a U enrichment cannot be maintained for longer than 10 years. Our simulation is able to quantitatively reproduce the observed behavior of the tailings.

Apatite Dissolution and Passivation
We simulated the fate of apatite grains of different sizes within a geochemical environment comparable to what is shown in Figure 5h. The evolution of apatite was simulated under water-saturated conditions and local diffusive transport for the sake of simplicity. The objective of such a simulation was to estimate the relative chemical stability of sulfate and phosphate solid phases and to take into account the porosity evolution induced by the dissolution and precipitation of these minerals. The results are shown in Figure 6, indicating that small grains can dissolve and be replaced by the precipitating UP crystals due to the lower solubility of uranyl-phosphate minerals. Conversely, due to the larger molar volume of uranyl phosphates, concomitant dissolution and precipitation yield a decrease in porosity, inhibiting diffusive transfer and thus the complete dissolution of the bigger apatite grains. These results are in good agreement with SEM observations (Figure 1e,f) and may explain the slow reactivity of apatite in the present system.

System Evolution after Burying the Tailings
The evolution of uranium mobility under the influence of evaporation was treated in Figure 4. After burial under fresh tailings, the local influence of evaporation decreases, and we simulate the fate of uranium mobility in these conditions. The results are given in Figures 7 and 8. Figure 7 shows the geochemical evolution of a primary gypcrete after burial. Without the influence of local evaporation, the ascending capillary flux progressively dissolves uranyl-sulfate evaporites. The released aqueous uranium is transported upwards, which slightly displaces the U-enrichment. At the primary gypcrete location, the U content in the solid phase decreases from 4000 to 2500 ppm, close to the measured value (Figure 1) in the paleocrust, with a slight enrichment above, also visible in Figures 7 and 8. During this transition, the uranium concentration and mobility are significantly reduced and controlled by the solubility of UP phases. Depending on the phosphate availability, the dissolved U concentration may vary. During this transition, the simulations suggest that porosity remains relatively stable because gypsum is stable throughout the column. This steady porosity therefore constitutes a paleomarker of atmospheric contact in evaporitic systems [22].

Discussion
The main objective of the simulations presented in this study was to quantitatively verify the hydrogeochemical scenario provided in [14,51], which explains the formation of surface U-enriched sulfated gypcrete. Furthermore, these observations in an evaporative system were used to test the numerical approach implemented in HYTEC for the fully coupled reactive flow system [47].
The applied methodology for simulation calibration followed the philosophy: as simple as possible, as complicated as necessary. The first and biggest challenge to set these simulations was the flow calibration. For simplicity, a smaller vertical column, focusing on the three first meters, was chosen to limit numerical instabilities. The ascending flow rate was chosen as representative of the evaporation rate and in order to not overestimate the total quantity of available initial water. With these simplifications, our simulations present some limitations, mainly regarding the temporal system dynamics. It is possible that the modeled processes occurred on a shorter timespan and were interrupted by the lack of available water. Furthermore, the interplay between the rare rainfall, water spray from mine operators and evaporation is challenging to constrain. All of these drawbacks should limit how representative the simulations are.
Second, a representative geochemical description was constrained as much as possible by thermodynamic data [13] and detailed observations from [14,51]. Selecting a relevant uranyl-sulfate phase in the database proved to be a challenging task due to the relatively large solubilities of uranyl-sulfate phases. Nevertheless, Sharifironizi et al. [63] described zippeite, a uranyl-sulfate mineral that had been identified in uranium mine waste, which was able to precipitate during evaporation in our geochemical system. Several minerals were then included in the simulation to ensure a pH buffering capacity and prevent the evaporation-induced acidification. The interactions between jarosite and ferrihydrite, as well as the apatite dissolution, were the main contributors to buffering the pore water. Furthermore, the precipitation of zippeite and the reactivity of gypsum and jarosite exerted significant control of the mobility of the major elements in the model, limiting the development of highly saline solutions.
Calibrated as such, the simulations yielded results close to what was observed in the field. Sensitivity analyses were performed on the initial porewater chemistry and mineral reactivity (apatite, uraninite) with limited impact on the simulated processes. Modeled porosity reductions were always between 8% and 12%, while uranium enrichment varied between 4000 and 5000 ppm. Therefore, the uncertainty of most model parameters would not yield significantly different results. One of the most impactful parameters is the reactivity of apatite grains, which controls the amount of dissolved phosphates and the subsequent precipitation of uranyl phosphates and also constitutes a buffer for pH. Considering the wide range of apatite grain sizes on the one hand and the potential observed and simulated passivation process on the other makes the reactivity description challenging.
Despite the limited sensitivity of the modeling results to changes in the initial mineralogy or porewater chemistry, changes in the initial dissolved sulfate or uranium concentration constitutes an exception. Too low of an estimate for the initial dissolved sulfate and/or the uranium concentration may prevent the saturation of uranyl-sulfate minerals (zippeite-type) upon evaporation and thus prevent the formation of uranium-enriched gypcrete. This could potentially explain the different behavior of another mine from Niger. Indeed, the SOMAIR mine [15] displayed no levels of U enrichment. Considering that these two mines were operated in a very similar way and that they are located close to each other, this different behavior will probably find its roots in the initial geochemical conditions. The simulations suggest that the initial dissolved uranium concentrations and residual uraninite reactivity could be the primary controls of the overall tailings' evolution.
These observations may provide interesting feedback regarding the choice of chemical processing used for the ore. In both cases (Cominak vs. Somair), the long-term uranium mobility is likely controlled by uranyl-phosphate phases.
From an operational perspective, these simulations allow confirming several interesting points for the long-term management of mill tailings. The gypcrete layer at the surface of the pile limits the aerial dispersion of tailings particles. Second, the transformation of sulfated uranyl minerals into phosphate uranyl minerals, known to exhibit lower solubilities even in acidic conditions, plays in favor of uranium immobilization. Last, the reduction in porosity observed within the tailings, which were well reproduced by our simulations, may limit radon exhalation.

Conclusions and Perspectives
This paper presented a reactive transport approach to model the complex coupled hydrogeochemical processes that control the evolution of mill tailings in arid regions. The mineralogy was constrained based on observations, and literature values were adopted to describe the transport and hydraulic properties. A simplified description that focused on the levels of interest was adopted to facilitate the convergence of simulations at a very low water content. Our model was able to quantitatively reproduce the main observed phenomena. In particular, the formation of a gypcrete enriched in uranyl sulfate was predicted with relatively good accuracy regarding the porosity decrease and solid U content. Furthermore, the evolution of such gypcrete after being buried matched the observations: steady porosity and the conversion of uranyl sulfate into uranyl phosphates. The role of the initial sulfate availability was investigated and showed to be a major controller of uranium mobility in this industrial mining environment.
From the point of view of mine operators, these simulations confirm that the long-term controls of uranium mobility are the low-solubility phosphate phases. In addition, these simulations provide solid thermodynamic and quantitative support to understanding how complex materials, such as mill tailings, may evolve in such arid environments.
The use of reactive transport modeling in systems encountering heavy evaporation is still limited, partially due to numerical challenges associated to the three-phase water conservation and the potential water disappearance. In this paper, we have shown that reactive transport simulators can be extended to such systems. With some additional work to properly take into account the wetting/evaporation cycles as well as the hydrogeochemical heterogeneities, this may open up new perspectives for reactive transport to model evaporitic systems, such as salars. Finally, this paper has shown that reactive transport modeling can be a powerful tool in the post-mine context.