Edinburgh Research Explorer The Effect of Grain Size on Porewater Radiolysis Citation for published version: DeWitt, J, McMahon, S & Parnell, J 2022, 'The Effect of Grain Size on Porewater Radiolysis', Earth and Space Science, vol. 9, no. 6, e2021EA002024, pp. 1-21. https://doi.org/10.1029/2021EA002024 Digital Object Identifier (DOI): 10.1029/2021EA002024 Link: Link to publication record in Edinburgh Research Explorer Document Version: Peer reviewed version Published In: Earth and Space Science General rights Copyright for the publications made accessible via the Edinburgh Research Explorer is retained by the author(s) and / or other copyright owners and it is a condition of accessing these publications that users recognise and abide by the legal requirements associated with these rights. Take down policy The University of Edinburgh has made every reasonable effort to ensure that Edinburgh Research Explorer content complies with UK legislation. If you believe that the public display of this file breaches copyright please contact openaccess@ed.ac.uk providing details, and we will remove access to the work immediately and investigate your claim. Download date: 17. Nov. 2022 manuscript accepted at Earth and Space Science 1 The Effect of Grain Size on Porewater Radiolysis J. DeWitt1, S. McMahon2,3, and J. Parnell4 1 Department of Physics, East Carolina University, Howell Science Complex, 10th Street, Greenville NC 27858, US 2 School of Physics and Astronomy, University of Edinburgh, James Clerk Maxwell Building, Peter Guthrie Tait Road, Edinburgh EH9 3FD, UK 3 School of Geosciences, University of Edinburgh, Grant Institute, Edinburgh EH9 3FE, UK 4 School of Geosciences, University of Aberdeen, Meston Building, King’s College, Aberdeen AB24 3UE Corresponding authors: Joel DeWitt (dewittjo@ecu.edu), Sean McMahon (sean.mcmahon@ed.ac.uk) Key Points: Porewater radiolysis in rocks and minerals yields molecular hydrogen (H2), an energy source for microbial life on Earth and perhaps Mars. Our Monte Carlo simulations of mineral radiation reveal how grain size controls H2 yield for a given radioisotope composition and porosity. For some realistic scenarios, clay-sized grains can produce up to an order of magnitude more H2 per unit time than sand-sized grains, other things being equal. manuscript accepted at Earth and Space Science 2 Abstract The radiolysis of porewaters by uranium, thorium, and potassium in mineral grains is a recognised source of molecular hydrogen in rock- and sediment-hosted fluids. This radiolytic hydrogen is of geomicrobiological interest as a potential energy source (electron donor) for microbial metabolism, especially in energy-limited settings such as the marine deep biosphere or the subsurface of Mars. Previous efforts to predict the production of radiolytic hydrogen from columns of rock and sediment have tended to rely upon analytic models that cannot account for the attenuation of mineral radiation by grains larger than ~30 microns. To address this, we have developed a Monte Carlo method to simulate the physics of mineral radiation and evaluate the production of H2 as a function of mineral grain size and radioisotope composition. The results confirm that grain size is a major control on radiolytic H2 yield. For example, using the standard geological classification of grain sizes, we find that clay can produce up to an order of magnitude more H2 per unit time than sand. The magnitude of this effect is illustrated using compositional data from real geological units in order to demonstrate the dependence of radiolytic hydrogen flux on natural radionuclide concentration and bulk porosity. 1. Introduction The oxidation of molecular hydrogen (hydrogenotrophy) provides energy for microbial growth in a wide variety of environments and may have sustained Earth’s earliest metabolising organisms (e.g., Kral et al., 1998; Schulte et al., 2006; Sleep et al., 2011; Lane and Martin, 2012; Preiner et al., 2020). In some subsurface habitats, naturally occurring H2 serves as an electron donor for methanogens, acetogens, sulfate-reducers and other anaerobic microorganisms forming part of the “deep biosphere”. This includes volcanic and siliciclastic aquifers as well as marine sediments low in organic matter (Stevens and McKinley, 1995; Chapelle et al., 2002; Lin et al., 2005; Blair et al., 2007). Somewhat analogous habitats could potentially support similar metabolisms in the deep subsurface of Mars and other planets and moons (Stevens and McKinley, 1995; Chapelle et al., 2002; Blair et al., 2007; Nixon et al., 2013). Although there is no direct evidence for life on these worlds as yet, methane has been detected repeatedly in the Martian atmosphere (most recently by the Curiosity rover and Mars Express Orbiter; Webster et al., 2015; Giuranna et al., 2019), and molecular hydrogen is also present in aqueous plumes venting from the subsurface ocean of Saturn’s moon Enceladus (Waite et al., 2017). Hypothetically, life on these worlds might be supplied with molecular hydrogen from several natural sources also known to be active on the Earth, including organic fermentation, reactions between water and silicate minerals (known as “serpentinization”; McCollom, 1999), diffusion from crystallographic defects (Freund et al., 2002), and the reduction of water during “mechanoradical” processes driven by frictional grinding on fault planes (Hirose et al., 2011; McMahon et al., 2016). Independently of these sources, H2 is produced by extremely rapid interactions between protons, electrons and radical species released from water molecules under alpha, beta and gamma irradiation issuing from radioactive elements in minerals (235U, 238U, 232Th, 40K and their respective short-lived radiogenic daughters). Radiolysis is independent of temperature, and can generate H2 from both porewater (including ice) and crystallographic water in hydrated minerals (Blair et al., 2007; Smetannikov, 2011). Water radiolysis also produces reactive oxygen species capable of oxidizing minerals to form electron acceptors such as SO42-; manuscript accepted at Earth and Space Science 3 radiolysis can thus support microbial life at both ends of the electron transport chain (Lefticariu et al., 2010). Several studies have attempted to model the radiolytic hydrogen (H2) flux from rocks and sediments beneath the seafloor and land surface under different conditions, including conditions that represent different planetary bodies. Lin et al. (2005) and Blair et al. (2007) considered crystalline continental rocks and marine shelf sediments, respectively, focusing on water- saturated pore spaces between grains, rather than fractures. Both studies used an intuitive analytical model derived from Aitken (1985) via Hofmann (1992). This model combines measured U, Th and K abundances, their activities, and reference data for the expected H2 yield from water per energy absorbed (Spinks and Woods, 1990). It accounts for the attenuation of radiation energy within the bulk medium by factoring in the radiation stopping power of water, the stopping power of the mineral grains (assumed to equal that of silica, SiO2, which has been experimentally determined for alpha, beta and gamma radiation), and the ratio of water to mineral matter by volume (i.e., porosity). Using this model, Lin et al. showed that radiolysis is the most likely cause of high H2 concentrations measured in groundwaters of the Witwatersrand Basin in South Africa. They also modelled radiogenic 4He production and compared the predicted H2:He ratio with the empirically measured dissolved gas concentrations in Witwatersrand groundwater; the difference showed that H2 was either produced less efficiently than the model predicted or was being consumed by a combination of abiotic reactions and hydrogenotrophic microorganisms. Onstott et al. (2006) adapted the model for Mars by substituting typical Martian values for crustal porosity and the concentration of parent nuclides. Assuming fully water-saturated basement, H2 production rates were found to be comparable to crystalline basement on Earth because of the higher porosity of Mars’ crust, despite the lower concentrations of uranium and thorium. Updated (higher) radionuclide concentrations were used by Tarnas et al. (2018, 2021) to estimate an even greater rate of radiolytic H2 production on Mars, using the same model. Bouquet et al. (2017) applied a similar model to Jupiter’s moon, Europa. Dzaugis et al. (2016) introduced a model of radiolytic H2 production from water hosted within planar fractures in radioactive rock, rather than intergranular pore spaces, and applied it to subseafloor basaltic basement. Dzaugis et al. (2018) showed that this fracture-based model yields radiolytic H2 production rates for Martian basaltic crust broadly comparable to those predicted for fine-grained Martian sediments by the model of Blair et al. (2017). One implication of these studies, in agreement with other considerations, is that the best place to seek extant life on Mars today is probably the subsurface (e.g., Tarnas et al., 2021). Hitherto, studies of porewater radiolysis have not properly accounted for the attenuation of each radiogenic particle or ray within its mineral grain (clast or crystal) of origin. This attenuation can be safely ignored when the grain diameter is much smaller than the stopping distance of the relevant particle in the mineral medium, as for example in clay. However, alpha particles in the kinetic energy range characteristic of mineral radiation (~5.6 MeV on average) typically have stopping distances in SiO2 (and similar materials) of ~23 µm (National Institute of Standards and Technology, 2021a). Consequently, if such an alpha particle is produced near the edge of a grain of diameter ?20 µm (i.e., coarser than a medium silt) and travels towards the interior of the grain, it will not have enough kinetic energy to escape from the other side, resulting in a null contribution to the radiolysis of the surrounding porewater. Similar considerations apply to beta particles—although their stopping distances are substantially larger—given a characteristic mean kinetic energy of ~0.33 MeV, equating to a stopping distance of ~452 µm in SiO2 (National Institute of Standards and Technology, 2021b). Furthermore, the manuscript accepted at Earth and Space Science 4 radiolysis contribution due to gamma (photon) radiation requires special treatment given their chargeless, massless character. Here, we quantitatively account for the variability of grain size in order to expand the range of environments in which radiolytic H2 production and its potential contribution to microbial productivity can be estimated. The method detailed in this paper is generalizable to the extent that it can be adapted for a variety of geological contexts and sediment or rock types with any level of radioactivity, whose constituent grains can be of arbitrary size and shape. In this particular study, we seek to model radiolytic H2 production from porous rocks and sediments with grain sizes ?20 µm for the first time (to our knowledge). One motivation for doing so is the fact that high concentrations of radioactive minerals on Earth commonly occur in relatively coarse sedimentary rocks, e.g., sandstone-type uranium deposits. In addition, relatively coarse granular materials can accommodate more microorganisms, which may tend to be excluded by the restrictive pore throat sizes of finer-grained materials, even at the same porosity (e.g., Krumholz et al., 1997). 2. Materials and Methods We have written a Fortran program that numerically simulates the physics of natural radioactive decay within a matrix of water and mineral grains using a Monte Carlo (MC) method. (The algorithmic foundation of this program is given in Supplementary Information S1.) MC models of alpha attenuation in mineral grains have previously been deployed in some geochronological studies (e.g., Farley et al., 1996; Hourigan et al., 2005). Although analytical models are sometimes easier to implement than numerical models, an MC approach offers maximum flexibility since it can ultimately be adapted for grains of any size and shape, with any distribution of radionuclides, any arrangement of pore space, etc. (Ketcham et al., 2011). In the present study, we aim to address the effect of grain size independent of other factors. For simplicity, grains are therefore represented as perfect spheres composed of pure SiO2. Real mineral grains typically deviate from a spherical geometry (especially clay minerals) but modelling the geometry in this way leverages the rotational symmetry of the system, which eases the conceptual and computational aspects of the MC method. The primary goal of our simulations is to score the total amount of kinetic energy transferred to the surrounding water layer due to the decay of radioactive isotopes within the SiO2 grains. Figure 1 shows the geometry of a spherical SiO2 grain and a composite of sample thorium decay events occurring within it. manuscript accepted at Earth and Space Science 5 Figure 1. A randomly selected SiO2 grain of radius 0.059 cm (590 ?m) showing a sample thorium decay against , and coordinate axes (black arrows with gradations). Alpha, beta, and gamma decay events register as red, green, and blue, respectively. These depth vectors locate exit points that correspond to solutions of the sphere . The program leverages the MC method and the physics of radioactive decay to drive a given simulation towards a well-averaged endpoint. The steps the simulation uses are as follows: 1. Choose a random SiO2 grain whose radius is within the desired interval of clay, silt, or sand (Wentworth, 1922), or a logarithmic distribution of all three; 2. Choose a random K, Th, or U constituent (weighted by the proportion of that element in the composition; see “Activity of Bulk Rock or Sediment” below); 3. Choose a random point ( , , ) within the grain at which the decay of the constituent plus its daughters will proceed; 4. Examine the entire decay chain in order, choosing random exit vectors, and applying the appropriate physics of alpha, beta, and gamma propagation in SiO2; 5. Repeat Steps 1 – 4 many times (typically 105 trials) in order to ensure the simulation converges while keeping a running total of the amount of kinetic energy (in MeV) deposited in the surrounding water layer; 6. Output radiolytic H2 production from the volume of rock or sediment; 7. Using Step 6 as a baseline, scale H2 production by volume-corrected porosity (see “Porosity” below). Athy’s law can be used to scale radiolytic H2 production in the case where porosity can be estimated for a column of rock or sediment as a function of depth. Point 2 assumes that all radioactive nuclides are uniformly distributed throughout each grain, that all grains are the same in U, Th and K content. This is a useful simplification for examining the effect of grain size as an independent variable, but in real rocks radionuclides may be heterogeneously distributed within grains and most concentrated in particular minerals (e.g., manuscript accepted at Earth and Space Science 6 zircon) that have a different grain size from the rest of the sediment or rock. For example, radioactivity is commonly associated with secondary diagenetic minerals such as iron oxides and clays that adsorb U and Th (Murray et al., 2015). Point 6 assumes that all radiolytic H2 is generated from water and not from the radiolysis of hydrated or hydroxylated minerals, a less well understood phenomenon (Smetannikov, 2011; Westbrook et al. 2015). These considerations are not taken further in the present study but suggest possible avenues for future study. The “length” of the simulation is controlled by the number of K, Th, and U decays + daughters that are tracked (read: trials), with the understanding that a large number of trials will enable the simulation to converge in a manner dictated by the physics of radioactive decay. For the purposes of this study, our simulations output 1) H2 production in units of moles per [Earth] year-cm3 as a function of sediment depth in meters, or 2) H2 production as a function of grain size (logarithmically distributed) in cm. 2.1 Radioactive Decay This section describes how the components of natural radioactivity interact with the medium through which they are propagating. This is also referred to as stopping power or linear energy transfer (LET). Our simulations treat alpha, beta, and gamma decay piecewise, wherein each decay daughter in a given chain is assigned a random exit vector. The entire decay chains for 238U, 235U, 232Th, and 40K were hard-coded into the program using data from the Interactive Chart of Nuclides provided by Brookhaven National Laboratory (2021). (The average energies required for each transition within a given decay chain were implemented rather than the “endpoint” or maximum energies, as previous studies have done.) These vectors are then used to propagate each decay product through a geometrically-computed depth of SiO2 material before depositing the remainder of its kinetic energy in the surrounding water layer. This strategy is advantageous since it avoids the assumption of alpha, beta, and gamma dominance through the use of lump energy sums (Blair et al., 2007). Overall, alpha, beta, and gamma emission resulting from natural radioactive decay will occur with average kinetic energies of 5.63 ± 1.68, 0.334 ± 0.267, and 0.227 ± 0.639 MeV, respectively. 2.1.1 Alpha Particles The interaction of a heavy charged projectile with the electrons in a stopping medium (absorber) is described by the Bethe-Bloch formula (Bethe, 1930). The amount of kinetic energy imparted to the target medium per unit path length is given as (1) , where is the coupling strength of electromagnetism, is the atomic number of the charged particle or ion projectile, is the rest mass-energy of the electron, is the electron density of the absorber with effective atomic number and effective atomic mass manuscript accepted at Earth and Space Science 7 (when referring to compounds and mixtures), and is the mean ionization potential of the target medium. The dimensionless number from special relativity. The last term in the square brackets is the density effect correction that accounts for the phenomenon that, in a dense medium, the field which perturbs electrons far from the projectile track is modified by the dielectric polarization of the atoms between the distant electrons and the projectile (Fermi, 1940). The quantity in Equation 1 is sometimes preceded by a negative sign to indicate that kinetic energy is transferred to the stopping medium. Alpha particles (helium nuclei, = 2) that are emitted during the course of natural radioactive decay possess an average kinetic energy that is relatively low, or ~5.6 MeV. This low kinetic energy necessitates a correction to the stopping power due to the phenomenon of electron capture, which accounts for the fact that the bare nuclear charge of the projectile is reduced and can thus be replaced by an effective projectile charge (Weaver and Westphal, 2002). Once corrected, Equation 1 is used to numerically propagate alpha particles through a known amount of SiO2 with characteristics and . From a radiolysis standpoint, the production of molecular hydrogen per unit path length is greatest at the end of an alpha particle’s range, a region commonly referred to as the Bragg peak in radiation physics. A charged particle or ion projectile with atomic number of rest mass significantly larger than that of the electron (charge ) is considered a “heavy” charged particle. This would include mesons, protons, alpha particles, and higher nuclei. Electrons and positrons are considered “light” charged particles; the physics of electron propagation through matter is discussed in the next section. 2.1.2 Beta Particles While the Bethe-Bloch formula can be used in principle to calculate the stopping power of any charged particle passing through matter, it must be modified in the case of beta particles (electrons). There are several reasons for this: one must always consider special relativity when handling electron-electron interactions; an electron that acts as a projectile can lose nearly all of its kinetic energy in a single interaction with a target electron; quantum mechanics states that we will not be able to distinguish between the projectile and the target electrons after the collision (they are identical particles); for a high energy electron that traverses a high stopping medium, radiative processes such as bremsstrahlung (from the German: “breaking radiation”) will be significant. The total stopping power of an electron as a function of kinetic energy in a specified target material is the sum of the collision plus bremsstrahlung contributions (Bethe and Heitler, 1934): (2) , manuscript accepted at Earth and Space Science 8 where is the stopping power due to electron collisions, , (3) and is the stopping power due to electron bremsstrahlung, (4) , where is the fine-structure constant. Equations 2-4 are used to numerically propagate beta particles through a known amount of SiO2 with characteristics similar to those described by Equation 1. The simulation of radioactive decay in a mineral medium requires that Equations 1-4 be applied to the propagation of alpha and beta particles inside spherical SiO2 grains. As a means of developing and validating our numerical method, it would be well if the algorithm given in Supplementary Information S1 were checked against a pure geometric argument. It has been found analytically that if particles with stopping distance s are emitted at random vectors from points inside a spherical grain of diameter , the proportion (or percent fraction) of these particles that escape from the grain rather than terminating (or stopping) within it is given by . (5) Equation 5 has been obtained previously (e.g., Farley et al., 1996); a derivation is provided in Supplementary Information S2. Figure 2 shows a graphical comparison of these two methods in the case of alpha particles, along with an extension to beta particles using the same algorithm given in Supplementary Information S1. manuscript accepted at Earth and Space Science 9 Figure 2. Escape efficiency of alpha particles with stopping distance = 19.6 ?m as a function of spherical grain radius, where the methods of MC simulation and geometric analysis are compared. The escape efficiency of beta particles is also shown using the same algorithm. The comparison shown in Figure 2 corresponds to an alpha particle test case with stopping distance = 19.6 ?m (kinetic energy: 4.93 MeV). These curves disagree on average by 0.2% over 100,000 trials per micron radius. The good agreement between the core algorithm of the simulation and our geometric analysis partially validates our MC approach in this study. 2.1.3 Gamma Rays Unlike energetic charged particles which have a finite range when traversing matter, gamma rays (photons) of a given energy do not possess such a discrete range, but rather are absorbed exponentially as a function of distance in the absorbing medium. This phenomenon is referred to as attenuation. (Alpha and beta particles, strictly speaking, do not attenuate; we use this term contextually to service the reader.) For a beam of monoenergetic photons, the rate of absorption of the photons as a function of depth in the absorber is manuscript accepted at Earth and Space Science 10 , (6) where is the absorption cross section of each atomic scattering center (or absorbing molecule), and is the volume density of atomic scattering centers. (The negative sign above represents absorption.) We define , where is referred to as the linear attenuation coefficient. We can then integrate Equation 6 with respect to distance to obtain . (7) Similarly, if we consider depth in terms of not distance, but of areal density (in g/cm2), we can define the mass attenuation coefficient such that (8) , where is the mass density (in g/cm3), is the atomic molar mass (in g/mol) of the absorbing medium, and is Avogadro's number (in molecules/mol), respectively. Equation 7 now becomes , (9) where the values of are tabulated. Using Equation 9, a characteristic attenuation (or absorption) length can be formulated as an analog to the stopping distance of alpha and beta particles; for a single photon, this is taken as the distance into a material when the probability of absorption has dropped to , or ~37%. Explicitly, this gives a propagation distance that is equal to the inverse of the mass attenuation coefficient, or . A typical attenuation length in this study is on the order of centimeters. There are a number of different mechanisms by which photons can interact with matter. However, only three of these mechanisms—the photoelectric effect, Compton scattering, and electron/positron pair production—dominate, while the remainder usually make only a negligible contribution. The total absorption cross section in Equation 8 is the sum of the cross sections for each type of photon interaction: , (10) manuscript accepted at Earth and Space Science 11 where is the cross section for photoelectric interactions, is the cross section for Compton interactions, and is the cross section for pair production. (A cross section in nuclear physics can informally be read as an energy-dependent probability for interaction.) Since from Equation 8, there is a corresponding mass attenuation coefficient for each cross section. 2.1.3.1 Photoelectric Effect Free electrons cannot fully absorb photons without violating conservation of energy and conservation of momentum, and for that reason, free electrons are not subject to the photoelectric effect. Indeed, the photoelectric effect should be thought of as interaction between an incident photon and an entire atom, not just one of its electrons (Einstein, 1905). Photoelectric interactions are more likely to eject electrons from the more tightly bound inner electron shells of the atom than the less tightly bound outer shells, since interactions involving inner shell electrons occur closer to the center of the atom where it is easier to transfer an incident photon’s linear momentum to the atom. This is because nearly all the mass of the atom is concentrated in its nucleus. The cross section for the photoelectric effect is: (11) , where is the fine-structure constant, is the atomic number of the absorbing medium, is the classical electron radius, the quantity is the rest mass-energy of the electron, and finally is the energy of the incident photon. Equation 11 should be applied if the energy of the incident photon is between 0.1 and 0.35 MeV, or where the photon energy is above the K-absorption edge. 2.1.3.2 Compton Scattering The interaction of photons with “nearly” free electrons is described by a number of mechanisms: Photon interactions with “truly” free electrons is referred to as Thomson scattering; Coherent (or elastic) photon interactions with “nearly” free electrons is referred to as Rayleigh scattering; Incoherent (or inelastic) photon interactions with “nearly” free electrons is referred to as Compton scattering. Conservation of energy for the Compton scattering interaction gives manuscript accepted at Earth and Space Science 12 (12) , where from special relativity, with and the energies of the incident and scattered photon, respectively. After obtaining the differential cross section for Compton scattering—which is specified by the Klein-Nishina formula (Klein and Nishina, 1929)—it can be shown that integration over all angles gives the cross section for the Compton interaction: (13) , where effectively compares the incident photon energy to the rest mass-energy of the electron, or . 2.1.3.3 Pair Production Pair production is essentially the reverse process of bremsstrahlung. When in close proximity to a heavy nucleus, a photon can be converted into an electron/positron pair. Conservation of energy for this process dictates , (14) where and are the energy of the electron and positron, respectively. Because the rest mass of the electron and the positron must come from the energy of the photon, there is a threshold energy of = 1.022 MeV below which this process cannot occur. The cross section for pair production originates from quantum electrodynamics, and is given by (Gould and Schréder, 1967): (15) , where from special relativity. In this case is actually the relativistic energy of the electron/positron in the center-of-mass frame of reference. This cross section is applied for the transitions 40K ? 40Ar in 40K (1.461 MeV), 208Tl ? 208Pb in 232Th (3.376 MeV), and 214Bi ? 214Po in 238U (1.476 MeV). manuscript accepted at Earth and Space Science 13 2.2 Activity of Bulk Rock or Sediment The activity per unit volume of bulk rock/sediment in units of Bq/cm3 is given by (16) , where is the mass density of the rock/sediment constituents in g/cm3, is the Avogadro constant, is the molar mass of the radioisotope, N?g/g is the concentration of the radioisotope in micrograms per gram (?g/g), and is the half-life of the radioisotope in seconds. Equation 16 is used to convert the concentration of a given radioisotope component (uranium, thorium, or potassium) to the volume-normalized activity in bulk sediment/rock. Table 1. A summary of radioisotopes relevant to natural radioactive decay in minerals. Radioactive Isotope Molar Mass, (g/mol) Half-life, (s) 238U (99.274%) 238.0508 1.4090E17 235U (0.720%) 235.0439 2.2195E16 232Th 232.0381 4.4308E17 40K (0.0117% of 40K in 39K) 39.0983 3.9452E16 2.3 Radiolysis Hydrogen production is evaluated according to the water radiolysis model of Blair and coworkers (Blair et al., 2007): , (17) where , , and are the H2 yields per unit kinetic energy produced in pure water by radiolysis, and , , and are the kinetic energies absorbed by water in bulk rock/sediment due to alpha, beta, and gamma decay, respectively. The values of , , and are 1.57E4, 5.30E3, and 4.50E3 molecules/MeV, respectively (Spinks and Woods, 1990; Blair et al., 2007); the values of , , and (in MeV) are determined by the simulation. These G values are affected very little by changes in temperature, pressure, and pH (Draganic and Draganic, 1971). They are, however, enhanced to some extent by salinity due to the reaction of inorganic solutes with hydroxide (Spinks and Woods, 1990; Lin et al., 2005). manuscript accepted at Earth and Space Science 14 2.4 Porosity A porosity correction is needed to account for the fact that only the kinetic energy deposited in the intergranular water (rather than other grains) will contribute to H2 production. The porosity of a compressible rock or sediment column decreases with increasing depth as the grains are compacted together. The simulations presented in this paper assume that all the porosity is occupied by water; we use Athy’s law to quantify porosity as a function of sediment (or rock) depth (Fowler and Yang, 1998): , (18) where is the porosity near the surface and is the compaction coefficient and has units of m–1. In sandstones on Earth, takes an average value close to 5.3E-4 m–1 (McMahon and Parnell, 2014); for geological features on Mars, this study uses a slightly lower compaction coefficient of 3.5E-4 m-1 (Clifford et al., 2010). Equation 18 represents a steady-state leading-order solution to equilibrium compaction. A volume correction factor is required in order to properly scale H2 production as a function of depth. This factor, 1 - ?z, accounts for the fact that in a given volume of rock/sediment, a large porosity will translate to a smaller proportion being occupied by radioactive mineral grains rather than water. The kinetic energy yield from a given volume of wet rock/sediment is equal to the energy yield of each grain (supplied by the simulation), multiplied by the total number of grains contained in the volume. This total grain number is equal to the volume of the rock/sediment multiplied by a factor of 1 - ?z, normalized by the volume of each grain (Blair et al., 2007). In other words, one must multiply Equation 17 by Equation 18 (Athy’s law), corrected by 1 - ?z = 1 - ?!? "#$. The reader should note that the complete porosity correction used to scale H2 production as a function of depth z is equal to ?z(1 - ?z), which contains a global maximum at ?z = 0.5. manuscript accepted at Earth and Space Science 15 Figure 3. A comparison of the predicted radiolytic hydrogen flux using the MC method (orange line) with the analytic model of Blair et al. (2007) using the radionuclide concentrations, grain sizes, and porosities reported for a marine sediment column sampled by the Ocean Drilling Program in the Peru Basin. The left panel shows the anticipated H2 production after 200 trials, while the right panel shows the same after 200,000 trials. The differences between the predictions of Blair et al. and this study are methodological and do not relate to grain size effects. Figure 3 compares the predicted H2 production given by the MC method with those calculated by Blair et al. (2007) from marine sediment column measurements as a function of depth. This figure demonstrates the compatibility of our simulations with the deep marine sediment column in which empirically-determined radionuclide concentrations, grain sizes, and porosities (measured by the Ocean Drilling Program in the Peru Basin) were provided by Blair and coworkers. Figure 3, in conjunction with Figure 2, further validates our MC methodology. The left panel in Figure 3 displays a comparison after 200 MC trials, while the right panel shows the same after 200,000 trials. Our simulations scale radiolytic hydrogen production using the measured porosity ?z and implemented with ?$(1 ? ?$) described above. In general, this figure shows that the MC method accurately reproduces the calculations of Blair et al. where confidence is high (i.e., standard error is low); any disagreement above ~55 meters is be manuscript accepted at Earth and Space Science 16 attributed to 1) the error in measured 232Th concentrations, and 2) the differing decay energy sums of the two approaches. (Whereas Blair et al. used the endpoint, or maximum, decay energies within a given chain, this study used instead the average to better represent the spectrum of possible outcomes.) Finally, it should be noted that the grain sizes, with an average diameter of 6.2 ?m over 120 meters, are too small in this case to attenuate the modeled flux. 3. Results 3.1 Effects of Grain Size Figure 4. Radiolytic hydrogen production, normalized to maximum in clay, as a function of sediment (or rock) depth for the test case of 7.4 ?g/g U, 2.8 ?g/g Th, and 1.3% K, with = 0.5 and = 5.3E-4 m–1. The SiO2 grain size is sampled over clay (0.06 - 3.9 ?m), silt (3.9 - 63 ?m), and sand (63 - 2,000 ?m). (In this simple model, porosity is independent of grain size.) Figure 4 shows modelled H2 production, normalized to maximum in clay, as a function of depth for the test case of 7.4 ?g/g U, 2.8 ?g/g Th, and 1.3% K (arbitrary but fairly typical values for geological materials) with = 0.5 and = 5.3E-4 m–1; grain size is randomly sampled over intervals corresponding to the grain sizes of clay (0.06 - 3.9 ?m), silt (3.9 - 63 ?m), and sand (63 - 2,000 ?m). Radiolytic hydrogen production in this figure is scaled by the factor ?$(1 ? ?$) described in Section 2.4, which for simplicity is treated as independent of grain size. Figure 4 manuscript accepted at Earth and Space Science 17 incorporates the understanding that the porosity of a rock or sediment column—and therefore overall water content—increases nearer the surface. This result shows that the effect of grain size on the H2 flux is substantial (purely as a result of mineral attenuation of radiation, and ignoring any effect of grain size on porosity). Assuming a mass density of 2.6 g/cm3 results in values of 0.24, 0.030, and 1.1 Bq/cm3, respectively, using Equation 16 and the data given in Table 1. Recognizing also that 235U contributes to the overall activity (0.011 Bq/cm3), the total activity per unit volume of the sediment is simply the sum of the individual activities, or ~1.4 Bq/cm3. From an MC perspective, it is helpful to convert the values above to algorithm-friendly percent fractions of about 18%, 1%, 2%, and 79% for 238U, 235U, 232Th, and 40K, respectively. Through this example we can see that a random number generator will select potassium predominantly over the course of a given simulation. Overall, Figure 4 demonstrates the extent to which radiolytic hydrogen production is larger in smaller grains, and so largest in clay when compared with silt or sand. If alpha radiation is used as a baseline with its comparatively large conversion coefficient ( , Equation 17), then this can be understood in terms of an average stopping distance that is at least an entire order of magnitude larger than a characteristic clay grain. By contrast, H2 flux is reduced in silt-sized grains, where characteristic grain sizes are comparable to the average stopping distance of alpha radiation (~23 ?m). In sand, the production of H2 is reduced by an entire order of magnitude compared to clay and silt for the porosity-depth relationship considered here (where porosity is itself not dependent on grain size). As before, this phenomenon can be understood in terms of characteristic grain size versus stopping distance and/or attenuation length (in the case of gamma emission). If the average stopping distance of beta radiation from natural radioactive decay is ~452 ?m, then there is the possibility that this decay mode—in conjunction with alpha—is also strongly attenuated, thereby allowing the H2 contribution due to gamma emission to dominate (attenuation length: ~3.2 cm or 32,000 ?m). It is clear from Figure 4 that the decrease in radiolytic hydrogen production as a function of depth is due in part to a higher order term in the solution to equilibrium compaction, or ??%!?"%#$. Quantitatively, normalized H2 fluxes near the surface of 100%, 61%, and 10% in clay, silt, and sand, respectively are reduced to 83%, 51%, and 9% at the maximum depth of 1 km. These values differ according to rock type such that, on average, porosity tends to be greater in coarser rocks; we do not explore this effect here (or the effect of compositional differences between clay-sized and coarser mineral grains typically encountered in siliciclastic rocks). manuscript accepted at Earth and Space Science 18 Figure 5. Radiolytic hydrogen production via alpha, beta, and gamma decay as a function of SiO2 grain size for the test case of 7.4 ?g/g U, 2.8 ?g/g Th, and 1.3% K. Top: total H2 production per unit volume of water-saturated sediment or rock. Bottom: percentage contribution of each type of radiation to the total H2 yield. Figure 5 compares the effect of grain size on H2 production due to alpha, beta, and gamma radiation for the test case presented in Figure 4. (Neither the porosity correction nor the volume correction is applied here since neither would qualitatively alter the comparison.) The top panel shows H2 production in units of pico (p, E-12) moles per [Earth] year-cm3 as a function of SiO2 grain radius in cm. Note that H2 production due to alpha decay is strongly attenuated as a function of grain size. This phenomenon is explained by the relatively short average stopping distance of alphas (~23 ?m or 2.3E-3 cm) compared to beta particles (electrons) and gamma rays (photons). Furthermore, this figure also shows that H2 production due to alpha decay mirrors the total amount over a logarithmic distribution of grain sizes from clay to sand. This is due to the preponderance of alpha particles in the radioactive emissions of the simulated material, and to the high radiolytic H2 yield per unit kinetic energy for alpha particles, which dominates that of beta particles and gamma rays by an entire order of magnitude (Spinks and Woods, 1990). Figure 5 (bottom) shows H2 production in terms of percent fraction as a function of SiO2 grain radius in cm. This figure (which would have to be adapted for simulations of materials with manuscript accepted at Earth and Space Science 19 different U, Th and K concentrations) is useful in detailing the relative contributions to H2 production as a result of alpha, beta, and gamma decay in clay, silt, and sand. In general, Figure 5 suggests that the total H2 production in this test water-SiO2 matrix is dominated by the influence of alpha decay. Beta and gamma decay, by contrast, make relatively modest contributions, where little or no attenuation is observed. In clay, H2 production from alphas in this simulation is approximately constant at 10.0 ± 0.1 pM/yr-cm3 where the grain size is smaller than the characteristic, or average, stopping distance as a result of natural radioactive decay. By contrast, H2 production due to beta and gamma decay is comparatively small, with estimates of 0.3377 ± 0.0003 and 0.5220 ± 0.0002 pM/yr-cm3, respectively. This is explained by the differing radiolytic H2 yields per unit kinetic energy of these radiation types, i.e., , , and factors. The percent contributions in clay due to alpha, beta, and gamma decay are 92%, 3%, and 5%, respectively. The contribution from alpha decay is thus compatible with Blair et al. (2007), who report an estimate of 90 ± 2 % in the case of clay-like grain sizes and similar radioactivity levels. These simulations also indicate that the main radioisotope contributor in this particular test material is 238U, with the remaining components of 235U, 232Th, and 40K producing only trace amounts of H2 irrespective of the decay mode. This figure shows that these percent contributions are roughly constant due to the comparatively weak attenuation of kinetic energy in the small grains characteristic of clay. The difference in grain size between clay and silt translates into a marked decrease in H2 yield by the latter (at a given porosity). This is not surprising given that a typical alpha particle undergoing radioactive decay in SiO2 will range on the order of E-3 cm. Figure 5 also shows that the degree of attenuation is substantial: The lower grain-size bound of silt due to alpha decay produces 9.6 pM/yr-cm3, while the upper grain-size bound produces 3.2 pM/yr-cm3. Similar to clay, H2 production in silt due to beta and gamma radiation are virtually unchanged due to their comparatively large stopping distance and attenuation length, respectively. The effect of increased grain size is more readily apparent in the percent fractions of H2 production due to alpha, beta, and gamma decay. At the lower grain-size bound of silt, alpha radiation is responsible for roughly 92% of H2 production. This reduces to 79% at the upper bound, with the percent fractions of beta and gamma increasing to 8% and 13%, respectively, as a result. In sand, the attenuation of alpha particles is more pronounced than that observed in silt. (It should be noted that the reduction of H2 production in this test water-sand matrix qualitatively resembles the escape efficiency of alpha particles seen in Figure 2.) The lower bound of sand due to alpha decay produces 3.2 pM/yr-cm3, while the upper bound produces 0.1 pM/yr-cm3; between 1 and 2 mm grain size, the production of H2 due to alpha and beta decay is comparable. The effect of alpha particle attenuation is readily apparent in sand, where the percent fraction contributions of beta and gamma radiation compensate proportionally. At the upper bound (~0.2 cm), alpha, beta, and gamma emissions contribute 14%, 18%, and 68% to H2 production, respectively. Furthermore, a local maximum in the curve corresponding to beta decay is observed at ~0.08 cm, suggesting that beta particles “range-out” in grain sizes characteristic of sand. This is consistent with the average beta-particle stopping distance of ~0.05 cm in natural radioactive decay (National Institute of Standards and Technology, 2021b). No such maximum is observed in the curve corresponding to gamma radiation, however, owing to the comparatively large average attenuation length (~3.2 cm) of this decay mode in rock/sediment, not to mention the differing physics of photon propagation in matter. Given the grain sizes explored in this study, we have naturally found that all gamma radiation escapes from the originating grain. manuscript accepted at Earth and Space Science 20 3.2 Application to Natural Systems Having validated our method for predicting the radiolytic hydrogen flux in clay, silt, and sand, in this section we consider how grain size could affect radiolytic yields from a variety of geologically interesting natural materials of known (or estimable) porosity and radionuclide concentration. Table 2 summarizes the values of key parameters (other than grain size) used here to simulate porewater radiolysis at Jezero Crater on Mars, the average upper continental crust of the Earth, monazite beach sediment found in Ullal (India), the upper half of the Peru Basin sediment column examined in Figure 3, and Upper Vredefort Crust of South Africa. Table 2. A summary of the input parameters used to generate Figures 6-10. The activities are computed using Equation 16 with a mass density of 2.7 g/cm3. The total activity per unit volume of rock/sediment is the sum of all radioisotope contributions. We assume that hard, crystalline rocks have depth-invariant porosity since they are effectively incompressible over the depth range sampled, i.e., the uppermost one kilometre of crust (Gleeson et al., 2016); for an alternative approach appropriate to larger depth intervals see Sherwood Lollar et al (2014). For the other geological settings, we set , the compaction coefficient, to 5.3E-4 m-1 on Earth and 3.5E-4 m-1 on Mars (Clifford et al. 2010). Geological (m- Reference Formation (Bq/ (%) (%) (%) (%) (% 1) cm3) ) Jezero Crater, 0.3 ~96 ?2 <<1 ?2 35 3.546 Dzaugis et al. Mars 1E-4 (2018); Clifford et al. (2010) for k Upper 2.2 ~90 ?5 <1 ?4 1 0 Aquilina et al. continental crust (2004) (Earth, granitic) Monazite beach 6.4 ~7 ~77 <1 ?15 40 5.3E- Radhakrishna sediment, Ullal, 4 et al. (1993), India except ?0 Peru Basin 2.3 ~86 ?7 <1 ?6 88 5.3E- Blair et al. sediment 4 (2007) column Upper Vredefort 3.3 ~92 ?5 <1 ?2 0.2 0 Lin et al. Crust 5 (2005) manuscript accepted at Earth and Space Science 21 The factor (decay/sec-cm3) listed in Table 2 is used in part to convert the H2 production of Equation 17 (molecules H2 / decay) to units of M/yr-cm3. The other activities listed in this table serve as “probability selectors” for the MC simulations. (For example, the radioisotope 40K, being relatively abundant in minerals, is commonly selected by the random number generator.) In general, these selectors act as proxies for the prevalence of alpha, beta, or gamma emission within the decay chains used by the simulations. The parameters and (surface porosity and compaction coefficient, respectively) are used in Athy’s law. For the upper half of the Peru Basin sediment column, we used the average of the porosity measurements above ~55 meters to extrapolate this geological formation to greater depths. manuscript accepted at Earth and Space Science 22 Figure 6. Radiolytic hydrogen production, normalized to maximum in clay, as a function of depth in a variety of geological settings. The SiO2 grain size is sampled over clay (0.06 - 3.9 ?m), silt (3.9 - 63 ?m), and sand (63 - 2,000 ?m). Similar to Figure 4, Figure 6 shows modelled H2 production, normalized to maximum in clay, as a function of sediment/rock depth for Jezero Crater on Mars, monazite beach sediment found in Ullal (India), and the upper half of the Peru Basin sediment column examined in Figure 3. As before, the grain size is randomly sampled over intervals corresponding to clay (0.06 - 3.9 ?m), silt (3.9 - 63 ?m), and sand (63 - 2,000 ?m). Data corresponding to the upper continental crust of the Earth and Upper Vredefort Crust of South Africa is not shown in this figure due to the null compaction coefficients assumed for this study, which yield constant radiolytic hydrogen fluxes as a function of depth. Relative to maximum production in clay, these fluxes drop to 67% and 68% in silt for the upper continental crust and Upper Vredefort Crust, respectively. In sand, these fluxes are reduced further to 17% and 19%. For clay, H2 production is reduced to 81% relative to maximum in Jezero Crater at a depth of 1 km; in monazite beach sediment—reduced in size to clay—a reduction to 75% is observed. The relatively high surface porosity of the Peru Basin sediment column sees an H2 production near the surface of 42% relative to maximum, which in this case occurs at approximately 1.06 km. Though not shown in Figure 6, this production will decrease at larger depths in accordance with ?z(1 - ?z). For silt and relative to maximum in clay, the production of molecular hydrogen is 72% and 65% near the surface in Jezero Crater and monazite beach sediment (size-reduced to silt), respectively. The Peru Basin sediment column displays a flux 28% near the surface by virtue of its high porosity. At maximum depth, H2 production is reduced to 58% and 49% in Jezero Crater and monazite beach sediment, respectively; the Peru Basin sediment column increases to 65% at 1 km. For sand and again relative to maximum in clay, H2 production is 27%, 8%, and 6% near the surface in Jezero Crater, monazite beach sediment, and the Peru Basin sediment column, respectively. Whereas these values are reduced at maximum depth to 22% and 6% in Jezero Crater and monazite beach sediment, respectively, by contrast H2 production increases to 14% in the Peru Basin sediment column. The reader should note the difference in H2 flux between Jezero Crater and monazite beach sediment when uranium enrichment is considered (2% vs. 15%, respectively). Since naturally occurring uranium decays by alpha emission, Figure 6 confirms the production of radiolytic hydrogen is comparatively large in clay versus sand. This observation could be relevant to uranium ore environments, in particular those found in the Witwatersrand Basin in South Africa. manuscript accepted at Earth and Space Science 23 Figure 7. Total radiolytic hydrogen production due to alpha, beta, and gamma decay as a function of SiO2 grain size in a variety of geological formations. Figure 7 shows total H2 production in units of pico (p, E-12) moles per [Earth] year-cm3 as a function of SiO2 grain radius in cm for the geological formations presented in Figure 6. This figure shows that the monazite beach sediment of Ullal, India produces an H2 flux that is three orders of magnitude larger than Jezero Crater on Mars. Consistent with the observation in Figure 6, the typical Earth group predicts comparable fluxes of radiolytic hydrogen owing to their similar overall radioactivity. This observation is largely due to their comparable concentrations of thorium and uranium, which decay primarily via alpha emission. If monazite beach sediment were reduced to a grain size in the range of clay, H2 production from all sources of radiation would be approximately 216 ± 2 pM/yr-cm3 (compared to 86 and 10 pM/yr-cm3 at the upper bounds of silt and sand, respectively). This comparatively large H2 flux stands in contrast to the typical Earth group, which displays an average of 10.5 ± 1.8 pM/yr-cm3. At the lower end, Jezero Crater on Mars gives an estimate of 0.57 ± 0.01 pM/yr-cm3. At the grain size of silt, the typical Earth group gives a reduced H2 flux of 4.7 ± 0.8 pM/yr-cm3 where grain sizes approach those characterizing sand (i.e., the upper bound of silt). Similarly, Jezero Crater reduces to 0.3 pM/yr-cm3. This phenomenon is explained by the tendency of alpha emissions to attenuate in silt where the characteristic range of these particles (~23 µm) is comparable to the grain size. The manuscript accepted at Earth and Space Science 24 slope of the decline in H2 flux with increasing grain size becomes shallower in sand. Specifically, at ~0.2 cm, radiolytic hydrogen production in the typical Earth group and Jezero Crater is reduced to 1.3 ± 0.3 and 0.1 pM/yr-cm3, respectively. Finally, since neither the porosity correction nor the volume correction is applied here, the reader should note that the ordering from highest to lowest in Figure 7 naturally differs from that of Figure 6. Figure 8. Radiolytic hydrogen production due to alpha decay as a function of SiO2 grain size in a variety of geological settings. Figure 8 (top) shows H2 production due to alpha decay in units of pico (p, E-12) moles per [Earth] year-cm3 as a function of SiO2 grain radius in cm for a variety of geological formations. Similar to Figure 5, Figure 8 (bottom) also shows H2 production in terms of percent fraction as a function of SiO2 grain radius in cm. When only alpha radiation is considered, our simulations predict that radiolytic hydrogen production would be ~207 ± 2 pM/yr-cm3 if monazite beach sediment were reduced to a grain size in the range of clay. The typical Earth group displays a middle ground of H2 flux due to alpha decay, with an average of 8.9 ± 1.7 pM/yr-cm3 in clay. Jezero Crater is expected to manuscript accepted at Earth and Space Science 25 produce a flux that is an order of magnitude below the typical Earth group, with an estimate of 0.42 ± 0.01 pM/yr-cm3. From the standpoint of alpha decay, our simulations show (in descending order) 95.81 ± 0.04 % (if reduced to clay), 88.5 ± 0.1 %, 85.1 ± 0.2 %, and 82.5 ± 0.2 % for monazite beach sediment, the Peru Basin sediment column, the average of the upper continental crust of the Earth, and Upper Vredefort Crust, respectively, with an average of 85.4 ± 3.0 % for the typical Earth group. This is explained by the 232Th and 238U parameters of Table 2, which are derived from the source thorium and uranium measurements of previous studies. (The result for Jezero Crater is slightly lower at 73.9 ± 0.3 % for the same reason.) All geological settings studied here attenuate in a predictable manner in silt-sized grains, with estimates of 76.9, 3.2 ± 0.6, and 0.1 pM/yr-cm3 for monazite beach sediment (if appropriately reduced in grain size), the typical Earth group, and Jezero Crater, respectively. Similar to Figure 5 (bottom), in silt the effect of increased grain size is visible on the percent fraction of H2 production due to alpha decay. The similar shapes of the curves observed in the typical Earth group are due to their analogous radioactive composition; at the upper bound of silt, the percent contribution from alpha decay is reduced to 67.7 ± 5.0 %. Bookending the typical Earth group are size-adjusted monazite beach sediment and Jezero Crater, which are reduced to 89.6% and 50.1%, respectively. In grain sizes characteristic of sand, Figure 8 shows that H2 production due to alpha decay is reduced substantially: At ~0.2 cm, this figure gives estimates of 2.59, 0.11 ± 0.02, and 0.005 pM/yr-cm3 for monazite beach sediment, the typical Earth group, and Jezero Crater, respectively. Similar to Figure 5 (bottom), alpha emissions from natural radioactive decay are further attenuated in sand. At the upper bound (~0.2 cm), alpha emission contributes 26.5%, 8.2 ± 1.8 %, and 4.2% to H2 production in monazite beach sediment, the typical Earth group, and Jezero Crater, respectively. As noted in Figure 5, the upper bound of sand in Figure 8 displays a tapering effect that is consistent with the escape efficiency of alpha particles seen in Figure 2. manuscript accepted at Earth and Space Science 26 Figure 9. Radiolytic hydrogen production due to beta decay as a function of SiO2 grain size in a variety of geological settings. Figure 9 (top) shows H2 production due to beta decay in units of pico (p, E-12) moles per [Earth] year-cm3 as a function of SiO2 grain radius in cm for a variety of geological formations. In contrast to Figure 8, where the attenuation of alpha particles is observed in silt, in Figure 9 this same phenomenon occurs in grain sizes that characterize sand. This is explained by the substantially larger stopping distance of beta particles that are routinely emitted via natural radioactive decay (~0.05 cm). Similar to Figure 8, this figure also shows that—across all grain sizes—the monazite beach sediment of Ullal, India produces the most H2 overall, a finding that is one order of magnitude above the typical Earth group. Figure 9 (bottom) shows H2 production in terms of percent fraction as a function of SiO2 grain radius in cm. Here a global maximum in sand is observed for all geological units. If monazite beach sediment were reduced to a grain size in the range of clay, our simulations predict that radiolytic hydrogen production would be approximately constant at 3.122 ± 0.002 pM/yr-cm3. In clay the typical Earth group displays an H2 flux due to beta decay that is one order of magnitude smaller: 0.5 ± 0.1 pM/yr-cm3. Jezero Crater on Mars is expected to produce a flux that is an order of magnitude below the typical Earth group, with an estimate of 0.04776 ± 0.00003 pM/yr-cm3. In descending order, our simulations show percent fraction manuscript accepted at Earth and Space Science 27 contributions due to beta decay of 8.37 ± 0.08 %, 4.84 ± 0.9 %, and 1.45 ± 0.01 % for Jezero Crater, the typical Earth group, and monazite beach sediment (in the form of clay), respectively. Similar to Figure 8, these results are explained by the 40K parameters of Table 2, where potassium decays the most frequently via beta emission. Figure 9 shows that monazite beach sediment, if reduced to grain sizes characteristic of silt, is predicted to produce an H2 flux due to beta emission that is only slightly diminished compared to that in clay, or 3.0 pM/yr-cm3 at the upper bound. What is more, the typical Earth group shows almost no variation both qualitatively and quantitatively (at ~0.006 cm: 0.5 ± 0.1 pM/yr-cm3). This is also observed in the case of Jezero Crater, or 0.04685 pM/yr-cm3 at the upper bound of silt. Examining Figure 9 (bottom), if monazite beach sediment were reduced to silt, the percent contribution due to beta decay would increase marginally to 3.5% at the upper bound. Furthermore, Jezero Crater and the typical Earth group also increase to 15.8% and 10.5 ± 1.4 %, respectively. Comparing these results to Figure 8 (bottom), we can understand this phenomenon in terms of the tendency of alpha particles to range-out in silt (i.e., beta emission compensates proportionally). At the upper bound (~0.2 cm), Figure 9 (top) gives estimates of 1.27, 0.20 ± 0.04, and 0.02 pM/yr-cm3 for monazite beach sediment, the typical Earth group, and Jezero Crater, respectively. Beta emissions are substantially attenuated in sand where the characteristic range of these particles (~0.05 cm) is comparable to the grain size. That is, beta particles range-out in sand in much the same way as alpha particles in silt. Moreover, it can be seen from Figure 9 that the peak contribution from beta emission occurs at about 0.04, 0.06, and 0.1 cm for Jezero Crater, the typical Earth group, and monazite beach sediment, respectively. manuscript accepted at Earth and Space Science 28 Figure 10. Radiolytic hydrogen production due to gamma decay as a function of SiO2 grain size in a variety of geological settings. Figure 10 shows H2 production due to gamma decay in units of pico (p, E-12) moles per [Earth] year-cm3 as a function of SiO2 grain radius in cm for a variety of geological formations. In contrast to alpha and beta emission, no attenuation is observed in radiolytic hydrogen production as a function of grain size. This is explained by 1) the differing physics of photon interaction with matter, and 2) the large attenuation length (~3.2 cm) characteristic of gamma decay in a natural setting. Similar to Figure 9 (bottom), this figure shows that Jezero Crater produces the most H2 across all grain sizes due to gamma decay, while monazite beach sediment produces the least, with the typical Earth group presenting a middle ground. These results are sensible given the 40K parameters of Table 2, which suggest the frequent emission of a residual gamma ray as radioactive potassium decays to 40Ar. Gamma rays are photons that lack charge and mass. Consequently they are not subject to electromagnetic slowing as in the case of alpha and beta particles. Thus, the H2 contribution of gamma rays is constant where the grain size is much less that the characteristic attenuation length found in natural radioactive decay. Our simulations confirm this phenomenon with a summary in Table 3. manuscript accepted at Earth and Space Science 29 Table 3. Average H2 production in units of pM/yr-cm3 due to gamma decay for a variety of geological settings predicted using the parameters given in Table 2. Jezero Upper Monazite beach Peru Basin Upper Crater, continental sediment, Ullal, sediment Vredefort Crust Mars crust (Earth) India column Average 0.10122 0.8509 5.928 0.9185 1.278 Standard 0.00004 0.0004 0.002 0.0004 0.001 Deviation In clay (0.06 - 4 ?m), Figure 10 shows percent fraction contributions due to gamma decay of 17.7 ± 0.2 %, 9.8 ± 2.2 %, and 2.74 ± 0.03 % for Jezero Crater, the typical Earth group, and monazite beach sediment (size-reduced), respectively, in descending order. As noted in Figure 9 (bottom), the sequencing observed here is attributed to the relative radioactivity of these geological formations. In silt (5 - 60 ?m), this figure shows that Jezero Crater, the typical Earth group, and monazite beach sediment (reduced in size to silt) increase to 34.1%, 21.8 ± 3.6 %, and 6.9%, respectively. As noted previously, this is a further consequence of the tendency of alpha particles to range-out in silt, where beta emission (Figure 9, bottom) and gamma emission (Figure 10) compensate proportionally. In sand (70 - 2,000 ?m), alpha and beta particles are shown to have a diminished effect on the percent contribution to H2 production. At this point gamma emission dominates radiolytic hydrogen production to the extent that Jezero Crater, the typical Earth group, and monazite beach sediment contribute 80.9%, 76.8 ± 2.0 %, and 60.6%, respectively, at the upper bound. Gamma rays characteristic of natural radioactive decay are expected to attenuate in grain sizes of between one and two orders of magnitude beyond those of sand. 4. Discussion and Conclusions This study has developed and validated a Monte Carlo method that simulates the physics of mineral radiation in order to predict the production of H2 as a function of mineral grain size and radioisotope composition. The results given here confirm that grain size is a major control on overall radioactivity, and consequently, on radiolytic H2 yield. The dominant mode of decay (i.e., as controlled by the radioisotope composition) also tends to influence the degree to which H2 flux is attenuated as a function of grain size; materials whose radioactivity is dominated by their potassium content will show a less severe attenuation of H2 flux with increasing grain size since 40K radiation is dominated by beta rather than alpha emissions. When our method is applied to a generic test case, the standard geological classification of grain sizes seems to suggest that clay- sized particles can produce up to an order of magnitude more H2 per unit time than sand-sized manuscript accepted at Earth and Space Science 30 particles (other things being equal). Furthermore, we have used the physics encoded in our simulation method to confirm that alpha particles range-out in silt, while beta particles range-out in grain sizes characteristic of sand. Additionally, all gamma rays, by virtue of their attenuation length, were found to escape to the surrounding water layer in our simulations, with the result that gamma radiation dominates H2 production for grain sizes larger than sand. Consequently, we do not expect gamma radiation to attenuate until grain sizes exceed those of sand by at least one order of magnitude. We have applied our numerical method to a variety of interesting natural systems as a means of exploring the magnitude of porewater radiolysis. This exploration was conducted using compositional data from real geological units in order to demonstrate the dependence of natural radionuclide concentration and bulk porosity on radiolytic hydrogen flux. A survey of our simulation results seems to suggest that H2 production in these geological units is modulated by alpha decay in rock/sediment grain sizes characteristic of clay and silt, while beta and gamma emissions take over as a secondary control in sand-sized grains. As these results were obtained for homogeneous spherical grains, it is important to note that this work has not addressed the relationships between size and shape obtaining in real geological materials. The importance of these relationships for radiolytic yields is not yet clear; in real rocks and sediments, grains of a size likely to attenuate alpha radiation considerably (e.g., sand) are also likely to be fairly spherical and equant, whereas highly non-spherical grains (e.g., platy clay crystals) are just as likely to be small and thus minimally attenuating. Equally, this work has not considered the broad correlations between grain size and porosity that exist in natural materials; these avenues could be explored in the future. Moreover, future work might also address the effect of surface-area-to-volume ratio rather than grain size on the modeled H2 production, and the highly heterogeneous distribution of porosity and mineral-water interfaces occurring in some rock types of particular astrobiological interest (e.g., finely fractured crystalline rocks or vesicular basalts; McMahon et al., 2013). It would also be interesting to examine the possible release of radiolytic H2 from hydrated or hydroxylated mineral grains themselves (which may tend to be more abundant in the clay size-fraction) rather than the porewater. Ultimately, our Monte Carlo method can be readily adapted for materials with constituent grains and pore spaces of arbitrary size and shape, of any composition, and on any planetary body. As a consequence, it should now be possible to estimate more accurately how radiolytic hydrogen might contribute to the habitability of microbial life in porous materials on Earth, Mars, and beyond. Acknowledgments The authors greatly appreciate the helpful discussions with Eric Benton at Oklahoma State University and Regina DeWitt at East Carolina University. S.M. also thanks Barry B. McMahon for helpful discussions about the geometry in Supplementary Information S2. Open Research Data Availability Statement manuscript accepted at Earth and Space Science 31 The radiolytic hydrogen production data used for the grain size analysis in this study are available via CC BY-SA 4.0 (DeWitt et al., 2022). The Fortran program used for the simulation of porewater radiolysis is available via CC BY-NC-ND 4.0 and developed openly at East Carolina University (DeWitt, 2022). References Aitken, M.J. (1985). Thermoluminescence Dating. Orlando: Academic Press. Aquilina, L., de Dreuzy, J.R., Bour, O., and Davy, P. (2004). Porosity and fluid velocities in the upper continental crust (2 to 4 km) inferred from injection tests at the Soultz-sous-Forêts geothermal site. Geochimica et Cosmochimica Acta, 68(11), 2405-2415. Bethe, H. (1930). Zur theorie des durchgangs schneller korpuskularstrahlen durch materie. Annalen der Physik, 397(3), 325-400. Bethe, H. and Heitler, W. (1934). On the stopping of fast particles and on the creation of positive electrons. Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character, 146(856), 83-112. Blair, C.C., D'Hondt, S., Spivack, A.J., and Kingsley, R.H. (2007). Radiolytic hydrogen and microbial respiration in subsurface sediments. Astrobiology, 7(6), 951-970. Brookhaven National Laboratory (2021). Interactive Chart of Nuclides (2021), http://www.nndc.bnl.gov/nudat2 Bouquet, A., Glein, C.R., Wyrick, D., and Waite, J.H. (2017). Alternative energy: production of H2 by radiolysis of water in the rocky cores of icy bodies. The Astrophysical Journal Letters, 840(1), L8. Chapelle, F.H., O'Neill, K., Bradley, P.M., Methé, B.A., Ciufo, S.A., Knobel, L.L., and Lovley, D.R. (2002). A hydrogen-based subsurface microbial community dominated by methanogens. Nature, 415(6869), 312-315. Clifford, S.M., Lasue, J., Heggy, E., Boisson, J., McGovern, P., and Max, M.D. (2010). Depth of the Martian cryosphere: Revised estimates and implications for the existence and detection of subpermafrost groundwater. Journal of Geophysical Research: Planets, 115(E7). DeWitt, J., McMahon, S., and Parnell, J. (2022). Figure data and analysis for "The Effect of Grain Size on Porewater Radiolysis" [Dataset]. Dryad. https://doi.org/10.5061/dryad.qfttdz0j1 DeWitt, J. (2022). Alpha: A Fortran program for simulating porewater radiolysis [Dataset]. Dryad. https://doi.org/10.5061/dryad.tdz08kq0z. Alpha v1.0.0 [Software]. Zenodo. https://doi.org/10.5281/zenodo.5553572 Draganic, I.G. and Draganic, Z.D. (1971). The Radiation Chemistry of Water, Academic Press, New York. manuscript accepted at Earth and Space Science 32 Dzaugis, M.E., Spivack, A.J., Dunlea, A.G., Murray, R.W., and D’Hondt, S. (2016). Radiolytic hydrogen production in the subseafloor basaltic aquifer. Frontiers in microbiology, 7, 76. Dzaugis, M., Spivack, A.J., and D'Hondt, S. (2018). Radiolytic H2 production in martian environments. Astrobiology, 18(9), 1137-1146. Einstein, A. (1905). Über einen die Erzeugung und Verwandlung des Lichtes betreffenden heuristischen Gesichtspunkt Annalen der Physik 17, 132. Farley, K. A., Wolf, R. A., and Silver, L. T. (1996). The effects of long alpha-stopping distances on (U?Th)/He ages. Geochimica et Cosmochimica Acta, 60(21), 4223-4229. Fermi, E. (1940). The ionization loss of energy in gases and in condensed materials. Physical Review, 57(6), 485. Fowler, A.C. and Yang, X.S. (1998). Fast and slow compaction in sedimentary basins. SIAM Journal on Applied Mathematics, 59(1), 365-385. Freund, F., Dickinson, J.T., and Cash, M. (2002). Hydrogen in rocks: an energy source for deep microbial communities. Astrobiology, 2(1), 83-92. Gleeson, T., Befus, K.M., Jasechko, S., Luijendijk, E., and Cardenas, M.B. (2016). The global volume and distribution of modern groundwater. Nature Geoscience, 9(2), 161-167. Gould, R.J. and Schréder, G.P. (1967). Pair production in photon-photon collisions. Physical Review, 155(5), 1404. Hofmann, B.A. (1992). Isolated reduction phenomena in red beds: a result of porewater radiolysis. In Y.K. Kharaka, and A.S. Maest (Eds.), Water–Rock Interaction (pp 503–506). Rotterdam: A.A. Balkema. Ketcham, R. A., Gautheron, C., and Tassan-Got, L. (2011). Accounting for long alpha-particle stopping distances in (U–Th–Sm)/He geochronology: Refinement of the baseline case. Geochimica et Cosmochimica Acta, 75(24), 7779-7791. Klein, O. and Nishina, Y. (1929). Über die Streuung von Strahlung durch freie Elektronen nach der neuen relativistischen Quantendynamik von Dirac. Zeitschrift für Physik, 52(11), 853-868. Kral, T.A., Brink, K.M., Miller, S.L., and McKay, C.P. (1998). Hydrogen consumption by methanogens on the early Earth. Origins of Life and Evolution of the Biosphere, 28(3), 311-319. Krumholz, L.R., McKinley, J.P., Ulrich, G.A., and Suflita, J.M. (1997). Confined subsurface microbial communities in Cretaceous rock. Nature, 386(6620), 64-66. Lefticariu, L., Pratt, L.A., LaVerne, J.A., and Schimmelmann, A. (2010). Anoxic pyrite oxidation by water radiolysis products—a potential source of biosustaining energy. Earth and Planetary Science Letters, 292(1-2), 57-67. manuscript accepted at Earth and Space Science 33 Lin, L.H., Hall, J., Lippmann?Pipke, J., Ward, J.A., Sherwood Lollar, B., DeFlaun, M., et al. (2005). Radiolytic H2 in continental crust: nuclear power for deep subsurface microbial communities. Geochemistry, Geophysics, Geosystems, 6(7). McCollom, T.M. (1999). Methanogenesis as a potential source of chemical energy for primary biomass production by autotrophic organisms in hydrothermal systems on Europa. Journal of Geophysical Research: Planets, 104(E12), 30729-30742. McMahon, S., Parnell, J., Ponicka, J., Hole, M., and Boyce, A. (2013). The habitability of vesicles in martian basalt. Astronomy & Geophysics, 54(1), 1-17. McMahon, S. and Parnell, J. (2014). Weighing the deep continental biosphere. FEMS Microbiology Ecology, 87(1), 113-120. McMahon, S., Parnell, J., and Blamey, N.J. (2016). Evidence for seismogenic hydrogen gas, a potential microbial energy source on Earth and Mars. Astrobiology, 16(9), 690-702. Murray, K. E., Orme, D. A., and Reiners, P. W. (2014). Effects of U–Th-rich grain boundary phases on apatite helium ages. Chemical Geology, 390, 135-151. National Institute of Standards and Technology (2021a), Stopping-power and range tables for helium ions, physics.nist.gov/PhysRefData/Star/Text/ASTAR.html National Institute of Standards and Technology (2021b), Stopping-power and range tables for electrons, physics.nist.gov/PhysRefData/Star/Text/ESTAR.html Nixon, S., Cousins, C.R., and Cockell, C. (2013). Plausible microbial metabolisms on Mars. Astronomy & Geophysics, 54, 1.13-1.16. Lane, N. and Martin, W.F. (2012). The origin of membrane bioenergetics. Cell, 151(7), 1406- 1416. Preiner, M., Igarashi, K., Muchowska, K.B., Yu, M., Varma, S.J., Kleinermanns, K., et al. (2020). A hydrogen-dependent geochemical analogue of primordial carbon and energy metabolism. Nature ecology & evolution, 4(4), 534-542. Radhakrishna, A.P., Somashekarappa, H.M., Narayana, Y., and Siddappa, K.A. (1993). New natural background radiation area on the southwest coast of India. Health Physics 65(4), 390– 395. Schulte, M., Blake, D., Hoehler, T., and McCollom, T. (2006). Serpentinization and its implications for life on the early Earth and Mars. Astrobiology, 6(2), 364-376. Sherwood Lollar, B., Onstott, T. C., Lacrampe-Couloume, G., & Ballentine, C. J. (2014). The contribution of the Precambrian continental lithosphere to global H2 production. Nature, 516(7531), 379-382. manuscript accepted at Earth and Space Science 34 Sleep, N.H., Bird, D.K., and Pope, E.C. (2011). Serpentinite and the dawn of life. Philosophical Transactions of the Royal Society B: Biological Sciences, 366(1580), 2857-2869. Smetannikov, A.F. (2011). Hydrogen generation during the radiolysis of crystallization water in carnallite and possible consequences of this process. Geochemistry International, 49(9), 916- 924. [Sonzogni and Shu, 2020] Sonzogni, A. and Shu, B. (2020), Nudat 2, www.nndc.bnl.gov/nudat2/ Spinks, J.W.T. and Woods, R. J. (1990). An Introduction to Radiation Chemistry. New York: Wiley. Stevens, T.O. and McKinley, J.P. (1995). Lithoautotrophic microbial ecosystems in deep basalt aquifers. Science, 270(5235), 450-455. Tarnas, J.D., Mustard, J.F., Sherwood Lollar, B., Bramble, M.S., Cannon, K.M., Palumbo, A.M., and Plesa, A.C. (2018). Radiolytic H2 production on Noachian Mars: implications for habitability and atmospheric warming. Earth and Planetary Science Letters, 502, 133-145. Tarnas, J.D., Mustard, J.F., Sherwood Lollar, B., Stamenkovi?, V., Cannon, K.M., Lorand, J.P., et al., 2021. Earth-like Habitable Environments in the Subsurface of Mars. Astrobiology, 21(6), 741–756. Taylor, J.R. (1996). An Introduction to Error Analysis: The Study of Uncertainties in Physical Measurements, 2nd Edition., Sausalito, CA: University Science Books. Waite, J.H., Glein, C.R., Perryman, R.S., Teolis, B.D., Magee, B.A., Miller, G., et al. (2017). Cassini finds molecular hydrogen in the Enceladus plume: evidence for hydrothermal processes. Science, 356(6334), 155-159. Weaver, B.A. and Westphal, A.J. (2002). Energy loss of relativistic heavy ions in matter. Nuclear Instruments and Methods in Physics Research Section B: Beam Interactions with Materials and Atoms, 187(3), 285-301. Webster, C.R., Mahaffy, P.R., Atreya, S.K., Flesch, G.J., Mischna, M.A., Meslin, P.Y., ... and Lemmon, M.T. (2015). Mars methane detection and variability at Gale crater. Science, 347(6220), 415-417. Wentworth, C.K. (1922). A scale of grade and class terms for clastic sediments. The Journal of Geology, 30(5), 377-392. Westbrook, M. L., Sindelar, R. L., and Fisher, D. L. (2015). Radiolytic hydrogen generation from aluminum oxyhydroxide solids: Theory and experiment. Journal of Radioanalytical and Nuclear Chemistry, 303(1), 81-86.