Location of Uranium‐Rich Brines Determines the Distribution and Grade of Unconformity‐Related Uranium Deposits

Recent fluid inclusion analysis indicates that uranium‐bearing brines were present in both the sandstone and the basement at the time of ore genesis in the Athabasca Basin, Canada. However, the question of how the location of uranium‐rich brines controls the formation of unconformity‐related uranium deposits remains unresolved. In this study, four reactive flow modeling scenarios are designed to address this outstanding issue. Our numerical results confirm that the basement‐hosting uranium‐rich brines are capable of forming high‐grade concentrated uranium deposits in the footwall region of the fault zone below the unconformity interface. In contrast, the sandstone‐hosting uranium‐rich brines lead to the formation of low‐grade sheet‐like uranium mineral precipitates above the sandstone layer rather than in the basement. Therefore, this study reveals that unconformity‐related uranium deposits tend to be formed when the basement unit serves as the main uranium source. Taking the sandstone layer as the main uranium source, no unconformity‐related uranium deposits can be formed since the resultant uranium precipitation has no spatial relation with the unconformity interface.

. In association with the formation of unconformity-related uranium deposits, previous numerical studies have largely concentrated on coupling fluid flow with heat transport and/or mechanical deformation (e.g., Cui et al., 2012aCui et al., , 2012bEldursi et al., 2020;Li et al., 2016Li et al., , 2017Pek & Malkovsky, 2016). Only very limited studies have coupled hydrothermal fluid flow with chemical reactions. For instance, Raffensperger and Garven (1995) simulated the formation of uranium deposits under equilibrium conditions with methane as a reducing agent. Aghbelagh and Yang (2014) set up more robust kinetic behavior for mineral dissolution and precipitation, and addressed the role of a faulted graphite zone. More recently, they examined the effect of fault dip angles and permeabilities on uranium mineralization (Aghbelagh & Yang, 2017). However, all the previous studies have assumed that uranium-rich brines are contained only in the sandstone layer, ignoring the fact that uranium-bearing brines are also present in the basement. Furthermore, they only emphasized the importance of fluid circulation within the sandstone layer, without addressing the role of fluid flow in the confining cover and basement unit. Therefore, this study aims to fill these knowledge gaps through evaluating the role of the location of uranium-rich brines in controlling uranium ore genesis.

Conceptual Model Development
As shown in Figure 1, our conceptual model is 5 km deep and 6 km wide, with the top boundary 3 km below the surface, which is based on some common features of typical deposits of this type in the Athabasca Basin. The solution domain is discretized uniformly by 160 cells vertically and horizontally, containing a confining cover, an intermediate sandstone layer, a basement unit, and an unconformity interface separating the sandstone and basement. Previous numerical studies have also employed similar layered models (e.g., Cui et al., 2012a;Li et al., 2016).
In this paper, it is assumed that aqueous methane is produced by dissolving graphite in fault zones at typical ore-forming brine temperatures, according to reaction 1: which is then used as a reducing agent to reduce uranium ore (Aghbelagh & Yang, 2014Raffensperger & Garven, 1995), following reaction 2: In the Athabasca Basin, the maximum total thickness of sedimentary rocks was 5-7 km, although due to erosion, the thickness is currently much less (Hoeve et al., 1981). Fluid inclusion analysis of the quartz indicates a homogenization temperature of 150-170°C and a salinity of 25 wt% NaCl equivalent for the basin brine (Jefferson et al., 2007). Reverse faults play a key role in the formation of unconformity-related uranium deposits, as they focused mineralizing fluids to deposition sites and provided reducing agents (i.e., CH 4 ) for uranium precipitation (e.g., Jefferson et al., 2007;Kotzer & Kyser, 1990;Kyser et al., 1989). The width of the fault zones varies from tens of meters to hundreds of meters (e.g., Bruneton, 1993;Li et al., 2016). They occur predominantly in the basement, with roots up to hundreds of meters deep, generally extending into the basinal sandstone for tens to hundreds of meters (e.g., Bishop et al., 2016;Derome et al., 2005;McGill et al., 1993). In addition, the dip of the fault zones also varies from very low to near vertical (e.g., Bishop et al., 2016;Finch, 1996;Hajnal et al., 2005;Hoeve & Sibbald, 1978). Based on these constraints on the fault zones in the Athabasca Basin, a faulted graphite zone is included in the model as illustrated in Figure 1. It slopes to the right at an angle of 59°, extends 625 m vertically (20 cells vertically deep), and is 128 m thick (4 cells horizontally wide). The fault zone is rooted mainly in the basement, with a vertical extension of 562.5 m below the unconformity and 62.5 m above. Similar fault zone dimensions have been used in previous studies (e.g., Eldursi et al., 2020;Li et al., 2016). Previous studies have reported that fault offsets of the unconformity interface are common in the Athabasca Basin, and the displacement ranges from tens of meters to hundreds of meters (e.g., Harvey & Bethune, 2007;Jefferson et al., 2007;McGill et al., 1993). Thus, in this study a displacement of 62.5 m (2 cells vertically thick) is included to the hanging wall of the fault zone to reflect the fault offset of the unconformity, as shown in Figure 1. To our best knowledge, no previous numerical study on unconformity-related uranium deposits has considered this parameter.
The confining cover consists of low-permeability shallow marine sedimentary rocks (mudstone, siltstone, and sandstone), mainly composed of calcite, kaolinite, quartz, dolomite, anhydrite, muscovite, and hematite. The intermediate sandstone layer is the main aquifer for fluid circulation, and its main minerals include quartz, muscovite, K-feldspar, anhydrite, hematite, and chlorite. The basement unit is almost impermeable and is composed of quartz, muscovite, K-feldspar, chlorite, anhydrite, hematite, and pyrite. The fault zone appears to act as a fluid conduit based on geological evidence that fault zones of this type were reactivated after the filling of the basin (Kotzer & Kyser, 1995) and remained conductive until recently (Hoeve & Quirt, 1984). Graphite is its main mineral together with muscovite, quartz, chlorite, pyrite, and kaolinite. Compiled from the parameters previously employed in related modeling studies and the physical properties obtained from in-lab measurements of the Athabasca group sandstone samples (e.g., Aghbelagh & Yang, 2014Cui et al., 2012a;King et al., 1988;Li et al., 2016;Raffensperger & Garven, 1995), Table 1 shows the density, porosity, permeability, and thermal conductivity of the four geological units, whereas Tables 2-5 tabulate the initial mineral volume fractions and aqueous component concentrations for each hydrostratigraphic unit. Table 6 shows the chemical parameters for the minerals involved, including the kinetic rate constant at 25°C K 25 (mol/m 2 s), activation energy E a (kJ/mol), and reactive surface area A m (cm 2 /g). Furthermore, the confining cover and the sandstone layer are assumed to be in oxidizing and acidic conditions, with log fO 2 = −14.8 and pH = 5.3, and log fO 2 = −22.8 and pH = 5.1, respectively, where fO 2 is oxygen fugacity. The basement unit and the fault zone are assumed to be in reducing and more acidic conditions, with log fO 2 = −46.8 and pH = 4.5, and log fO 2 = −51.3 and pH = 4.1, respectively.
As shown in Table 1, the sandstone layer is assumed to be highly permeable with a permeability of 3 × 10 −13 m 2 . This is justified by the following: (a) In-lab measurements of the physical properties of the Athabasca group sandstone samples (King et al., 1988). (b) The fluid inclusion study by Richard et al. (2016) indicates that the large-scale fluid circulation was across an entire 3-6-km-thick sedimentary pile and at least several hundreds of meters below the unconformity in the basement in the Athabasca Basin. (c) In previous ore-forming flow modeling studies, the permeabilities assigned for the sandstone unit at similar depth all have the same order of magnitude (e.g., Cui et al., 2012aCui et al., , 2012bOliver et al., 2006;Raffensperger & Garven, 1995). In particular, Li et al. (2018) and Li et al. (2021) assigned the same permeability of 3 × 10 −13 m 2 for the sandstone layer in the Athabasca Basin.
The top and bottom boundaries are set at fixed temperatures of 90 and 240°C, respectively, following a typical geothermal gradient of 30°C/km. The side and bottom boundaries are assumed to be impermeable to fluid flow, but the top boundary is set at a fixed fluid pressure of 30 MPa (e.g., Cui et al., 2012a). As for the boundary conditions of the chemical domain, the top and bottom are assumed to have fixed mineral volume fractions and aqueous component concentrations, equal to those of their respective units. For the side boundaries, the normal gradients of the volume fractions and concentrations are set to zero (e.g., Aghbelagh & Yang, 2017). The initial temperature distribution is calculated using the same geothermal gradient, and the initial fluid pressure is determined on a basis of hydrostatic conditions.

Numerical Modeling Method
Numerical computations were performed by using the commercial software package TOUGHREACT (Xu et al., 2011), which is capable of simulating chemically reactive nonisothermal flow of fluids in porous and fractured media, as well as interactions between mineral assemblages and fluids.
Governing equations include coupled nonisothermal fluid flow, solute transport, and reactive geochemistry. The flow and transport equations can be derived from the principle of mass and energy conservation. The chemical transport equations are written in terms of total dissolved concentrations of chemical components that are concentrations of their primary species plus their associated aqueous secondary species (Steefel & Lasaga, 1994;Walter et al., 1994). Advection and diffusion processes are considered for chemical transport, and diffusion coefficients are assumed to be the same for all aqueous species. The fluid density and viscosity are calculated in TOUGHRE-ACT through an Equation of State EOS1, which is particularly applicable to hydrothermal fluid flow problems. By choosing EOS1, all fluid properties (density, specific enthalpy, viscosity, and saturated vapor pressure) are calculated from the steam table equations as given by the International Formulation Committee (1967) (Pruess et al., 1999;Xu et al., 2004).
In this study, aqueous complexation, acid-base, redox, and cation exchange are assumed to be in local equilibrium. Mineral dissolution and precipitation are assumed to proceed under kinetic conditions, with the exception of anhydrite and calcite where an equilibrium approach is employed due to their fast reaction rate when reacting with aqueous species (Xu et al., 2011). The original chemical data in TOUGHREACT (Xu et al., 2011) did not include uranium aqueous species and minerals, so were modified by adding uranium aqueous species (e.g., UO 2

2+
, UO 2 Cl + , and UO 2 Cl 2 0 ) and uranium minerals (e.g., uraninite and rutherfordine). Refer to Yang (2014, 2017) for the details of the modified thermodynamic data. More details on the numerical modeling approach and the geochemical system can be found in previous publications by Xu et al. (2004Xu et al. ( , 2011 and Yang (2014, 2017 Note. Based on previous studies by Raffensperger and Garven (1995) and Aghbelagh and Yang (2017).

Table 2 Initial Mineral Volume Fractions and Aqueous Component Concentrations for the Confining Cover
To represent a geochemical system in TOUGHREACT, it is required to select a subset of aqueous species as primary species. All other species are called secondary species, including aqueous complexes and minerals. Secondary species are produced from primary species through chemical reactions (Xu et al., 2004). The aqueous uranium species in the Athabasca Basin are in the form of UO 2 2+ and uranyl chloride complexes (e.g., UO 2 Cl + , and UO 2 Cl 2 0 ) (e.g., Migdisov et al., 2018). Notice that UO 2 Cl + and UO 2 Cl 2 0 are produced from UO 2 2+ via the following two reactions: UO 2 2+ + Cl − = UO 2 Cl + and UO 2 2+ + 2Cl − = UO 2 Cl 2 0 (Migdisov et al., 2018). Therefore, in this study UO 2 2+ is selected as the primary species, and UO 2 Cl + and UO 2 Cl 2 0 are selected as the secondary species that are considered as important forms in transporting oxidized uranium in ore-forming fluids via the above-mentioned two reactions, refer to Appendix 2 in Aghbelagh and Yang (2014). Richard et al. (2012Richard et al. ( , 2016 performed fluid inclusion analysis of quartz veins in barren samples coexisting with several uranium deposits in the Athabasca Basin, indicating that the concentration of aqueous uranium UO 2 2+ in the basement brines ranges from 1.0 × 10 −6 to 2.8 × 10 −3 mol/L and the average is 1.0 × 10 −4 mol/L. A recent analysis of fluid inclusions in barren sandstones in the Athabasca Basin (Chi et al., 2019) indicates that UO 2 2+ concentrations range from 2.2 × 10 −6 to 9.9 × 10 −5 mol/L, with an average of 2.5 × 10 −5 mol/L. Based on these data sets, this study assumes that the UO 2 2+ concentrations in the confining cover and fault zone are 1.0 × 10 −6 and 1.6 × 10 −6 mol/L, respectively. For the sandstone layer and basement unit, we consider the following four scenarios in order to address whether the basement-hosting uranium-rich brines or the sandstone-hosting brines favor the ore genesis based on the observed uranium concentrations.

Results and Discussions
In Scenario 1, the sandstone layer accommodates the uranium-rich brine with a UO 2 2+ concentration of 1.0 × 10 −4 mol/L, while the basement UO 2 2+ concentration is 2.5 × 10 −5 mol/L. Scenario 2 is the opposite of Scenario 1, using the measured average UO 2 2+ concentrations for the basement unit and sandstone layer. That is, the basement accommodates a UO 2 2+ concentration of 1.0 × 10 −4 mol/L, while the sandstone has a UO 2 2+ concentration of 2.5 × 10 −5 mol/L. In Scenario 3, the basement unit has a UO 2 2+ concentration of 2.8 × 10 −3 mol/L (the measured maximum basal concentration), and the sandstone UO 2 2+ concentration is Note. Based on previous studies by Raffensperger and Garven (1995) and Aghbelagh and Yang (2017). 2.2 × 10 −6 mol/L (the measured minimum basinal concentration). In Scenario 4, the sandstone layer accommodates a UO 2 2+ concentration of 2.8 × 10 −3 mol/L, and the basement UO 2 2+ concentration is 1.0 × 10 −6 mol/L (the measured minimum basal concentration). All other conditions remain the same in these four scenarios. Figure 2a illustrates the fluid flow vectors at 300,000 years. The max flow rates in the confining cover, sandstone layer, basement unit, and fault zone are equal to 1.83 × 10 −3 , 1.50, 1.14 × 10 −3 , and 1.52 m/year, respectively. Figure 2b shows an enlarged view surrounding the fault zone, indicating that the left half of the fault zone carries a downflow and the right half carries an upflow. Flow vectors in the cover and basement unit are too small to be identified. To better visualize them, the original flow vectors in Figure 2a are sparsely sampled. We also used a constant vector length to allow direction of flow to be seen and used a color scale to represent differences in the magnitude of flow velocity, as shown in Figure 2c. It can be seen from Figures 2a and 2c that there are four convection cells developed in the sandstone layer, with two upwelling flow zones along the side boundaries flanked by two downwelling zones. The dominant upwelling zone is shifted to the left by about 500 m from the fault zone, which is due to the "downslope" flow over the fault offset of the unconformity that is directed toward the left, see Figure 2b. As such, the clockwise convection cell on the right is greater in size than the anticlockwise convection cell on the left.

Scenario 1
The basinal fluid travels upward through the three upwelling zones. When it reaches the top of the sandstone layer, most of it diverges to the left and right and then flows downward through the two downwelling zones. The rest of the basinal fluid penetrates across the sandstone-cover interface to mix with the cover fluid, which also moves upward but at a significantly lower rate (max 1.83 × 10 −3 m/year in the cover vs. 1.50 m/year in the sandstone). The cover fluid then diverges and flows downward to join the downwelling flow zones in the sandstone layer. On the other hand, as the basinal fluid descends through the downwelling zones to the bottom of the sandstone layer, most of it diverges and flows parallel to the unconformity interface and then merges with the upwelling zones. The rest seeps Note. Based on previous studies by Raffensperger and Garven (1995) and Aghbelagh and Yang (2017).

Table 4 Initial Mineral Volume Fractions and Aqueous Component Concentrations for the Basement Unit
across the sandstone-basement interface to mix with the basement brine, which also moves downward but at a much lower rate (max 1.14 × 10 −3 m/year in the basement vs. 1.50 m/year in the sandstone). The basement brine then diverges and flows upward to join the upwelling zones in the sandstone layer. It is also noticed from Figure 2c that almost all basement brine of the entire basement is concentrated in the footwall of the fault zone below the unconformity, except for some shallow brine near the side boundaries that discharges to the sandstone layer. Figure 2d illustrates the temperature distribution at 300,000 years, which has an asymmetric mushroom shape due to the aforementioned fluid circulation within the sandstone layer. The fluid flow in the cover and basement has a negligible effect on the temperature distribution, as evidenced by the nearly evenly spaced isotherms in these two units. However, both the cover flow and basement flow are important for uranium mineralization, although their rates are nearly three orders of magnitude lower than the rate of the basinal flow.
Based on the oxygen fugacity fO 2 assigned to each geologic unit in Section 2, two oxidation-reduction fronts occur. One is the sandstone-cover interface and the other is the sandstone-basement interface (i.e., the unconformity). At early time, uranium precipitation occurs along these two fronts when the basinal fluid percolates across the sandstone-cover interface and reacts with the more oxidizing cover fluid and when the basinal fluid infiltrates the basement and reacts with the reducing basement brine, respectively. Figure 3a shows the precipitated uraninite at 5,000 years, which is characterized by a horizontal sheet-like shape immediately above the sandstone-cover interface and below the unconformity with a max volume fraction of 4.42 × 10 −7 . Notice that the uraninite precipitation above the sandstone layer is stronger than that below the unconformity. The difference in amounts of uranium mineralization may be caused by the permeabilities assigned. As listed in Table 1 Note. Based on previous studies by Raffensperger and Garven (1995) and Aghbelagh and Yang (2017). Note. Based on previous studies by Raffensperger and Garven (1995) and Aghbelagh and Yang (2017). assigned a permeability that is slightly higher than that of the basement, which facilitates more basinal fluid percolation into the cover. Figures 3b, 3c and 3d illustrates the precipitated uraninite at 100,000 years, 200,000 years, and 300,000 years, respectively. As time progresses, more basinal uranium-rich fluid seeps into the overlying cover and the underlying basement, resulting in increased uranium deposition. The max volume fractions of the precipitated uraninite in Figures 3b, 3c and 3d are 3.73 × 10 −6 , 3.75 × 10 −6 , and 4.05 × 10 −6 , respectively. Starting from 100,000 years, the precipitated uraninite in the confining cover is mainly distributed just above the sandstone-cover interface between the upwelling and downwelling zones, and it is also distributed along narrow bands in the downwelling zones close to the side boundaries, representing the converging belts of the percolated basinal fluid with the more oxidizing cover fluid, see Figure 2c. As mentioned above, most of the basement brine is focused in the footwall region of the fault zone below the unconformity. As a result, starting from 100,000 years the uranium deposition in the basement is concentrated in the footwall. Again, the precipitated uraninite in the basement grows over time, see Figures 3a-3d. However, because the uranium-rich brine is hosted in the sandstone layer, the volume fraction and size of the precipitated uraninite in the basement are always smaller than in the cover. Figure 4a shows the precipitated uraninite at 5,000 years, indicating that a wide range of uraninite precipitation occurs immediately below the unconformity with a max volume fraction of 2.52 × 10 −7 . Uraninite mineralization also occurs above the sandstone-cover interface but at a much lower volume fraction because the uranium-rich brine is now hosted in the basement. The depletion of UO 2 2+ from uraninite precipitation results in a decrease in its concentration in the basement near the unconformity. Thus, at late time the basement brine becomes more important as it drives deep aqueous uranium to the shallow region. As shown in Figure 2c, the basement brine is focused dominantly in the footwall region of the fault zone, allowing deep aqueous uranium to rise and react with shallow oxidizing fluid that seeps into the basement from the overlying sandstone layer. Figure 4b illustrates the precipitated uraninite at 100,000 years. Now uranium deposition is concentrated in the footwall with a max volume fraction of 3.17 × 10 −6 . Thus, the footwall region seems to represent a favorable structural trap for uranium deposition. As time increases, more of the basement brine is concentrated into the structural trap, bringing more aqueous uranium from the depths. This increases the size and volume fraction of the precipitated uraninite. Figure 4c shows the precipitated uraninite at 200,000 years with a max volume fraction of 8.99 × 10 −6 , and Figure 4d illustrates the precipitated uraninite at 300,000 years with a max fraction of 1.08 × 10 −5 . The uranium deposit is located in the same structural trap but the size and volume fractions increase slightly from 200,000 to 300,000 years.

Scenario 2
Comparing Scenario 2 with Scenario 1 reveals that the basement-bearing uranium-rich brine enables the formation of more concentrated, higher grade uranium deposits in the footwall region. In contrast, the sandstone-bearing uranium-rich brine leads to sheet-like uranium mineralization above the sandstone-cover interface but at a much lower grade.

Scenario 3
Figures 5a and 5b illustrates the precipitated uraninite at 200,000 and 300,000 years, with max volume fractions of 2.75 × 10 −4 and 3.06 × 10 −4 , respectively. As in Scenario 2, the uranium deposit is again formed in the same structural trap but the volume fraction of precipitated uraninite is increased by about 30 times, as evidenced by comparing Figures 5a and 5b with Figures 4c and 4d. This is because the UO 2 2+ concentration in the basement fluid is now 28 times higher than in Scenario 2. The precipitated uraninite in the confining cover has a negligible volume fraction in comparison with that in the basement unit because the sandstone layer has a lower UO 2 2+ concentration of 2.2 × 10 −6 mol/L in this scenario compared with Scenario 2.

Scenario 4
Figures 6a and 6b illustrates the precipitated uraninite at 200,000 and 300,000 years, respectively, with a max volume fraction of 1.03 × 10 −4 . As in Scenario 1, the predominant uraninite precipitation occurs immediately above the sandstone-cover interface and along the converging belts in the downwelling zones in the cover. The max volume fraction is increased by about 30 times, and this is because the UO 2 2+ concentration in the sandstone layer is 28 times higher than in Scenario 1. The volume fraction of the precipitated uraninite in the footwall region is reduced substantially because the basement brine is now assigned a lower UO 2 2+ concentration of 1.0 × 10 −6 mol/L in comparison with 2.5 × 10 −5 mol/L in Scenario 1.

Conclusions
Reactive mass transport modeling has been conducted to evaluate the role of the source location of uranium-bearing brines in controlling uranium ore genesis. Our results from four numerical case studies show that the fluid circulation in each geological unit is important for temperature distribution, transport of aqueous components, and uranium deposition. Flow rates in the confining cover and the basement unit are about three orders of magnitude lower than those in the sandstone layer and fault zone. However, fluid flow in the cover and basement is critical to the genesis of uranium deposits. At early time, uraninite precipitation initially occurs immediately above the sandstone-cover interface and below the unconformity, but at a very low volume fraction. At later time, different locations of uranium-rich brines result in different distributions and grades of precipitated uraninite. When the sandstone layer accommodates the uranium-rich brine, the predominant uranium precipitation occurs in the confining cover rather than in the basement unit since it is relatively easier for the basinal fluid to percolate across the sandstone-cover interface than across the sandstone-basement interface. As a result, the precipitated uraninite takes on a sheet-like shape directly above the sandstone layer and between the upwelling and downwelling zones in the cover. When the basement unit holds the uranium-rich brine, the predominant uranium precipitation occurs in the footwall of the fault zone below the unconformity because the basement reducing brine is mainly focused in such a favorite trap to react with the basinal oxidizing fluid being percolated from the overlying sandstone layer. As a result, uraninite deposition is concentrated in the footwall region with a much higher volume fraction. Therefore, our numerical results reveal that the scenarios with the basement unit serving as the main uranium source favor the formation of unconformity-related uranium deposits, while the scenarios

Data Availability Statement
Reactive flow modeling was conducted by using the commercial software package TOUGHREACT. Details of the software and modeling methodology can be found in published journal articles by Xu et al. (2011) and Yang (2014, 2017). All data used in this research are tabulated in Tables 1-6. No other data were used nor created for this research. The data tables are also available at Xu et al. (2022).