Lags in Desorption of Lunar Volatiles

Monte Carlo simulations of gas motion inside a granular medium are presented in order to understand the interaction of lunar gases with regolith and improve models for surface-boundary exospheres, a common type of planetary atmosphere. Results demonstrate that current models underestimate the lifetime of weakly bonded adsorbates (e.g., argon) on the surface by not considering the effect of Knudsen diffusion, and suggest that thermal desorption of adsorbates should be modeled as a second-or-higher-order process with respect to adsorbate coverage. An additional discrepancy between present models and outgassing from a realistic porous boundary is found for surface-adsorbate systems containing a distribution of activation energies (e.g., water). In that case, the mobility of adsorbates between desorption events (i.e., surface diffusion), not considered in global models of the exosphere, controls their surface residence time via transitions between sites of low and high binding energy. Without mobility the equatorial surface retains more water over a lunar day because sites of low binding energy are not repopulated by motion along the grain surface when depleted. The effects of Knudsen and surface diffusion apply to other volatile species and help us partly understand why measurements of lunar exosphere constituents appear to indicate stronger bonding of gas with the lunar surface than measured in some laboratory experiments.


Introduction
Because the exosphere of the Moon is a macroscopic manifestation of interactions between gas atoms and the surface (e.g., Stern 1999), exospheric measurements can be used to extract information about the gas-surface bond. With some exceptions (H 2 , He, Ne), gas particles accumulate, or adsorb, on the cold surface at lunar night and desorb from the warmer dayside soils, with the local time of the density peak constraining the temperature dependence of the desorption rate and the strength of the bond. Exosphere models are tasked with retrieving the strength of gas-surface interaction from such measurements by repeatedly simulating adsorption and desorption events until loss of these atoms and molecules by photoionization or dissociation. They tend to overestimate the desorption rate because they do not account for the porosity of the powdery surface, and must compensate by using a higher effective activation energy for desorption than measured in laboratories to fit exospheric measurements (e.g., Grava et al. 2015;. This has been interpreted as indicating that experiments conducted in laboratories do not represent the pristine conditions found on the Moon. While this conclusion may be fair, we demonstrate that an improved treatment of the porous surface in exosphere models reduces these discrepancies. Furthermore, we show that accurate knowledge of the distribution of adsorption energies from experiments is a necessary but not sufficient condition for high-fidelity exosphere models. Additional uncertainty is introduced by the role of surface diffusion, whose effect on desorption has been ignored by previous exosphere models. We used two exemplary species, argon and water, to illustrate the range of outgassing delays caused by diffusion. Argon bonds very weakly with the surface (e.g., Grava et al. 2015;Kegerreis et al. 2017, and references therein), while the interaction between water and the lunar surface is stronger (Poston et al. 2015;Jones et al. 2020a). Argon gas is a product of radioactive decay of 40 K in the lunar soil and has repeatedly been measured in the lunar exosphere, both during the Apollo days and more recently by the Lunar Atmosphere and Dust Environment Explorer (LADEE; Benna et al. 2015;Hodges & Mahaffy 2016). To interpret the existing data Kegerreis et al. (2017) implemented a simple model to account for diffusion of argon adsorbates. They assumed that some argon atoms penetrate randomly into the ground at night, only to reemerge during the daytime with a uniform distribution in delay times when the regolith warms. The argon density peak at sunrise, and the day to night ratio near sunrise, could be simulated if the trapping fraction of a few parts per thousand atoms was assumed on every bounce. They concluded that diffusion provides an alternative to assuming very high activation energies for desorption Hodges & Mahaffy 2016). A sophisticated model of a granular surface is used here to quantify how Ar moves between the subsurface and vacuum.
Desorption of adsorbed water has been presumed to be the cause of slope changes in reflectance spectra measured over a lunar day by the Lunar Reconnaissance Orbiter (LRO; Hendrix et al. 2019). Water gas may be generated by recombinative desorption of solar wind-implanted hydroxyl (Jones et al. 2018), and by meteoroid impacts, as shown by the detection of watergroup products in the lunar exosphere during meteor showers (Benna et al. 2019). The gas form, water or hydroxyl, was not determined from these LADEE measurements. If water, the LADEE measurements imply that the water does not accumulate on the nightside and it lacks signatures of repeated adsorptiondesorption cycles, perhaps because of dissociative adsorption to the regolith (Gun'ko et al. 1998). Whether or not water is continuously present in the lunar exosphere, in this study we use water as an example of an adsorbate that interacts with the lunar surface with a wide distribution of binding energies (Jones et al. 2020a). On such heterogeneous surfaces the residence time of different molecules on the surface will differ considerably. By analogy to volume diffusion controlling the timescale of hydrogen release from grain rims when they contain a distribution of crystal defects (e.g., Farrell et al. 2017), surface diffusion has been hypothesized to slow thermal desorption of adsorbates from heterogeneous planetary surfaces (Killen et al. 2018). This effect has not been considered in exosphere models and is quantified for the first time in this study.

Model Formulation and Kinetic Parameters
A detailed description of the processes taking place on an isolated astrophysical grain is provided by He & Vidali (2014). Additionally, when a grain is enveloped in a grain pile, some desorption events lead to readsorption and some lead to diffusion deeper into the powder. And, finally, surface diffusion, i.e., the thermal migration of adsorbates between adjacent sites, or local minima of the adsorption potential on a grain, enables particles to find grain contact points and slowly travel onto other grains (e.g., Sarantos & Tsavachidis 2020).
Three random sphere-packings were computer-generated (Kloss et al. 2012) to simulate the disordered structure of lunar regolith. Each packing had width of 1 × 1 mm 2 , depth of ∼5-7 mm, and consisted of 30,000 spherical grains. The grain size distribution for each packing was based on three Apollo and Luna samples that were representative of size distributions for mature (72141), submature (24109), and immature (67481) lunar soils. 3 The samples were truncated at one standard deviation, with a most likely particle diameter of ∼29 μm, a mean particle size of 44 μm, a maximum grain size of 374 μm, and a void fraction of 0.5 for the sample of intermediate maturity. Although piles of spherical grains provide better representation of the granularity of lunar regolith for gas calculations than the current state-of-the-art models, most lunar regolith particles are not spheres, and future work should include more appropriate particle geometries.
An obstructed random walk of 20,000 randomly distributed test particles was initiated in these packings to simulate gas transport. The spaces between grains in these simulations provided realistic diffusion paths. Similar approaches have previously been used to study gas transport in porous media under atmospheric pressures (e.g., Zalc et al. 2003;Hlushkou et al. 2013), but under lunar conditions the test particles collide only with grains or escape to vacuum. Particles were either deposited all at one time for isothermal simulations (Figure 1), or with a time profile for temperature-programmed desorption (TPD) simulations (Figures 2-3). Particles in the continuous dosing scenario (TPD runs) were uniformly distributed in time from sunset to sunrise to simulate adsorption during the long lunar night. For time-dependent runs the temperature was changed for every time step in order to simulate the heating rate experienced by patches of the lunar surface. Continuous linear heating or cooling from a lunar thermal model (Hurley et al. 2015) was used, but with a more gradual temperature ramp-up near the terminator due to realistic shadowing (Williams et al. 2017). No lateral or vertical temperature gradients were assumed, and their additional effects on the lifetime of adsorbates (e.g., Davidsson & Hosseini 2021) should be assessed in future work. Particles were removed when they escaped the medium, and particles remaining in the grain pile were permitted to undergo diffusion. Desorption events were recorded. We thus calculated a realistic initial condition for the distribution of adsorbates with depth as the surface element arrives at the dawn terminator.
On every time step we separate the adsorption, desorption, and surface diffusion events. At the beginning of the step we insert new particles (for the time-dependent dosing simulations) and we test for each particle whether desorption is to take place. If desorbed, we mark with ray tracing whether the particle hits another grain, including perhaps many successive collisions with other grains until it comes to a stop (Knudsen diffusion), and if it escapes to vacuum it is removed from the simulation. If no desorption occurs, we check the probability for surface diffusion, and for those particles that will undergo hopping events we update the binding energy and location. When the sites have different adsorption energies, we randomly select the binding energy after each sticking event or surface diffusion hop by importance sampling of the binding energy distribution. This means that there is no correlation between activation energies of successive jumps. The distance traveled along the spherical grain from an adsorbateʼs current position within a time step of 1-4 s is simulated probabilistically as described in Sarantos & Tsavachidis (2020) (see also Hołyst et al. 1999).
Sticking is one of the parameters enhancing the retentiveness of soil. Higher sticking probability per collision of a desorbed atom with intervening grains leads to readsorption onto adjacent grains, prolonging the time to escape from the granular bed. The adsorption probability during single collisions with grains is described by a sticking coefficient, S, while the probability of reflection is 1-S. Values of S for gases on lunar regolith have not been reported, so we used values for simple metal oxides. The sticking probability of Ar on MgO decreases exponentially with gas incident energy, E, for gas speeds >250 m s −1 (Dohnalek et al. 2002). Using their empirical relation S(E) = exp −a(E-E0) , where a = 0.202 mol kJ −1 and E 0 = 2.08 kJ mol −1 , and integrating over a Maxwell-Boltzmann distribution for thermally accommodated gas, we find population-averaged sticking probabilities of S ≈ 1 at T = 50 K, S ≈ 0.76 at T = 200 K, and S ≈ 0.55 at T = 400 K. The high sticking probabilities used in our simulation are consistent with the efficient energy accommodation of lunar argon as indicated by LADEE measurements of the scale height (Hodges & Mahaffy 2016). The sticking probability of water on CaO single crystals is constant at temperatures 120-300 K and near 0.85-0.9 (Seifert et al. 2021). Similarly, a temperature-independent sticking coefficient near unity has been measured for low-coverage water on MgO single crystals at temperatures 100-250 K (Stirniman et al. 1996). In the results presented here we assumed S = 0.85 for single collisions of water with the lunar surface at all temperatures.
The adsorption time of a test particle between desorption events is calculated as −ln(1−ξ)/R des , where R des is the desorption rate and ξ a random number that is uniformly distributed between zero and one. This numerical scheme generates exponentially distributed waiting times for desorption with mean residence time τ = 1/R des . The desorption rate increases with grain temperature, T, R des = v 0 × exp(−E des /KT). Desorption parameters for Ar were adopted from the work of Grava et al. (2015), who simulated the decay of argon density from sunset to sunrise using Apollo data.
The surface interaction was empirically described by an unusually low pre-exponential factor v 0 = 2 × 10 9 s −1 and E des = 0.281 eV, values adopted here, although other values of v 0 and E des will give similar results (Kegerreis et al. 2017). The kinetic parameters describing the interaction of water with lunar grains were obtained by Jones et al. (2020a). A pre-exponential factor of v 0 = 10 13 s −1 and a wide distribution of binding energies between 0.6 and 1.9 eV, with the most likely value 0.63 eV, were derived from TPD experiments on Apollo Highlands sample 14163, and were adopted in these simulations.
The surface diffusion of adsorbates is simulated as a random walk between neighboring sites of potential energy with a jump rate, R hop = v 0 × exp(−E diff /KT). The barrier to surface diffusion, E diff , is uncertain and is usually expressed as E diff = αE des , where α < 1 is a surface corrugation ratio. For complex surfaces and weakly adsorbed species α ∼ 0.7-0.8 (e.g., Perets et al. 2007). Values around 0.5-0.6 are reported for gases adsorbing on glass (Gilliland et al. 1974). Theoretical estimates of Ar tracer diffusion on oxides, carried out by Riccardo & Steele (1996), suggested α ≈ 0.45, the value assumed for the argon simulations here. As to water, its mobility on complex oxide surfaces is uncertain. On the one hand, high mobility of adsorbed molecules along the surface is expected from the finding that their entropies near desorption are 2/3 of their gas phase entropies, implying that only motion perpendicular to the grain surface is restricted (Campbell & Sellers 2012;Weaver 2013). Theoretical calculations on some model oxide surfaces confirm this expectation for water. For instance, the motion of an isolated water molecule adsorbed on MgO was estimated by McCarthy et al. (1996), indicating that adsorbed water molecules are very mobile on this surface (α ≈ 0.25). On the other hand, water molecules on many oxides are partially dissociated, and the dissociation products are shown to be immobile (e.g., Matthiesen et al. 2009). Yet for dissociation to occur, some mobility of water molecules is necessary to even reach the active sites (Schaub et al. 2001). For these reasons, a wide range of values (α = 0.25-1) for water surface diffusion have been assumed in our work. When there is a distribution of sites this ratio, α, was assumed to be constant at all sites regardless of the depth of the adsorption well.
The model of gas-solid interaction adopted here is suitable to extended regions of the Moon. Adsorbate-adsorbate repulsive or attractive interactions can be neglected if the average separation between adsorbates is large. Hendrix et al. (2019) inferred that less than 0.01 of a water monolayer (1 ML ∼ 10 15 cm −2 ) exists on the lunar surface. The adsorbate abundance, C ads , in flux balance with the argon gas density, n g , can be estimated as C ads = n g × v g × t res . To support a gas density n g = 2 × 10 4 cm −3 (Benna et al. 2015) for a thermally accommodated gas of mean speed v g with the Grava et al. (2015) parameters for the surface residence time, t res , the surface coverage of adsorbed argon must be C ads ∼ 10 12 -10 13 cm −2 , i.e., 1 in every 100-1000 sites near sunrise are occupied. In many gas-surface systems abundances higher than a tenth of a monolayer must be inferred before coverage-dependent kinetic parameters must be considered (e.g., Zhdanov 1991). Furthermore, gas atoms and molecules are deposited not only on the top 1 or 2 grains of the porous regolith, but up to 10 grains into the subsurface (Figure 1(a); see also Kulchitsky et al. 2018). Thus, the grains are essentially bare over large portions of the lunar surface, and the noninteracting test particle approach is valid.

Adsorbed Argon
The left panel of Figure 1 shows how the argon reservoir evolves with time during an isothermal simulation at T = 115 K. Because of the gaps between grains and because some atoms do not stick the first time they meet a grain, deposited particles (black line, t = 0) can reach significant depths compared to the mean grain size. About 10% of all argon test particles are initially deposited at a depth exceeding 200 μm. With diffusion between grain voids this fraction rises to ∼15% within one day as desorption events occur. This population represents a "fat tail" of slowly desorbing particles that, as we will show next, do not follow exponentially distributed waiting times for desorption from the grain pile. Due to this motion the gas removal rate from a powder is slower than an exponential decay, and is better approximated by a second-order process.
The right panel of Figure 1 shows the cumulative probability of argon desorption (desorbed fraction) with time. The black line shows the outgassing from our simulation of the porous medium. A first-order model ( µ dN dt N), where N is the number of test particles remaining in the packing, was used to simulate the outgassing from a single grain (blue line). Defining lag as the time it takes to remove 63.2% of the adsorbate reservoir, the porous soil takes 6 to 7 times longer to reach full effect than the exponential decay with the same microphysical parameters. A reduction of the rate by a factor of 4 when a grain is enveloped in a pile of grains somewhat extends the applicability of the exponential decay to longer timescales (red line). The sum of two exponential decays, one fast and one slow, did not improve the fitting. Thus, thermal desorption is not an exponential decay: because of diffusion, the release is decelerating. A second-order model ( µ dN dt N 2 ) better estimates the average residence time of adsorbates, although it fails for the last ∼30% of particles at depth (yellow line), while a third-order model (magenta line) predicts too slow a depletion of the surface reservoir. In conclusion, a kinetic model with order between two to three with respect to coverage should be adopted to model thermal desorption from a powder. Figure 2 shows the desorption of sequestered argon as the porous surface approaches and crosses sunrise. Time-dependent simulations were performed at two different latitudes, equatorial (left) and 45°north (right). In black is the outgassing from the grain pile under the temperature heating profile shown in gray. Also shown are the TPD curves from a smooth surface (exponential decay with Eb = 0.281 eV, blue), and a model with reduced rate as established in the isothermal simulation (red), which both show the asymmetric peak shape typical of first-order TPD curves. The firstorder model with a reduced rate correctly predicts the timing of the peak at sunrise, but overestimates its amplitude (the slow ramp-up in outgassing from the powder). Additionally, a porous surface outgasses for up to 10 hr longer than a smooth surface. Not only is the symmetric TPD curve from the grain pile qualitatively consistent with second-order desorption (magenta line), but also we confirmed quantitatively with an Arrhenius plot of ln(R/N n ) versus 1/T (De Jong & Niemantsverdriet 1990), where n is the order of desorption, that second or even higher desorption order is justified. The effect of increasing the sticking coefficient is an increase in the order of desorption.
What is the mechanism causing slow desorption? In order to appreciably diffuse inward during the long lunar night via the slower surface diffusion mechanism, an argon atom would require a low activation barrier for surface diffusion. However, in time-dependent simulations that better approximate the expected depth distribution of adsorbates at sunrise, we found only small differences in the outgassing profile when α 0.25. This finding suggests that the initial distribution of particles with depth, followed by Knudsen diffusion, were the main cause of the delay in desorption. The increase of Knudsen diffusion into porous media at the onset of desorption has been demonstrated experimentally (Ballinger et al. 1989). Outward diffusion of an adsorbate through MgO powders with nanosized interparticle voids showed both a fast and a slow diffusion channel, suggestive of bonding on two types of sites, and thermal desorption followed a symmetric TPD curve unlike a first-order curve (Kim et al. 2009, their Figure 4). We suggest that these experimental results are similar to our findings that desorption from a powder is slower than an exponential decay of the adsorbate reservoir even for well-defined sites.
While Kegerreis et al. (2017) correctly identified diffusion as a rate-limiting process, we do not confirm their implementation of diffusive effects. Our recommendation is a change in the adopted desorption order. Second-order kinetics, which are usually associated with recombinative desorption, describe the time profile of diffusion-controlled gas release from a grain pile even for monoatomic gases. Furthermore, with this formulation we retrieve meaningful desorption energies (i.e., the single-grain attachment of 0.281 eV), whereas if we treated desorption as a first-order process, we would erroneously infer a distribution of binding energies for this gas-surface system from the TPD curve (see Sarantos & Tsavachidis 2020). Figure 2. Time-dependent simulations of argon, accumulated during the lunar nighttime onto a computer-generated powder, illustrate the slow release of sequestered argon from porous regolith at sunrise (black). Simulations of exponential decay (blue and red), the current standard used in exosphere models, fail to reproduce the long tail of the regolith outgassing. A second-order model (magenta) better describes the diffusion-controlled removal of gas near the terminator as the regolith is heated (right axis). Figure 3 illustrates the effects of porosity and diffusion on the timing of the desorption of deposited water. For these simulations we used the activation energy of water release from a Highlands sample 14163 (Jones et al. 2020a). For the surface diffusion barrier we assumed α = 0.25-0.99. To highlight the porosity effects we also used the equations in Barrie (2008) to illustrate the expected TPD curve from a flat surface with this binding energy distribution under the assumptions of first-order desorption and no surface diffusion. The latter are the typical assumptions in exosphere models.

Adsorbed Water
The simulated TPD spectrum, recorded as the total number of desorption events every 1000 s, is shown in Figure 3(a). For rapid surface diffusion values (α 0.5), shown in blue, we get a TPD curve similar to argon. In this limit shallow sites are immediately upon desorption replenished by diffusion from the active sites, accelerating desorption. For the same distribution of desorption energies but intermediate values of surface diffusion (α = 0.7), shown in black, the desorption shifts to higher temperature and exhibits a tail of long-lived adsorbates. Our result confirms the analysis by Xia et al. (2007), who predicted that the desorption rate in a TPD experiment from surfaces with a broad energy distribution should peak at lower temperatures when surface diffusion is considered. A case with equal diffusion and desorption barriers is shown in red, and the case for first-order desorption from a flat surface with no diffusion is shown in green (Barrie 2008). Contrasting these four cases, it can be concluded that high migration barriers (α > 0.5) and readsorption reduce the outgassing rates.
The cumulative desorbed fraction (Figure 3(b)) shows a sigmoid curve characteristic of deadtime, while the powder retains different amounts of water depending on the surface diffusion assumption. With no surface diffusion (red curve), the surface begins to outgas at lower temperature but overall retains more molecules over a full thermal cycle. Up to 9% of the deposited water test particles were retained for an entire lunar day at equatorial latitudes when no surface diffusion and the rate reduction due to porosity were simulated, compared to 2% retention from a single grain and/or smooth surface (green). (This high retention efficiency reflects the fraction of sites with initial binding energies >1.5 eV in the Jones et al. 2020a data, readsorption because of the high sticking coefficient, and Knudsen diffusion deeper into the grain pile.) For diffusion barriers lower than the desorption barrier (blue and black lines), there is an additional deadtime of ∼3 × 10 4 s (∼10 hr, or ∼6°in local time) before desorption commences as shallow sites populate active sites via surface diffusion at early morning temperatures. Although this migration delays the onset of desorption, mobile atoms, which would be trapped throughout a lunar cycle if immobile, can at higher temperatures migrate to shallow sites and desorb repeatedly, like an avalanche, overcoming loss to readsorption. Hence, over a lunar day the surface is less retentive at equatorial and midlatitudes for this gas-substrate system when the surface diffusion barrier is much lower than the desorption barrier. On the other hand, for latitudes experiencing maximum dayside temperatures less than 190 K the initial deadtime from surface diffusion results in a surface being more retentive for this distribution of binding energies. We note that the desorption curves are approximately independent of the assumed grain size distributions (right panel of Figure 3), although their asymptotic content may differ.
Our conclusions on the role of surface diffusion apply regardless of the assumed distribution of binding energies provided that there is a tail of chemisorption sites. This statement was confirmed by repeating these simulations using the energy distributions of Jones et al. (2020a) from the mare sample 10084, which produced similar sigmoid curves (not shown here). Furthermore, the absence of high-temperature desorption when high adsorbate mobility is assumed is similar to that observed by Jones et al. (2020a) from the Mare sample, Figure 3. Porosity and surface diffusion complicate the desorption of adsorbed molecular water from regolith simulants. All cases adopted the same distribution of desorption energies, E des , but varied the barrier for surface diffusion, E diff . Solid lines show desorption from the submature soil, whereas the dotted and dashed lines show the effect of decreasing and increasing maturity through the particle size distribution. suggesting that differences in surface diffusion be considered as one of the processes that contribute to the differences between Mare and Highlands soils (Poston et al. 2015;Jones et al. 2020a).

Conclusions
We find two sources of gas lag corresponding to two forms of diffusion. When all sites are equivalent, the desorption lag is due to Knudsen diffusion followed by readsorption. The time dependence of the adsorbate reservoir beyond the half-life, and hence the exosphere longevity, is not correctly described by an exponential decay because thermal desorption and Knudsen diffusion form a negative feedback loop inside a powder (Figures 1-2). To describe the removal of gas and derive meaningful activation energies for desorption without explicitly treating the multiple types of diffusion covered here, expressing the gas release as a second-order reaction with respect to surface concentration is a reasonable approximation in global exosphere models.
We suggest that the thermal desorption flux be modeled as F out (t) ≈ 0.25 × R des × σ × C ads (t) 2 , where F out is the instantaneous outgassing flux corresponding to an adsorbate coverage C ads , R des is the single-grain desorption rate, and σ = 10 −15 cm 2 is the cross section of an adsorption site. The factor 0.25 denotes the reduction of the yield from flat surfaces due to readsorption (Figures 1-2). This parametric change will describe how volatile adsorbates dissipate after sudden changes (e.g., crossing the terminator, following a natural or manmade impact, or following the descent of a lunar lander).
Surface diffusion, another diffusion mechanism ignored in exosphere models, has measurable macroscopic effects. For strongly bound gases, such as alkalis at the Moon, surface diffusion assists adsorbates to escape microscopic shadows, enabling more adsorbates to be degassed by photons but limiting the photodesorption rates (Sarantos & Tsavachidis 2020). For volatiles the effect of surface diffusion is more pronounced when a wide distribution of binding energies is assumed. Thermal hopping of adsorbates between physisorption (shallow) and chemisorption (active) sites, i.e., diffusion in the energy dimension, affects the desorption rate in a complicated way, suppressing desorption at lower temperatures (He & Vidali 2014), but enhancing desorption over a lunar day depending on the maximum dayside temperature of the soil. With this redistribution mechanism, desorption at sunrise is slower and occurs at later local times for species like water (Figure 3). More generally, when surface diffusion has a lower barrier than desorption, the occupancy of desorption energies is not evolved correctly with time in exosphere models, and the modeling error (e.g., local time of the peak density) grows as a wider distribution of binding strengths is assumed. Kinetic parameters for adsorbate mobility and desorption are equally needed from experiments.
The demonstrated effects of Knudsen and surface diffusion apply to obtaining upper limits for all lunar volatile gases. Volatiles such as CO, methane, and CO 2 , which adsorb weakly with the nightside soil, will experience slower outgassing and peak at later local times. Given binding energies and surface mobility, the modeled output is less sensitive to remaining parameters (sticking coefficient, grain size distribution). Only the timeline and magnitude of the gas delays change with the strength of the gas-surface bond, not the general trends described here (higher-than-first-order desorption from granular medium, nonlinear effect of mobility on outgassing from heterogeneous surfaces). For instance, the second-or-higher-order desorption effect under diffusive conditions also appears in the simulation of alkali gas removal from loosely packed beds at Mercuryʼs surface temperatures (Figure 3 of Sarantos & Tsavachidis 2020). Although the strength of that bond, 1.85 eV, is much stronger than the argon-surface bond, 0.28 eV, the granular medium produces the same slow desorption. For these reasons, the effects presented here do not merely describe specific gas-surface systems, but can be generalized to all exosphere-surface interactions involving thermal desorption.