Surface Exospheric Interactions

Gas-surface interactions at the Moon, Mercury and other massive planetary bodies constitute, alongside production and escape, an essential element of the physics of their gravitationally bound exospheres. From condensation and accumulation of exospheric species onto the surface in response to diurnal and seasonal changes of surface temperature, to thermal accommodation, diffusion and ultimate escape of these species from the regolith back into space, surface-interactions have a drastic impact on exospheric composition, structure and dynamics. The study of this interaction at planetary bodies combines exospheric modeling and observations with a consideration of fundamental physics and laboratory experimentation in surface science. With a growing body of earth-based and spacecraft observational data, and a renewed focus on lunar missions and exploration, the connection between the exospheres and surfaces of planetary bodies is an area of active and growing research, with advances being made on problems such as topographical and epiregolith thermal effects on volatile cold trapping, among others. In this paper we review current understanding, latest developments, outstanding issues and future directions on the topic of exosphere-surface interactions at the Moon, Mercury and elsewhere.

Schematic diagram of gas-regolith interactions, showing possible paths taken by exospheric gas species upon arrival at a surface regolith, including sticking to regolith particles, surface diffusion, vaporization, sputtering, and inter-grain hops. These processes together determine the retention, abundance, distribution and diffusion of exospheric gas within the surface regolith, and thereby impact the structure and dynamics of a planetary body's surface bounded exosphere

Introduction
In the gravitationally bound exospheres of massive airless planetary bodies like the Moon and Mercury, atoms and molecules undergo far more collisions with the planet surface than with each other. For this reason gas-surface atomic/molecular interactions, including (1) thermal accommodation, (2) sticking, (3) volatile diffusion within the regolith, (4) cold trapping, and (5) impact and radiation induced diffusion and desorption, have an outsized influence on the distribution, composition, and time variability of these atmospheres. The properties of an exosphere are therefore determined to a much greater degree by the exosphere's interaction with the surface than with itself.
Since the Apollo era the observations of spacecraft and earth-based telescopes have revealed spatial and temporal variations of exospheric gas and surface frost distributions over many orders of magnitude in response to diurnal and, at some bodies, seasonal changes of surface temperature. Examples include the distribution of argon in the lunar atmosphere and on its surface, first investigated by Apollo and later by the NASA Lunar Atmosphere and Dust Environment Explorer (LADEE) spacecraft, e.g., (Benna et al. 2015;Hodges 2018;Hodges and Mahaffy 2016), and the distribution of Sodium in Mercury's exosphere in relation to this planet's surface 'cold poles' (Cassidy et al. 2016) (Sect. 8). Across the solar system, massive airless planetary bodies tilted to the ecliptic, including Saturn's moons Dione and Rhea (Teolis and Waite 2016), Pluto's moon Charon (Grundy et al. 2016;Teolis et al. 2022), and likely the massive Uranian moons -exhibit extreme seasonal changes in their exospheres on much longer timescales, driven by bi-annual condensation and evaporation of exospheric gases frozen at their poles. Yet remarkably, as we discuss in this chapter, even the Moon's argon exosphere, despite its minor 1.5 • angle to the ecliptic, exhibits some seasonal variability likely owing to seasonal condensation and evaporation of frozen argon to and from polar cold traps. In general bodies like the Moon and Mercury appear to represent the limiting case of a weakly seasonal exosphere, in which a day-night adsorption-desorption cycle is the primary effect, while the preponderance of permanent over seasonally shadowed polar cold traps act as ultra-long term exospheric sinks.
Other processes, such as impact and or radiation enhanced diffusion and desorption of semi-volatile species frozen to the surface regolith, and the chemical activation of surface regolith by impacting ions, e.g. from the solar wind, are only just beginning to be explored. Finally, surface microstructure, diffusion of adsorbed gas between regolith grains (Fig. 1), and surface topographical effects, including permanent, seasonal, and diurnal cold trapping, all influence the residence time of exospheric atoms and molecules on the surface material. The total surface interaction/residence time of exospheric species within a granular planetary surface, including the combined effect of surface diffusion and multiple accommodation, adsorption and desorption events between regolith grains (Fig. 1), is a critical parameter controlling the exospheric distribution between the day and night hemispheres, as well as any seasonal variability.
In this chapter we discuss the physics, open questions and current knowledge of three aspects of the surface-exospheric interactions at large airless bodies of the inner solar system: (1) accommodation, sticking and desorption, (2) regolith diffusion, and (3) surface impact/radiation processing and induced diffusion/ejection of adsorbed material.

Thermal Accommodation & Sticking
Energy accommodation is a fundamental and ubiquitous process common to surfacebounded exospheres. Accommodation constitutes the first step in the surface interaction of exospheric atomic or molecular species, and determines the local energy and altitude distribution of the exosphere above the surface. In simplest terms, any surface bounded exosphere consists of both a non-accommodated and an accommodated contribution. The non-accommodated supra-thermal component represents the exospheric source; i.e. neutral species freshly injected into the exosphere by processes such as impact vaporization, photodesorption, sputtering or radiogenic decay/outgassing from the surface, or chargeexchange at altitude. Some atoms from these high-energy processes have sufficient speeds to escape the planet's gravity, or at least to reach high altitudes before falling back to the surface (e.g., Killen et al. 2018). The accommodated component, by contrast, interacts with the planetary surface to acquire an approximate Maxwell-Boltzmann velocity distribution at a temperature close to that of the local surface. The accommodated component has two sources: (1) exospheric neutrals that fall back to the surface under gravity and interact with the regolith, and (2) exospheric source atoms/molecules initially emitted from the surface in an accommodated state; e.g. off-gassing radiolytic or radiogenic species, or sputtered energetic species that impact other regolith grains and lose energy before ever entering the exosphere.
The supra-thermal and accommodated components may segregate into high and low altitude layers, as can be seen in Fig. 2 showing the predicted altitude profile for multiple exospheric species in Europa's exosphere, and for sodium and argon in the Hermean and lunar exospheres, respectively. Compared to the Moon's fully accommodated argon exosphere, energetic release processes at Mercury produce a more extended exospheric sodium altitude profile (Fig. 2) modified to a lesser degree by surface temperature-dependent sticking and accommodation (e.g., Burger et al. 2010;Leblanc and Johnson 2010;Tenishev et al. 2013). At both rocky and icy planetary bodies sufficiently volatile exospheric species may exhibit (1) an extended high-altitude exosphere from freshly sputtered molecules, and (2) a much denser low-altitude thermally accommodated layer (see Fig. 2; Europa's O 2 and H 2 O) with molecules lacking the kinetic energy to attain high altitudes in the exosphere. Hence the exosphere at high-altitude gives a more representative sampling of species immediately after ejection from the local surface below, whereas the low-altitude exosphere may be dominated by species which have interacted and accommodated with the surface many times. Surface gas densities may in fact be greatest at night owing to the colder surface there, which cools the exosphere, and reduces the gas thermal speed and exospheric scale height such as to  Smyth and Marconi (2006). Sufficiently volatile species like O 2 and H 2 O exhibit a low altitude thermally accommodated layer, but at high altitudes the dominance of supra-thermal (non-accommodated) molecules newly produced/ejected from the surface results in a shallow slope of density versus altitude. Right: Lunar Ar and Na gas densities at sub-solar noon, from the models used by Grava et al. (2015) and Tenishev et al. (2013) to fit LACE Ar and remote Na observations. Na density taken from Tenishev et al. (2013); their case 1 that assumes full Na sticking to the surface. Ar density versus altitude is extracted from the Grava et al. (2015) model. Lunar Ar is fully thermally accommodated resulting in a steep slope versus altitude. Lunar Na is freshly ejected largely by photon-stimulated desorption. Therefore Na is not accommodated, producing a shallow slope concentrate the gas there. Examples of this interaction regime are some noble gases at the Moon (Benna et al. 2015). However, the opposite may happen for species such as lunar argon (Sect. 7) that can stick to the surface at the given surface temperature; surface gas densities will be lower at night if sufficiently cold for exospheric atoms or molecules to freeze to the surface.
The degree to which exospheric species accommodate with the surface on a single collision is quantified by the energy accommodation coefficient, which is the ratio (0 to 1) of energy transferred in a surface collision relative to that required for full accommodation. The accommodation coefficient depends on the mass, kinetic energy, and internal degrees of freedom (e.g. molecular rotations and vibrations) of the impacting species, and the magnitude of the interaction potential or binding energy with the surface (Bonfanti and Martinazzo 2016). Another major consideration is regolith porosity, which enhances the 'effective' accommodation of exospheric species undergoing multiple inter-grain collisions within the regolith before returning to the exosphere (Fig. 1).
The angular dependence of the velocity of a gas thermally accommodated to the planetary surface follows a Knudsen cosine-law (Comsa 1968;Knudsen 1909), which implies that the direction in which an atom or molecule rebounds from a surface is independent of the direction at which it approaches. The Knudsen cosine law can be understood in terms of the dispersing micro-geometry of the solid surface (Feres and Yablonsky 2004), but the same decorrelation is expected from thermal accommodation during contact. Armand (1977) determined the desorption velocity distribution based on lattice vibrations, and arrives at the Maxwell-Boltzmann flux distribution, which obeys the Cosine Law.  Kolasinski (2012) The interaction of gas phase molecules with surfaces can involve vastly complex chemical physics and depends on not only the molecule and the surface, but also the presence of co-adsorbates, surface temperature, incoming molecular energy (translational, rovibrational), and the crystal structure (defects, step edges, crystallographic plane, terraces). There are numerous high-level texts with more in-depth coverage of the material (e.g. Billing 2000;Ibach 2006;Kolasinski 2012;Morrison 2013;Zangwill 1988). Molecular beam experiments are reviewed in e.g. Arumainayagam and Madix (1991). The following is meant as only a general introduction to adsorption.
Energy transfer is a pre-requisite for adsorption. For trapping or sticking to occur, the incoming molecule must lose energy, with the molecule finally transferred into a bound adsorbed state. Adsorption is further divided into two categories, activated and non-activated. An easy way to visualize this is to utilize a Lennard Jones potential as shown in Fig. 3. If an impinging molecule has too much energy, the molecule will reflect as in a typical elastic collision (a). A molecule or atom may lose some energy after reflecting off the potential barrier where it then proceeds down well to a chemisorbed state (b). Alternatively, a small barrier may exist such that a slow molecule will trap in the shallow well resulting in a physisorbed state (c). The height of this barrier is what controls the activated or non-activated adsorption process.
At 'low' kinetic energies of order 0.1 eV or below, typical of most simple exospheric species at or below a few km/s speeds, gaseous species will transiently stick to the surface, with ample 'sticking time' (picoseconds or longer) in the surface interaction potential well to fully accommodate before thermally desorbing. At these low kinetic energies (which are albeit still up to thousands of Kelvins) the gas species are fully adsorbed during their interaction with the surface, with near unity accommodation and sticking coefficients. By contrast, sputtered species ejected with 'high' energies just shy of the planetary escape speed, or species accelerated by radiation pressure, may fall back to the planetary surface with eVs of energy, depending on the planet's gravity (or escape speed) and species mass. At eV energies or above the impact of gaseous atomic species with surfaces may be better approximated as binary collisions between the impactor and a single target atom. Backscattered impactors may lose energy to the surface atoms elastically and inelastically in accordance with the back-scattering angle and relative projectile and target masses. Whereas light species (such as H, H 2 or He) may have relatively low accommodation coefficients at eV energies, heavier impactors with mass similar to the target material atoms may exhibit accommodation coefficients closer to unity.
A caveat is that atomic or other reactive species, such as Na at the Moon or Mercury, may have high eV binding energies owing to chemisorption with the surface. Even at eV kinetic energies these reactive species may transiently stick and accommodate with the surface irrespective of their mass. Another consideration for molecular species at eV energies is the absorption of some of the collision impact energy into internal interatomic motions within the impacting molecule, such that the molecule loses translational kinetic energy in the collision with the surface by becoming internally excited.

Adsorption Models
Numerous adsorption models have been developed over the years that have resulted in analytical expressions of the adsorption process that include the energy dependence of accommodation. The first theoretical treatment of energy exchange/accommodation came from Baule (1914) and later by Zwanzig (1960), who treated the surface as harmonic oscillators that are excited by collision allowing the incoming molecules to lose energy elastically. Under these conditions the accommodation coefficient α = 4μ (1+μ) 2 , where μ is the ratio of the incoming projectile mass to that of the surface species μ = mp ms . A variation of this follows the hard cube model, where the incoming molecule is accelerated by an attractive potential (V ) resulting in a simple modification to the Baule equation α = 4μ+V /ε (1+μ) 2 , where ε is the collision energy. These models assume that the surface is fixed at 0 K. Taking into consideration the average thermal distribution of surface velocities results in α = 4μ(ε+V −kT /2)/ε (1+μ) 2 (Bonfanti and Martinazzo 2016). Other classical models of energy accommodation as a function of gas and surface temperature were described by Fan and Manson (2010). For a quantum mechanical description, see texts by Billing (2000) and Zangwill (1988). Accommodation coefficient have often also been measured in the laboratory (e.g., Haynes et al. 1992;Persad and Ward 2016).
The well-known Langmuir model (Langmuir 1918), applicable for sub-monolayer coverage, uses a single binding energy and assumes that adsorption and desorption are at equilibrium and reversible processes without an activation barrier. Using these assumptions, the fractional coverage is derived as θ = Keq p 1+Keq p , where K eq = k ad k des is the ratio of the rate of adsorption to desorption and p is the partial pressure of the adsorbate. As evident from this equation, the sticking coefficient (S) is not dependent on the incoming incident energy, but simply the available surface coverage S ∼ (1 − θ). The Langmuir model was a big leap forward in thinking about gas-surface interactions, but it deviated from experimental results more often than not (Becker and Hartman 1953;King and Wells 1974;Singh-Boparai et al. 1975). Moreover, the Langmuir model does not take into consideration the energy of the incoming particle.
Ideally, if the adsorption process does not have an activation barrier, the fraction of atoms or molecules that adsorb to the surface will decrease with increasing incoming energy. Here, the sticking probability for a non-activated process will be unity until the incident energy reaches a critical value, at which point it rapidly diminishes to zero. In reality, this apparent step function is smoothed out by thermal effects from the surface. The sticking probability for an ideal surface without an activation barrier can be modelled treating the surface atoms as hard cubes (Sipkens and Daun 2017). The cube has an attractive square well which accelerates the molecule as it approaches the surface to a velocity u n . As an approximation only the cube's motion in the direction normal to the adsorption surface is considered, with the velocity dictated by the surface temperature T . The normal component of the particle velocity v after collision with the hard cube has a modified Maxwell-Boltzmann distribution,  Bowker (2016), showing processes of physisorption, chemisorption and diffusion over occupied (extrinsic) and unoccupied (intrinsic) adsorption sites. P desi and P dese denote the desorption probabilities from intrinsic and extrinsic sites, while P ai and P diff e are the probabilities of adsorption and of surface migration/diffusion to an intrinsic site 2kT v 2 ) dv. Sticking occurs if the normal velocity of the molecule after collision falls below a critical value which cannot escape the well. The sticking coefficient is derived from the total fraction of molecules below this critical velocity, v c , and can be analytically expressed as S = 1 + erf(bv c ) + e (−b 2 v 2 ) bunπ 1/2 , where b 2 = m 2kT . In the case of an activated adsorption process, the system will exhibit low sticking coefficients at low energies and low coverages. In this instance, a low energy molecule will reflect off the barrier and not adsorb. Here, the sticking coefficient is practically the inverse of the non-activated case and will be zero until it approaches the barrier height where it should approach unity and then drop off again as it approaches some critical value in energy where it approaches the elastic regime. Quantum tunneling and thermal effects will smooth this step function.
A follow-up to the Langmuir model was established by Kisliuk (1957), whose model allowed for the understanding of atypically high adsorption rates at high surface coverage (Kisliuk 1957). This model is best known as the precursor adsorption model. The theory behind the precursor adsorption model is relatively simple and is outlined in Fig. 4. Here, there are two types of precursors, intrinsic (adsorbate over an empty site allowing chemisorption) and extrinsic (adsorbate physisorbed over a filled site). The intrinsic molecules have a probability of adsorbing, desorbing, or migrating to an empty site, while the extrinsic molecules have a probability of desorbing or migrating. The sum of these probabilities yields an analytical expression of the sticking coefficient, S = S 0 (1 + θ 1−θ K) −1 , where K is the precursor parameter and is given by K = S 0 P des P ads with S 0 being the sticking constant at zero coverage, P des the desorption probability for the extrinsic molecule, and P ads the probability for intrinsic chemisorption; both of which are products of the general Arrhenius equation with activation energy (E), prefactor f , and the average residence time spent at a particular site (τ ) i.e. P = f e ( −E kT ) τ . Other widely used adsorption models are the Freundlich isotherm (with multiple binding energies), the Brunauer-Emmett-Teller (BET) isotherm (for multilayer adsorption), and the Zeta Adsorption isotherm (Narayanaswamy and Ward 2020).

Volatile Diffusion in the Regolith
Gas atoms and molecules can travel in the voids between grains and penetrate centimeters into the ground (Reiss et al. 2021). Furthermore, adsorbed atoms can also travel along a grain and onto other grains following a slower diffusion mechanism called surface diffusion ( Fig. 1). Although slow, surface diffusion can be the dominant mechanism of diffusion for species that bind strongly with the surface such as sodium (Sarantos and Tsavachidis 2020). We consider as an example the case of water vapor in the following discussion, to elucidate the physics of volatile migration and diffusion in a porous regolith.

Surface Bonding and Vapor Migration
Water molecules and other volatile species can migrate in the porous regolith via a repeated sequence of jumps. At room temperature the diffusion time is approximately the sum of times for each molecular flight, but at low temperature, migrating molecules spend most of their time residing on grain surfaces. In either case, the migration process can be thought of as a random walk, and it is thus described by a diffusion equation. The main uncertainty in quantifying the mass flux lies with the surface-vapor interactions.
At low temperature H 2 O can be in crystalline or amorphous form, or adsorbed at monolayer or sub-monolayer thicknesses onto a substrate. All naturally occurring ice and snow on Earth's surface is crystalline with a hexagonal lattice structure. Although many other crystal structures are known for water ice, few are stable at low pressure. Amorphous ice is solid (condensed) H 2 O that lacks a crystal structure, and forms at very low temperatures (<∼ 140 K at laboratory time scales; Sack and Baragiola 1993). It takes about a dozen monolayers of H 2 O until the properties of macroscopic ice are reached (Cadenhead and Stetter 1974). Lunar grains have a high specific surface area on the order of 500 m 2 /kg as indicated by measurements of grain surface areas in lunar soil samples (Heiken et al. 1991). The areal density of H 2 O molecules for a monolayer of solid ice is ϑ = (ρ/μ) 2/3 ≈ 10 19 molecules/m 2 , where ρ is the density of ice and μ the mass of a water molecule. Multiplied, a monolayer of H 2 O corresponds to ∼ 5 × 10 21 molecules adsorbed per kg of lunar regolith, which equates to about a hundred parts per million by mass. Whereas most of the lunar surface is extremely desiccated, much higher water concentrations have been measured at some locations in the mid-latitudes (Honniball et al. 2021) and lunar polar regions (Colaprete et al. 2010;Feldman et al. 2001).
Vapor pressure, sublimation, and molecular residence times are closely related. The sublimation rate of ice in vacuum F ice is given by the Hertz-Knudsen formula (Persad and Ward 2016;Watson et al. 1961b) with p v the equilibrium vapor pressure, k B the Boltzmann constant, T temperature, and μ the mass of an H 2 O molecule. When the vapor is in equilibrium with ice, the evaporation rate equals the condensation rate, and the above formula is the flux of a dilute ideal gas at the equilibrium vapor pressure. The flux of water molecules from the ice surface is the same with or without incoming water molecules, so this equation describes the sublimation rate into vacuum. A dependence F = F (θ, T ) of the sublimation rate on the surface coverage θ is also typical for sub-monolayer coverage. Conversely, when θ is much more than a monolayer, F = F ice (T ), and the equilibrium vapor pressure of a solid is of the form, where H subl is the sublimation enthalpy (the heat of sublimation). The average sticking time τ of an atom or molecule on the substrate surface can also be derived from the inverse of the Polanyi-Wigner equation where ν corresponds to the vibrational frequency of the bond between the adsorbate and the substrate surface, typically 10 12 -10 14 Hz, and Q is the energy of adsorption (Adamson 1982;Atkins 1986). Not all incident molecules stick to the surface. A fraction (1 − S) of incident molecules are reflected from the surface (Persad and Ward 2016). For ice, the sticking coefficient at 40-120 K temperatures is in the range 1-0.7 (Haynes et al. 1992). The binding energy of molecules adsorbed to the regolith or on amorphous ice differs from that in crystalline ice (Speedy et al. 1996). The sublimation enthalpy of crystalline ice is 51 kJ/mol, while adsorbed water has a 60-180 kJ/mol range of desorption energies on lunar regolith simulants (Hibbitts et al. 2011;Poston et al. 2013) and lunar samples (Jones et al. 2020;Poston et al. 2015) at sub-monolayer abundances. de Leeuw et al. (2000) report adsorption energies for H 2 O on forsterite surfaces of 100-172 kJ/mol. Other atoms such as sodium and potassium bind much more strongly with the silicate surface, with binding energies in excess of 174 kJ/mol (or 1.8 eV) (Yakshinskiy et al. 2000). Irradiation can increase the effective degree of bonding by creating defects and vacancies, while diffusion in powders may also increase the effective binding energy by slowing desorption (Sarantos and Tsavachidis 2020). Finally, surface microstructure prolongs the residence time of adsorbates on the surface by reducing the desorption yields via readsorption on adjacent grains (e.g., Cassidy and Johnson 2005; Sarantos and Tsavachidis 2020). Figure 5 shows the residence time of water molecules (harmonic mean) as a function of temperature. Below roughly 250 K, the molecules spend most of their time residing on the surface and a comparatively short time in flight (assuming a 1 mm inter-grain flight distance). Additionally, the likely mobility of water adsorbates on surfaces with a distribution of desorption energies further changes the residence time by enabling motion between low and high energy sites. Sarantos and Tsavachidis (2021) have recently suggested that surface diffusion affects the desorption rate in a nonlinear manner, suppressing desorption at lower temperatures yet enhancing desorption at higher temperatures.
Similar conclusions can be drawn for less volatile elements like sodium, for which it can be shown that the residence time of atoms against thermal desorption is ∼ 1000 years at maximum lunar surface temperatures, whereas on average several days elapse between desorption events for the more efficient photon-stimulated desorption mechanism (Sarantos and Tsavachidis 2020 and references therein). These timescales are much shorter at Mercury due to the higher surface temperatures. The surface residence times, along with the escape ratio, control the evolution of the atmosphere following episodic brightening events such as meteor showers (Colaprete et al. 2016).

Models of Subsurface Migration
Mathematically, the migration process can be described at three levels: • Molecular hops (random walk; discrete formulation) • Diffusion equation (continuum formulation; partial differential equation) • Boundary-value problem (for stationary solutions) The first type of model is a statistical simulation of individual molecules (random walk). Molecular random walk leads to diffusive migration where molecules reside on grain surfaces and then hop a distance , which is the molecular free inter-grain path within the regolith (Schorghofer and Taylor 2007). The mean grain size in lunar soil samples is typically 45-100 µm (Heiken et al. 1991), which may be taken (Schorghofer and Taylor 2007) as a crude approximation of .
The migration process can also be described by a continuum equation. The net flux J of molecules is given by differences in sublimation rates between two grain surfaces, and therefore where t and z are time and depth. The gradient in sublimation or desorption rate F is due to gradients in T or θ . Transport can be caused by differences in surface concentrations or by differences in temperature. At constant temperature, only concentration-driven migration occurs. Note that the time average of J , the net flux, is given by the gradient of the time average of F .
For situations that are stationary (time-independent, periodic, or quasi-periodic) the boundary-value formulation leads to insightful results. An example is the determination of the loss rate of ice buried by a layer of regolith of thickness z. Mass conservation dictates that after a transient period, the flux J is constant between an ice table and the surface. Thereafter, the value of J is determined by the difference between F on the surface (which is essentially zero, because θ is almost zero) and F ice at the ice-regolith boundary. The loss rate is given by The flux is reduced relative to sublimation from an exposed ice surface by a factor of the order / z. This equation is equivalent to the Knudsen flow through a porous medium. The Knudsen diffusion coefficient is proportional to the mean free inter-grain path , and the vapor density is related to p v , and therefore to F ice , via the ideal gas law. Schorghofer and Taylor (2007) obtained two more solutions to the vapor migration transport equations for constant temperature. For slow continuous water delivery to the surface, there is an optimum temperature for which the amount of adsorbed subsurface water is at a maximum. At low temperature little H 2 O accumulates, because the migration is too slow. At temperatures exceeding the optimum the loss from the surface is so fast that few source molecules reside on the surface at any time. Similarly, for an ice cover at constant temperature, there is an optimum temperature which maximizes the subsurface H 2 O content. In both cases, the grain coverage is limited to about a monolayer of water, because, at constant temperature, differences in surface concentration are necessary for net transport of additional molecules. For time-varying temperature, subsurface vapor migration can lead to a "pumping effect". Transport of H 2 O molecules can result from differences in sublimation rates with depth, even without change of mean temperature with depth. This phenomenon is described in more detail in Sect. 5.
The effect of inward migration on desorption of other adsorbates was studied with a three-dimensional model by Sarantos and Tsavachidis (2020). They used random walk simulations of desorption inside a granular medium consisting of spherical particles of different size. Using kinetic parameters that represent alkali, argon, and water adsorbates, it was demonstrated that surface and Knudsen diffusion, processes competitive to desorption, reduce the desorption rates from a porous medium beyond the rate reduction due to re-adsorption. Sarantos and Tsavachidis (2021) pointed out that thermal desorption from a granular medium is a decelerating process because it initiates Knudsen diffusion which in turn slows down desorption. These authors proposed that thermal desorption of adsorbates from a powder is not a first-order (i.e. constant rate) process, but is instead better approximated by a square or even higher order dependence of the desorption rate on the amount of adsorbate when diffusion is considered. This finding means that models consistently underestimate how long it takes for regolith to outgas volatiles such as CH 4 , CO and CO 2 adsorbed weakly on the lunar night side because Knudsen diffusion is not properly accounted for. It was also demonstrated by Sarantos and Tsavachidis (2020) that diffusion leads to a non-linear and non-monotonic relationship between the surface temperature and alkali photodesorption rates from regolith.
It should be noted that bulk diffusion, i.e., diffusion within a solid structure, controls the source of fresh atoms for many of the atmospheric species found around Mercury and the Moon. For example, Killen et al. (2004) have described how diffusion of sodium in the grain rims controls the release rate for photon-stimulated desorption, thermal desorption and ion sputtering, processes that, unlike meteoroid impacts, can extract atoms only from the top nm of grains. Continuous diffusion of argon through the pore space between rocks in the top 25 km of lunar crust has been considered as a source for the lunar argon atmosphere (Killen 2002). And, finally, bulk diffusion followed by re-combinative desorption is responsible for converting the solar wind protons into the lunar H 2 atmosphere (e.g., Tucker et al. 2019).

Cold Trapping
At low temperature, sublimation rates into vacuum become negligible. Figure 6 summarizes sublimation rates for a variety of volatiles (Zhang and Paige 2009). For example, the H 2 O loss rate at 115 K is around 0.1 m/Gyr, and the ice is said to be "stable" over geologic time periods. This is the phenomenon of cold trapping. Among volatiles that are abundant in comets (e.g. H 2 O, CO 2 ) or volcanic outgassing (CO 2 , SO 2 , H 2 O) the one with the lowest sublimation rate is H 2 O (Watson et al. 1961b). It has long been realized that craters near the poles of Mercury and the Moon can be permanently shadowed and volatiles could have accumulated in these eternally cold regions (Thomas 1974;Urey 1952). Recently, cold traps have also been predicted and then observed on Ceres (Hayne and Aharonson 2015;Platz et al. 2016;Schorghofer et al. 2016).
The concept of "stability" (negligible sublimation rate) can be extended from surface cold traps to the subsurface. However, whereas surface cold traps can be supplied with water directly through an exosphere, processes that emplace ice at depth are rarer, e.g., burial by impact ejecta. Maps of subsurface ice stability of the lunar polar regions have been published by Paige et al. (2010) and Schorghofer and Williams (2020).
Another form of trapping for less volatile elements is provided by the surface microstructure. Sarantos and Tsavachidis (2020) demonstrated that at lunar temperatures, conditions Summary of sublimation loss rates into vacuum for various species. Adapted from Zhang and Paige (2009) for which thermal desorption of alkalis has a mean lifetime of thousands of years, about half the adsorbates never participate in photodesorption because they hide in shadows the size of the lunar grains (i.e. the underside of grains). This is equivalent to assuming in models that every atom has a 50% probability of reaction between desorption events. On the other hand, if surface diffusion is fast due to mobility of these adsorbates along the grain surface, all atoms can dislodge from the microshadows and participate in desorption but at reduced rates.

Vapor Pumping by Temperature Cycles
Under suitable conditions, H 2 O molecules can be "pumped down" into the regolith by periodic (day-night) temperature cycles, leading to an enrichment of H 2 O in excess of the surface concentration. The amplitude of temperature oscillations quickly decays with depth, and the strongly nonlinear dependence of molecular residence times on temperature leads to vertical drift processes, see Fig. 7. The concept of an ice pump (i.e., vapor pumping that leads to the accumulation of macroscopic quantities of ice in the subsurface) was first developed in the context of Mars (Mellon and Jakosky 1993), where pumping can occur from a humid atmosphere into a porous subsurface. For the tenuous atmosphere surrounding the Moon, volatile H 2 O molecules on the surface are the source of the pump. If the surface concentration is high enough and the temperature amplitude is significant, H 2 O will be sequestered (Schorghofer and Aharonson 2014).
The simplest form of an ice pump involves pumping from an ice cover. Crucial for the pumping mechanism is the damping of the temperature amplitude with depth. Since the vapor pressure at the ice is the saturation pressure, which has a convex shaped temperature dependence, larger temperature amplitude implies a larger vapor pressure, thus creating a gradient in vapor pressures that preferentially moves molecules downward. If the concentration is sufficiently high that the nonlinear temperature amplitude effect on the top surface compensates for the reduced saturation pressure of ice at depth, pumping occurs. A weaker form of the lunar ice pump occurs when it does not produce ice but only excess adsorbate.
If the time averaged sublimation rate of adsorbed water on the surface is sufficiently high, it can balance that of pure ice at depth. A pumping differential may be defined as ) where the "ice table" corresponds to the shallowest depth with ice. If F > 0, then downward pumping exceeds the upward loss. If F < 0, the pumping is too weak and any subsurface ice experiences net loss.
The surface temperature and H 2 O surface concentration vary over time. Water molecules can accumulate at night (by delivery from an exosphere), but are lost when the surface warms Fig. 7 Schematic illustration of a subsurface temperature profile (solid line = instantaneous, dashed lines = minimum and maximum). Any volatile water molecule has a certain probability to hop up or down. A water molecule on the surface (upper dot) has a different mobility than a molecule at depth (lower dot) where temperature and temperature amplitude are different. In the long term, this leads to a net vertical flux of water molecules. When that flux is downward, this acts as an "ice pump". The nonlinear dependence of the sublimation rate on temperature is also illustrated (gray line). From Schorghofer and Aharonson (2014) in the morning. Quantifying the strength of an ice pump requires a quantitative model of the population of water molecules on the surface over time. Schorghofer and Aharonson (2014) have carried out model calculations for the surface population of adsorbed water molecules on the Moon. In their model, the temperature varies sinusoidally with time for half a solar day, which mimics daytime, and it is constant for the other half of the solar day, which mimics nighttime. Figure 8 shows the pumping differential as a function of mean and peak temperature. This phase diagram can be divided into three nearly complementary regions. At peak temperatures below about 120 K there is very weak pumping, and this parameter region nearly coincides with the temperature conditions for classical cold trapping. At peak temperatures larger than about 120 K and mean temperatures lower than about 105 K, significant pumping occurs. When the mean temperature is above ≈ 110 K, neither pumping nor cold trapping occurs. Classical cold trapping and strong pumping are nearly complementary. Areas of strong pumping exhibit large variations in surface water concentration over one lunation.
On bodies without an atmosphere, pumping is an inefficient process due to the rapid loss of surface molecules to space relative to their downward diffusion. For every molecular hop from the surface downward there is probabilistically one hop upward and thus one molecule lost. The column-integrated subsurface ice density is smaller than the time-integrated supply of water by at least a factor of / z 1. A rocky surface layer with large pore spaces (large ) favors fast diffusion, whereas dust (small ) is an inefficient medium for pumping. At most a few percent of the H 2 O delivered to the surface could have accumulated in the near-surface layer. Fig. 8 Pumping differential (defined in the text) color coded, according to model calculations for the Moon that assume a supply rate of 1 m/Ga, a space weathering rate of 1 m/Ga, and prescribed temperature variations. Below the dash line, maximum temperature is lower than 120 K, corresponding to the classical definition of a cold trap. There are three nearly-complementary regions: weak pumping and classical cold trapping (red), strong pumping (blue), and net loss (grey). From Schorghofer and Aharonson (2014) Maps of the pumping differential have been published by Schorghofer and Aharonson (2014) and Schorghofer and Williams (2020). The total area where substantial pumping occurs (a pumping differential larger than half the supply rate) is estimated to be more than five times the area of surface cold traps. Typically pumping occurs on pole facing slopes in polar areas, but within a few degrees of each pole the equator facing slopes are preferred.

Epiregolith Thermal Gradients
The surface layer on airless bodies is heated via solar radiation. To first order the distance from the sun, latitude, and albedo determine the average temperature at depth. In more detail, the thermophysical properties of the surface (e.g. thermal inertia, density, particle size, packing density) and rotational rate determine the depth at which this equilibrium is achieved and the magnitude of the thermal gradients that develop between the surface and subsurface. This 'surface thermal gradient' is time-varying, with the subsurface warmer than the surface at night and cooler than the surface during the day. In addition, for airless body surfaces dominated by fine particulates, the upper few particles radiate disproportionally to space relative to deeper particles and create an additional, persistent 'epiregolith thermal gradient' within the uppermost few 100 µm at times of day when the near-surface is sufficiently warm that radiation (rather than conduction) is the dominant heat transfer mechanism (e.g. Hale and Hapke 2002;Jakosky 1994, 1997;Logan and Hunt 1970;Logan et al. 1973). The surface and epiregolith thermophysical properties influence the stability of volatile reservoirs and the time-varying nature of exosphere transport.

Thermophysical Properties: Moon, Mercury & Asteroids
Studies of lunar surface temperatures and physical properties benefit from an unprecedented amount of data compared to other airless bodies, including telescopic, orbital, in situ, and analyses of returned samples. Properties on small scales were determined from laboratory measurements and Apollo drill cores, and showed the lunar regolith to have a highly insulating upper layer and a more conducting lower layer (e.g. Carrier et al. 1991;Jones et al. 1975;Langseth et al. 1976). Analytical and numerical models based on these characteristics predict the lunar surface to be extremely insulating and capable of maintaining temperature Fig. 9 Recent thermal models of the lunar surface have been revised using orbital observations from the Diviner lunar radiometer and use continuously increasing density and conductivity rather than one or two distinct layers. The new models (dotted line) are a better fit for Diviner nighttime data than a single layer (dashed line) or two-layer (solid line) models. From Vasavada et al. (2012) extremes on short distances to the extent that shadowed regions near the lunar poles could harbor water-ice deposits (e.g. Vasavada et al. 1999;Watson et al. 1961a).
High resolution temperature maps are provided by the Lunar Reconnaissance Orbiter Diviner Lunar Radiometer and confirmed the extreme nature of the lunar thermal environment Vasavada et al. 2012;Williams et al. 2019Williams et al. , 2017. Diviner also provides inferred physical properties of the surface including rock abundance (Bandfield et al. 2011), surface roughness (Bandfield et al. 2015), and regolith thermal inertia (Hayne et al. 2017). Critical to the interpretation of this dataset was a revision to the predominant two-layer lunar thermal model (e.g. Vasavada et al. 1999) to a continuously varying density and conductivity model (Hayne et al. 2017;Vasavada et al. 2012) that more accurately represents the regolith structure over depth scales of order 1 m (Fig. 9).
As described above, the thermal and diffusive properties of the lunar regolith afford significant protection to buried volatile reserves and thus increase their stability on geologic time scales. The same models used to derive surface properties from Diviner temperature observations and lunar topography can be used to predict subsurface temperatures (e.g. Hayne et al. 2017;Paige et al. 2010). In addition, sub-millimeter instruments (e.g. Chang'e 1 and 2 Microwave Radiometer) offer a direct temperature measurement of depths up to ∼ 2 m. These data predict an area ∼ 10× larger where ice is stable in the near subsurface than would be stable on the surface (Fig. 10)  . A diurnal or seasonal connection between this potential reservoir and any putative surface reservoir or exosphere has not Fig. 10 Near surface ice stability is modeled to be possible over vastly larger areas than strict permanently shadowed regions (white areas). The colored regions represent different depths to reach ice stability in the top 1 meter of regolith. From Paige et al. (2010) been established. However, the subsurface reservoir can be accessed via impacts, as was demonstrated by the LCROSS mission (e.g. Colaprete et al. 2010;Hayne et al. 2010).
The surface thermophysical properties of Mercury are currently poorly constrained. As such, the estimation of subsurface temperatures is uncertain. Numerical thermal modeling of Mercury was used to investigate temperatures associated with bright and dark regions suspected to be associated with surface and subsurface ice deposits Paige et al. 2013). However, these models assumed modified lunar thermophysical parameters as inputs. The Mercury Radiometer and Thermal infrared Imaging Spectrometer (MERTIS) instrument on BepiColombo will directly measure global surface temperatures at multiple times of day and allow accurate constraints on Mercury surface thermophysical parameters for the first time (Hiesinger and Helbert 2010).
Asteroids can differ significantly from the Moon with regards to regolith formation, age and gardening and thus a wide range of surface thermophysical properties are expected on different types of asteroids. Inferred thermal inertias of asteroids range from < 20 to 700 Jm −2 K −1 s −1/2 (Hanuš et al. 2018).

Temperature Gradients Inside the Epiregolith
Early laboratory measurements of returned Apollo soils showed significantly different emission behavior when samples were measured in typical laboratory conditions and under a simulated lunar environment (Logan and Hunt 1970;Logan et al. 1973). Since those first measurements, experiments have grown more sophisticated and spectral effects associated with varying sample composition, albedo, temperature, porosity, and roughness have been investigated (e.g. Donaldson Hanna et al. 2016;Henderson et al. 1996). Several environmental and sample characteristics explain these discrepancies. First, fine particulates under Fig. 11 Numerical models of the epiregolith consistently show thermal gradients but differ in their predicted amplitude. The magnitude of the thermal gradient also varies as a function of solar incidence angle and/or local time. Top: Noontime equatorial illumination of basalt for different regolith grain sizes from Henderson and Jakosky (1997). Bottom: Noontime and predawn equatorial illumination of a simulated regolith from Hale and Hapke (2002). Here kE T has units of J/m 2 /s/K, with k the solid state thermal conductivity and E T the thermal radiation extinction coefficient vacuum are extremely insulating. Second, when exposed to a cold shroud the surface particles rapidly cool off. Third, a broadband lamp deposits heat to depths greater than a few particles. And fourth, emission in the thermal infrared is depth-dependent on wavelength. In total, these simulated lunar conditions produce an emission spectrum where some radiation comes from the warmer interior (at the Christiansen feature emissivity maximum) and some radiation comes from the cold exterior (at the Reststrahlen Band absorption bands). Diviner orbital measurements have confirmed the spectral effects of these thermal gradients on the lunar surface (e.g. Donaldson Hanna et al. 2016;Greenhagen et al. 2010).
Efforts to model the epiregolith thermal gradients have demonstrated the underlying physics, including a dependence on scattering properties, particle size, and packing density. There remains some uncertainty as to the actual magnitude of near-surface thermal gradients on the Moon, Mercury and other airless bodies. Using similar but different computational methods, Henderson and Jakosky (1997) predict gradients of 40-50 K per 100 µm under lunar-like conditions, whereas Hale and Hapke (2002) predict gradients closer to ∼ 10 K per 100 µm (Fig. 11). Comparing laboratory measurements of quartz under simulated lunar conditions to models, Millán et al. (2011) infer thermal gradients similar to those calculated by Henderson and Jakosky (1997). Epiregolith temperature varies with time of day (Fig. 11, bottom) and latitude. Models indicate that even during the day, the lunar epiregolith is cooler at the surface than at depth, owing to the (wavelength-dependent) optical transmis-sion of thermal infrared radiation through the upper few 100 µm of regolith to space. As a result, temperature rises relatively steeply within the topmost hundreds of microns of regolith before gradually decreasing to approach the diurnal average temperature at ∼ 1 m depth (Fig. 11). At cooler, night-time temperatures, radiative cooling is less important and epiregolith temperature is nearly constant with depth.
The ability of Mercury and asteroids to produce epiregolith thermal gradients of a magnitude large enough to significantly affect the exosphere is unknown. Certainly, some asteroids have coarse particle regolith that would inhibit the formation of a thermal gradient in the epiregolith. In addition, even fine-grained bodies further from the sun would experience lower maximum temperatures where the magnitude of the epiregolith thermal gradient is largest and potentially reduce the effects of this phenomenon. Mercury, on the other hand, is expected to have fine-grained regolith, rotates slowly, and is near to the sun. However, the albedo of Mercury is darker than the Moon and this is expected to reduce epiregolith thermal gradient development.

Epiregolith Effects on Surface Reservoirs
Epiregolith thermal gradients on the Moon, Mercury or asteroids have the potential to be significant drivers of exospheric transport. Remote sensing reports in the UV and NIR often consider surface temperature when discussing the potential for diurnally varying (or invariant) surface hydration (e.g. Hendrix et al. 2019;Li and Milliken 2017). Variations of temperature with depth have been used to estimate the "pumping" of water vapor into the regolith as described above. Additionally, an increase of temperature with depth in the first mm of soil has been shown in simulations to keep sodium adsorbates near the top of the surface, to decrease inward diffusion, and to increase desorption flux (Sarantos and Tsavachidis 2020).
However, more important than a generalized 'surface' temperature, is the distribution of temperatures of the actual grains from which the light is scattered. Due to epiregolith thermal gradients, these temperatures must be lower than the temperatures predicted by typical thermal models (e.g. Hayne et al. 2017;Vasavada et al. 2012) that do not resolve the epiregolith. Observations of thermal emission at short wavelengths (∼ 3 microns) may also be biased towards higher temperatures due to surface anisothermality effects (e.g. Bandfield et al. 2018). As a result, a wide range of different 'thermal corrections' to these types of data have been presented in literature (e.g. Bandfield et al. 2018;Li and Milliken 2017;Wöhler et al. 2017). The significance of this variance on exosphere transport and surface volatile reservoirs is currently unknown and represents an area of ongoing research.

Lunar Exosphere
The Moon's 40 Ar exosphere, first studied over a half century ago during the Apollo missions as one of the first species detected in the lunar exosphere, is one of the oldest and most spectacular examples of a planetary atmosphere dominated by adsorption and desorption to and from the surface. The effect is best illustrated in Fig. 12, which shows the diurnal profile of exospheric surface 40 Ar gas density as a function of time measured during two different lunations by the surface mass spectrometer LACE (Hoffman et al. 1973). The argon gas density decreases from dusk to dawn as a result of the adsorption into the progressively colder nighttime surface (dayside LACE measurements were deemed unreliable). This implies that the dominant regolith interaction is adsorption-desorption (Hodges and Mahaffy 2016). The Fig. 12 40 Ar exospheric surface densities over a lunar day, measured by LACE from dusk to dawn over two different lunations. The density declines over night as 40 Ar increasingly sticks to the cooling surface, before rapidly rising near dawn as sunrise prompts 40 Ar to desorb back into the exosphere. Also shown is an exospheric model fit (Grava et al. 2015) that assumes a pre-existing 40 Ar exosphere augmented by a prompt 40 Ar release before the 24 Mar-7 Apr lunation. From Grava et al. (2015), data from Hodges (1975) dawn surge in the exospheric density (Fig. 12) is consistent with desorption of 40 Ar as the cold night surface is rapidly heated by the rising sun.
The shape of the density decline during the night is controlled by the strength of the interaction of the argon atoms with the surface. The residence time (time spent by an argon atom on a surface grain) has an exponential dependence on its heat of adsorption Q (Sect. 3.1), which is a still poorly constrained value that determines the shape of the density decline at night. The exospheric simulations performed by Hodges (1980) to explain the argon density decline revealed that the heat of adsorption was 0.26 eV (Grava et al. 2015 derived a similar value of 0.28 eV). This value is much higher than the value derived by adsorption experiments of 40 Ar on glass (0.16 eV; DeBoer 1968). The very high value of Q for 40 Ar inferred from LACE may be a radiation effect; e.g. trapping of 40 Ar atoms in radiation defects. Hodges (1980Hodges ( , 2018 suggested that 40 Ar's high heat of adsorption may be a result of the exceptional cleanliness of the lunar soil grains unattainable in laboratory, due to the removal of water molecules from high-energy adsorption sites by prolonged exposure to solar wind or serpentinization (sequestration of water molecules by reaction with olivine).
According to Hodges and Mahaffy (2016) this adsorption of argon on cold lunar surfaces enables 40 Ar to be trapped in Permanently Shaded Regions (PSRs) which never receive direct sunlight, or seasonal cold traps shielded from direct sunlight only during winter. Regular seasonal deviations of the exospheric 40 Ar abundance from a steady state exospheric model were observed by LADEE/NMS as shown in Fig. 13, consistent with adsorption of argon in seasonal cold traps. The 1.5 • tilt of Moon's rotation axis gives rise to "seasons" in the argon atmosphere of the Moon, where argon atoms freeze in the seasonal cold traps near one pole, and then desorb back into the exosphere later around lunar equinoxes. In Fig. 13, the peak exospheric density measured by LADEE was in January 2014, about 1 month after the lunar vernal equinox (in December), when the seasonal cold traps at the north pole begin to release Fig. 13 The deviation of argon exospheric surface density measured from LADEE/NMS (circles with 1-σ error bars) from a steady state model based on thermal desorption (gray line). From Hodges and Mahaffy (2016) the adsorbed atoms. This delay may be accounted for by variations in transport time through the exosphere and in thermal inertia (Hodges and Mahaffy 2016). Similar behavior may be expected for other condensable species such as H 2 O. It's also possible that sub-monolayer deposits of argon exist together with other condensable species on the lunar PSRs perennially, making the measurement of 40 Ar valuable for inferring the possible behavior of other condensable species such as water that are more difficult to measure.
A dawn surge (Fig. 12) was also detected in the exospheric methane density (Hodges 2016) by LADEE/NMS, indicating that this species also condenses to the cold lunar night side. However the CH 4 sunrise density maximum is shifted from dawn by almost 1 hour into the day (7:00 am local time). Argon, too, has a delayed sunrise density maximum, although the delay is shorter (at 6:40 am). The delay was interpreted by Kegerreis et al. (2017) to be due to temporary sequestration of argon atoms that migrate deep into the regolith during the night. This argon reservoir is at first unable to re-migrate upwards after sunrise as the sub-surface remains for a time colder than the surface. As a result, the argon atoms are released later, possibly as late as mid-morning. Evidence for this interpretation was obtained by Sarantos and Tsavachidis (2021) with a three-dimensional simulation (Sects. 3.1-3.2) of argon transport inside piles made of spherical grains of porosity ∼ 0.55. They found that inward diffusion does not occur at night, but rather as the regolith warms at dawn, with no temperature gradient required. If that is the case, then other condensable species can display the same behavior. Figure 14 shows the exospheric density measured at the lunar surface by LACE for several species: 40 Ar and 36 Ar, CO 2 , and CO or N 2 . Of these, only 40 Ar exhibits a clear sunrise exospheric density maximum. The lack of a sunrise density maximum at mass 28 is understandable for CO, which is less condensable than 40 Ar (Zhang and Paige 2009), and for N 2 , which is not expected to condense. However the lack of a sunrise surge for mass 44, if it is CO 2 , is less understandable, since this less volatile species may remain adsorbed up to substantially higher temperatures than 40 Ar. More measurements from orbit or from the lunar surface are needed to understand this discrepancy and to reveal the true identity of these detections.
Examples of the type of variation expected from gases that only weakly interact with the surface are provided by the detections of lunar He and Ne, which are illustrated in Fig. 15. Instead of the Ar vapor condensation that follows the night time cooling of the surface from dusk to dawn, the surface exospheric gas densities for the more volatile He and Ne atoms are Fig. 14 Surface gas densities measured near lunar sunrise from LACE mass channels 28 (N 2 or CO 2 ), 36 ( 36 Ar), 40 ( 40 Ar), and 44 (CO 2 ). "T" is the time at which the dawn terminator crossed the landing site, while "S" is the actual sunrise time at the landing site. The 40 Ar density rises as dawn approaches, as exospheric 40 Ar atoms "hop" across the terminator from the dayside and then freeze to the colder night side surface. The lack of a similar pre-dawn rise in CO 2 is surprising, as CO 2 is expected to also stick at the Moon's night side surface temperatures. However the 36 Ar signal is heavily contaminated by instrumental HCl, and all four mass channels rise on the dayside owing to outgassing of hydrocarbon contaminants (Hoffman et al. 1973). Therefore, more measurements are needed to verify LACE's CO 2 detection and the contribution of instrumental effects to the data.  . 15 Surface gas densities of exospheric species that due to weak gas-surface interactions desorb even at the coldest night side temperatures of the Moon and exhibit peak densities at night. From Benna et al. (2015) most dense on the night side, showing an anti-correlation to surface temperature T (which decreases continuously during the night due to thermal inertia), and follow approximately a law n ∝ T −5/2 at flux balance (Hodges and Johnson 1968).  Schleicher et al. (2004). North is up, dawn is to the left. Right: Simulated Na tangential column density with a model of successive sodium adsorption-desorption cycles, from Mura et al. (2009)

Mercury's Exosphere
Although numerous species have been observed in Mercury's exosphere as reviewed in the companion chapters by Grava et al. and Leblanc et al., sodium is a species that exemplifies circulation as a manifestation of gas-surface interactions discussed in this chapter. Mercury's sodium exosphere has been observed since its discovery in 1985 (Potter and Morgan 1985) and is produced by photon-stimulated desorption, and thermal and impact vaporization (Killen et al. 2007), with a relatively minor contribution from ion sputtering (McGrath et al. 1986). Since this first detection, several further observations from ground telescopes were performed, with Mercury at different True Anomaly angles and phase angle, with the aim of determining the temporal (i.e. function of True Anomaly Angle) and spatial (i.e. function of local time) variability of exospheric sodium. These observations are summarized in the companion paper by Leblanc et al. Almost twenty years after the first detection, two particular observations by Potter et al. (2002) and Schleicher et al. (2004) shed new light on the sodium cycle on Mercury and represented a milestone in the study of the its exosphere.
The observation by Potter et al. (2002) showed a long tail of escaping sodium neutral particles. This observation revealed how the radiation pressure of the Sun increases both the probability of sodium escape and the net sodium transport from the day side to the night side. Since the radiation pressure depends on the Doppler shift between neutral atoms and solar spectrum, the radiation pressure is not constant, but varies over the course of the Mercury year (Smyth and Marconi 1995) -and also varies, with some positive feedback effect, with the acceleration of the particles towards the night side.
Shortly thereafter, Schleicher et al. (2004), published the first sodium observation obtained during the 2003 transit of Mercury in front of the solar disc. Transit observations have the peculiarity of immediately providing information about the distribution of sodium at the terminator, thus revealing any dawn-dusk asymmetries and other latitudinal effects. The observations (Fig. 16, right) actually showed an excess of sodium at dawn, while sodium seemed virtually absent on the dusk side. Furthermore, it was noted that sodium increases in correspondence with the polar regions. Mura et al. (2009) reproduced these observations by means of an exospheric model with a rotating surface (Fig. 16, right), starting from the evidence, already noted by Sprague (1992), that the major processes of release occur on the day side, while the nocturnal surface necessarily fills with sodium, thanks to the exospheric migration and the radiation pressure. The sodium particles that precipitate on the night side are released only as the night side reaches the terminator, i.e. at dawn. This nocturnal filling and dawn release mechanism can explain both the asymmetry and, almost quantitatively, the column density of sodium at dawn observed from ground. Deviations between the observations and the model can easily ascribed to other replenishment mechanisms such as volume diffusion (Killen et al. 2004) or precipitation of micrometeorites.
In the observations of Schleicher et al. (2004), another feature that reveals much about the sodium cycle on Mercury is the presence of two bulges close to the north and south poles. Because of the observation geometry, it is not possible to tell the exact region of emission of these populations but it looks quite obvious that they may be ejected from the cusp regions, i.e. where the solar wind enters the magnetosphere and precipitates onto the surface (e.g. Mangano et al. 2013).
To test this hypothesis, the solar wind conditions at the time of observation were modeled to estimate the flux of the solar wind precipitating onto the cusps, and the resulting augmentation of exospheric sodium over these regions (Sarantos et al. 2008. Such simulations also provide the correct N/S asymmetry (modulated by the x component of the Interplanetary Magnetic Field), thus increasing the testability of the proposed mechanism. Finally, measurements of sodium's energy distribution from observations of the Doppler shift in its emission spectrum (of neutral sodium in the radial direction), and of the exospheric scale height, provide constraints on the Na exospheric temperature. Such Doppler measurements of the velocity distribution of the Na particles are reproduced only if the ejection source is assumed to have a temperature around 1000 K-1500 K. Hence, the simplest explanation is that the source is photon-stimulated desorption, modulated by proton precipitation through the cusps that induces a more intense diffusion (Sarantos et al. 2008 and, eventually, release from the surface. In principle, ion sputtering (IS) could work as well, but it would be necessary to review its energy distribution with respect to known laboratory measurements (Wurz and Lammer 2003).
Regarding the annual variability, and the persistence of the dawn-dusk asymmetry at different True Anomaly Angles (TAA) angles, Leblanc and Johnson (2010) have extensively studied the annual trend of Mercury's exosphere and, according to their model, the asymmetry is substantially permanent. Note that dawn, technically, changes position -from 90 • east of the sub-solar point, to 90 • west -for a short time close to perihelion, because in this period Mercury has a slightly retrograde motion. This work suggests that current models fail to reproduce the seasonal variability as established from existing ground-based observations. More detail about these models is provided in the companion paper by Leblanc et al. At the time of the 2006 transit of Mercury, new observations were performed by Potter et al. (2013), who did not observe the dawn/dusk asymmetry. This would suggest a strong variability as a function of the TAA. This effect has been studied extensively (Leblanc and Johnson 2010;Mangano et al. 2013;Milillo et al. 2020). However, the main issues with ground-based observations are that they suffer from the inevitable discrepancies between the calibrations of different telescopes, from the phase angle's high variability, and from the difficulty in disentangling the trend with true anomaly. Hence, before the MESSEN-GER (Mercury Surface, Space Environment, Geochemistry and Ranging) mission, it was thought that the day-night cycle (and the resulting dawn/dusk asymmetry) was a fairly robust paradigm, since the effect that brings the sodium particles onto the night side is basically permanent. The first global observations by MESSENGER (Cassidy et al. 2016, Fig. 2, right) have substantially modified the previous understanding of the sodium day/night cycle, or at least implied the presence of some effects not yet considered. In fact, such observations indicate that at certain true anomalies the sodium atmosphere peaks in the afternoon. A peak in the exospheric sodium column density persists over the cold planetary longitudes, known as "cold poles".
The surface temperature plays an important role in modelling the sodium exosphere. First, thermal desorption produces a thermalized, low altitude population of sodium particles quickly evaporating from the dayside surface. Surface temperature also controls the diffusion rate (Killen et al. 2004) which brings fresh sodium atoms to the uppermost layer of regolith grains, ready for photon-stimulated desorption to release them. In simplified models, the Mercury's surface temperature is assumed to depend only on some power of the solar zenith angle in the dayside (e.g., 1/4), and assumed to be uniform on the night side (Leblanc and Johnson 2003;Vilas et al. 1988;Wurz and Lammer 2003). As such, the sub-solar point and night side temperatures depend only on the distance from the Sun. However, due to the 3:2 spin/orbit resonance, the longitudes at 90 • and 270 • are at the subsolar point always in correspondence with the aphelion, causing these regions to receive, on a two-yearly cycle, less heat than any other equatorial region. Cassidy et al. (2016), observing that the sodium exosphere increases above the longitudes of the cold pole of Mercury, deduced that the surface at these longitudes functions as reservoirs/sinks for Na particles. The average temperature being lower there, over long timescales, causes less Na release on average, which may result in the formation of two sodium reservoirs at 90 • and 270 • . According to these authors, the reservoirs are filled when the anti-sun sodium transport is maximum and when the cold poles are close to the terminator.
This observation requires very complex modelling to be reproduced. In principle, to account for the presence of cold poles, one could just improve the surface temperature models, using for example more sophisticated ones. Unfortunately, this approach does not work well, neither by applying the Mura et al. (2009), nor the Leblanc and Johnson (2010) models. In general, the issue is that when a cold pole region is on the night side, it isn't a "sink" more than any other element at other longitudes (any sodium particle falling onto the night side remains close to where it falls). On the other hand, when the "cold pole region" is on the day side, it cannot emit sodium and act as a sink at the same time. In other words, the models are able to accumulate sodium at the cold poles, but, when Mercury is at TAA = 250 • (bottom of Fig. 17, right panel) such a "modelled" cold pole has already crossed through the whole dayside, and presumably has lost most of the surface sodium. However, on average, a model like that by Mura et al. (2009), if fed with an updated temperature map, is able to predict an average excess of sodium on the surface at 90 • and 270 • (on a 2-year/full solar day average). If we assume that such a cold pole enhancement has worked for ages in the past, causing a permanent excess of sodium at 90 • and 270 • , then we are able to better reproduce the observations (Fig. 17, right panel). More modelling is needed to fully understand these observations, and detailed steps for improving these simulations are offered in the companion paper by Leblanc et al. (this journal). Since it is difficult to separate annual variability from local time by means of ground-based observations, MESSENGER data will remain the only useful dataset for such studies until BepiColombo data becomes available in 2026.

Surface Features & Exospheric Replenishment
An important consideration in the analysis of Mercury's surface-exospheric interaction is the presence of geological features with unique surface compositions, which may supply the exosphere with inhomogeneously distributed ejected/sputtered or outgassed source species. These localized 'primary' gas sources act together with 'secondary' surface-exospheric interactions to determine Mercury's exospheric structure and composition. It is therefore essential to consider the spatial variations of the exosphere and their possible relationships with the surface composition and its geological features (Wurz et al. 2010). Mercury's surface is characterized by numerous geological features that are evidence for volcanic and tectonic activity in the past (e.g., Byrne 2014; Strom et al. 1975), when the presence of volatiles certainly played a major role in crustal evolution. Widespread effusive volcanic plains cover the planet and several vents of pyroclastic origin are present (Thomas et al. 2014). Thousands of scarps, which are the surface expression of the planet's global contraction (Strom et al. 1975), often form long prominent alignments (e.g., Byrne 2014), some of which likely formed under fluid or gas overpressure (Galluzzi et al. 2019). Such a geological asset certainly caused gas migration and seepage during Mercury's early history, at least up to 1.7 Ga, when the pyroclastic activity, the youngest of these processes, likely ended (Thomas et al. 2014). Today, just a few small fault-scarps are thought to be active (Watters et al. 2016), and albeit too small to be sources of high degassing, they are proof of a planet which is still contracting today. Volatile ascent is known to be particularly effective on a cooling planet, and the surficial layer of Mercury's crust is today a permeable net of fractures caused by fault patterns and impact craters (up to 4.5 km deep) that could easily permit outgassing.
Mercury's exospheric replenishment by surface feature activity likely takes place at several scales today. At the regional scale, Mercury is characterized by large geochemical terranes (Weider et al. 2015) mapped by using the data acquired by the X-Ray Spectrometer . The most extensive terrane is the high-Mg region covering a large part of Mercury's surface (Fig. 18a). This same region also appears to be characterized by elevated percentages of Ca , although the available map coverage is limited due to the possibility of measuring the abundance of Ca during solar flares only (Weider et al. 2015). MESSENGER measurements highlighted a correlation between the Mg-rich surface and the Mg exosphere (Merkel et al. 2018).
At the local scale, hollows are peculiar features that could contribute to the exosphere. Hollows are shallow irregular and rimless flat-floored depressions with bright interiors and halos, often found on crater walls, rims, floors and central peaks (Fig. 18b) (Blewett et al. 2018;Thomas et al. 2014). These features, from tens of m to tens of km, are scattered all over the surface. Due to their fresh appearance, the hollows are believed to form via a mech-  (Schaible et al. 2020), chemical sputtering, micrometeorite impact vaporization and pyroclastic volcanism. To understand if the hollows volatiles could contribute to the exosphere, laboratory measurements on the PSD of CaS powder (an analogue for oldhamite), which could be a predominant component of hollows identified within craters, have been performed. These measurements have been used to model the release of Ca by PSD from the Tyagaraja crater, which is characterized by hollows on its floor. The modeling showed that the localized neutral micro-exosphere produced from Ca PSD can be substantial even if only 1% CaS is assumed in the hollows field located in the crater (Bennett 2016). It is possible that hollows-forming volatile material could contribute to the Hermean exosphere, although the existing measurements do not have sufficient resolution to detect such a correlation.
The observed correlation between the Mg in the exosphere and the surface regions having the highest Mg/Si ratios (Merkel et al. 2018) demonstrates how surface composition data may be important in understanding the exospheric processes and variation. At the Moon the same conclusion was reached by the finding that the potassium exosphere is densest over the potassium-rich KREEP soils (Colaprete et al. 2016;Rosborough et al. 2019). MESSENGER was able to derive information on Mercury's surface composition with a spatial resolution of tens or hundreds of km and a low spectral sampling. Other refractory elements, like Ca, may exhibit a similar correlation, but it is necessary to know the surface composition with better accuracy. Improved measurements will help in distinguishing the contribution of the different source mechanisms, such as meteoroid impact vaporization which is assumed to be relevant for these elements.

Future Work
With the vast and complex array of gas-surface interactions expounded upon in this chapter, atmospheric science at the Moon, Mercury, and at other rocky and icy planetary bodies with surface-bounded exospheres is now at a crossroads. With an ever increasing understanding of the energy and material exchange between the atmospheres and surfaces of planetary bodies with surface-bounded exospheres, the field is now positioned for the development of a new generation of exospheric models. The most sophisticated exospheric models of the Moon and Mercury take into consideration solar energetic particle driven regolith grain diffusion as suggested by McGrath et al. (1986), or the surface sticking time of exospheric volatiles onto an assumed smooth surface, and treat cold trapping with a simple stochastic approach (e.g. Grava et al. 2015;Leblanc and Johnson 2010;Mura et al. 2009). Details including local temperature variations, versus location and depth, and the position and time dependence of topographical shadowing are now beginning to be implemented into exospheric models of the Moon (Prem et al. 2018;Schorghofer et al. 2017). Higher-fidelity models at the Moon, Mercury, and elsewhere would also take into account improving knowledge of the effects of porosity on outgassing from regolith Tsavachidis 2020, 2021). New estimates of the binding energies and ejection energy distributions of sputtered atomic species such as Na from surfaces of different mineral compositions are being furnished by molecular dynamics simulations (Morrissey et al. 2022), which will be highly useful in improving the accuracy of future exospheric models. Epiregolith thermal gradients should also be considered when modeling exosphere transport and surface volatile distributions. Such simulations would have the capability to map the time-evolution of surface 'frosts', in the shadows of mountains or craters, or in PSRs, with the changing solar position, taking into considering both (1) the diffusion of gas into the regolith, including the ice 'pump' phenomenon and epiregolith gradient effects discussed in this chapter, and (2) the removal/addition of volatiles to and from the exosphere.
A detailed understanding of the migration of frozen exospheric volatiles across the surface, and within the regolith, would empower future remote and in situ observers to understand and pinpoint to an unprecedented degree the probable locations and abundances of ices and other volatile repositories on the Moon and Mercury. The answers to these questions is highly dependent on the degree to which these 'ices' are sputtered or otherwise radiation processed; an understanding of which requires knowledge of the distribution versus position and depth of these materials. Exospheric simulations which simultaneously track the detailed surface thermal and frozen volatile distribution, including topographical effects , require significant computing resources, but ever increasing computing capabilities have for the first time made such models practical and feasible to carry out. The combined power of a new generation of models and experiments to inform and interpret spacecraft remote and in situ observations portends a new phase in the exploration of planetary atmospheres not only at the Moon and Mercury, but at the icy bodies of the outer solar system and beyond.
The effort to understand gas-surface interactions throughout the inner solar system will soon be assisted by a new generation of measurements from solar system exospheres. The renewed emphasis on lunar exploration will provide the international community data from commercial lunar landers equipped with a variety of spectrometers. Constraining the true identity of exospheric detections, such as the mass 44 signal ostensibly attributed to CO 2 by LACE (Sect. 7), will be a high priority for such instruments, and essential to understand the degree to which these species stick to and diffuse on the lunar surface. In situ determination of activation energies for desorption can be expected from the upcoming Volatiles Investigating Polar Exploration Rover (VIPER) mission. These measurements will eventually be complemented by packages to be deployed by humans during the Artemis missions. BepiColombo's SERENA instrument package (Orsini et al. 2021) will map the density, distribution and composition of Mercury's exospheric gas species through in situ measurements (with Strofio) and remote sensing of charge-exchanged species (with ELENA). Future spacecraft measurement of exospheric Na, Ca, K and S and their spatial relationship with the hollows will contribute to understanding not only of the relationship between the surface and exosphere, but also to the evolution and formation of the hollows. The BepiColombo MERTIS instrument will also provide critical observations to determine the magnitude of epiregolith thermal gradients at Mercury (Hiesinger and Helbert 2010). Additionally, the spacecraft's SIMBIO-SYS remote sensing instrument package (Cremonese et al. 2020) will improve the accuracy of composition data, providing global coverage of Mercury's surface at high (480 m/pixel) spatial resolution. These measurements to be carried out during the coming decade will improve on current understanding from MESSENGER mission data of Mercury's exospheric-surface interactions.