Using Io's Sulfur Isotope Cycle to Understand the History of Tidal Heating

Stable isotope fractionation of sulfur offers a window into Io's tidal heating history, which is difficult to constrain because Io's dynamic atmosphere and high resurfacing rates leave it with a young surface. We constructed a numerical model to describe the fluxes in Io's sulfur cycle using literature constraints on rates and isotopic fractionations of relevant processes. Combining our numerical model with measurements of the 34S/32S ratio in Io's atmosphere, we constrain the rates for the processes that move sulfur between reservoirs and model the evolution of sulfur isotopes over time. Gravitational stratification of SO2 in the upper atmosphere, leading to a decrease in 34S/32S with increasing altitude, is the main cause of sulfur isotopic fractionation associated with loss to space. Efficient recycling of the atmospheric escape residue into the interior is required to explain the 34S/32S enrichment magnitude measured in the modern atmosphere. We hypothesize this recycling occurs by SO2 surface frost burial and SO2 reaction with crustal rocks, which founder into the mantle and/or mix with mantle‐derived magmas as they ascend. Therefore, we predict that magmatic SO2 plumes vented from the mantle to the atmosphere will have lower 34S/32S than the ambient atmosphere, yet are still significantly enriched compared to solar‐system average sulfur. Observations of atmospheric variations in 34S/32S with time and/or location could reveal the average mantle melting rate and hence whether the current tidal heating rate is anomalous compared to Io's long‐term average. Our modeling suggests that tides have heated Io for >1.6 Gyr if Io today is representative of past Io.


Introduction
Io is the most volcanically active body in our solar system.Its volcanism is driven by tidal heating, powered by the orbital resonance of Io, Europa, and Ganymede around Jupiter (Peale et al., 1979).The high rates of volcanism lead to high resurfacing rates (and therefore a young surface), leaving little surface record of Io's past (Johnson et al., 1979).Additionally, Io's day-night temperature variations fluctuate around the freezing point of SO 2 (∼113 K, <10 9 bar), resulting in at least the lower atmosphere collapsing to form surface frosts and then subliming to re-release vapor every ∼41 hr (and additionally during eclipse; e.g., Tsang et al., 2016).This dynamic nature of Io makes the history of tidal heating difficult to constrain (de Kleer et al., 2019).For instance, the Laplace resonance is thought to have started quickly after the formation of these satellites (e.g., Peale & Lee, 2002), but it is not observationally confirmed how long the orbital resonance has been active.Cyclic behavior of the tidal heating has been proposed, with quantitative models indicating that the heat flow varies by an • Rayleigh distillation, where gravitational stratification drives isotopic fractionation, is a good approximation for Io's sulfur isotope cycle • Efficient mixing and recycling between Io's interior and atmosphere occurs by crustal sequestration, volcanic frost remobilization, and burial into the mantle • The difference between mantle and crustal frost 34 S/ 32 S decreases as mantle melting rate increases; atmospheric variability could measure this Supporting Information: Supporting Information may be found in the online version of this article.order of magnitude over oscillation periods of ∼100 Myr (e.g., Hussmann & Spohn, 2004;Ojakangas & Stevenson, 1986).However, this is not observationally constrained, and it is unknown whether the current heating rate is anomalous or typical of Io's long-term average (e.g., Bierson & Steinbrügge, 2021).
Stable isotopes offer a potential window into the past for systems such as Io (e.g., de Kleer et al., 2019), especially isotopes of elements found in constituents of the atmosphere and surface frosts that are measurable via remote and/or in situ techniques.Such measurements must be interpreted through isotope-enabled models that describe the sources, sinks, and distributions of these species (e.g., Donahue et al., 1997;Hunten, 1973;Jakosky, 1991).Recent examples include models constructed to interpret isotope measurements of Titan's modern atmosphere.Such models have been used to constrain the initial 15 N/ 14 N ratio for N 2 and the D/H ratio for CH 4 , as well as the time-scale for methane outgassing from the interior (Mandt et al., 2009(Mandt et al., , 2012)), and to place constraints on the timing and source of nitrogen in the atmosphere (Erkaev et al., 2020).For Pluto, modeling nitrogen isotopes highlighted the importance of condensation and aerosol trapping for the composition of HCN in the atmosphere (Mandt et al., 2017).On Mars, the D/H ratio of the atmosphere today can be explained by atmospheric escape modulated by sequestering water via crustal hydration (Scheller et al., 2021), whilst coupled CO 2 -N 2 -Ar isotope modeling was used to constrain the size and composition of the ancient atmosphere (Thomas et al., 2023).A combination of isotopic fractionation in impact-driven hydrothermal systems with photochemistry and chemical reactions in the atmosphere can explain the wide range of sulfur isotope ratios measured in Gale crater on Mars (Franz et al., 2017).
On Io, sulfur (mostly as SO 2 ) is the dominant volatile species; it covers the surface as frosts (e.g., Nelson et al., 1980), forms the main constituent of the atmosphere (e.g., Lellouch et al., 1990), is emitted from the many volcanoes (e.g., McGrath et al., 2000), and escapes from the atmosphere to populate Io's orbit as plasma (e.g., Broadfoot et al., 1979).Additionally, the natural isotopic abundances of the two most common stable isotopes of sulfur are sufficiently abundant ( 32 S ∼ 95% and 34 S ∼ 5%) that millimeter and/or spacecraft observations/sampling can detect and measure the isotopologues of sulfur-bearing species (e.g., de Kleer et al., 2024;Franz et al., 2017;Moullet et al., 2013).Hence, measurements of sulfur isotopes in the frost, atmosphere, and/or volcanic plumes of Io could be used to further our understanding of the geological history of Io.
de Kleer et al. (2024) measured the isotope ratio of sulfur ( 34 S/ 32 S) in Io's atmosphere using the Atacama Large Millimeter/submillimeter Array.Two observations of Io were made covering a frequency range containing multiple rotation transitions of SO 2 , SO, their isotopologues, as well as other species.A radiative transfer model was used to determine the 34 S/ 32 S ratio by fitting the observed emission lines.Although spatial data were collected, a single ratio and accompanying error was calculated for each of the two observations (values are consistent within 2σ) and both observations together.For both observations together, the 34 S/ 32 S ratio of SO 2 in Io's lower atmosphere is 0.0595 ± 0.0038 (equivalent to +347 ± 86‰ δ 34 S VCDT , where δ 34 S VCDT is defined in Equation 1 of Table 1).The signal-to-noise ratio for SO data was too low to fit a 34 S/ 32 S ratio of SO independently to the 34 S/ 32 S ratio of SO 2 .However, assuming that the 34 S/ 32 S ratio of SO was the same as the 34 S/ 32 S ratio of SO 2 was consistent with the measured SO data.
Assuming an initial δ 34 S VCDT of ∼0‰, de Kleer et al. ( 2024) argued that at steady state (i.e., atmospheric source and sink rates are equal, and the atmospheric source has a δ 34 S VCDT of 0‰) the atmosphere would be +83‰ δ 34 S VCDT if gravitational stratification and atmospheric escape are the primary causes of isotopic fractionation.
One explanation of the discrepancy between measurements and expectations is that the exceptional 34 S enrichment of the atmosphere largely reflects extensive time-integrated escape from a planetary sulfur reservoir (i.e., a physical domain that contains a finite amount of sulfur that is isotopically well-mixed) that actively exchanges with the tenuous atmosphere where isotopic fractionation occurs (de Kleer et al., 2024).Assuming Rayleigh distillation (i.e., an infinitesimal amount is removed from a reservoir and the residue is perfectly mixed after each removal), and that gravitational stratification followed by loss from the upper atmosphere is the cause of isotopic fractionation, requires ∼94%-99% of the initial sulfur in the reservoir to have been removed (de Kleer et al., 2024).This suggests that either much of Io's sulfur is in the core or has been lost in the past (likely during an early period of more rapid volcanic outgassing), or the initial 34 S/ 32 S ratio of the reservoir was much higher than assumed.
However, this first-order interpretation presumes that geological cycling of sulfur on Io can be approximated by Rayleigh distillation of the planetary reservoir through atmospheric escape.This proposition must be examined by a more detailed consideration of sizes, distributions, and isotopic compositions of Io's sulfur reservoirs; rates and fluxes of relevant processes involving transfer of sulfur from one reservoir to another; and the possibility that isotopic fractionations associated with processes other than gravitational stratification and atmospheric escape contribute to the overall sulfur isotope evolution of the atmosphere.
In this paper, we present a numerical model for the sulfur cycle on Io and describe constraints on the sizes and sulfur isotopic ratios of the different reservoirs, as well as the rates and isotopic fractionation factors for the different processes (Section 2 and Sections S2-S4 in Supporting Information S1).We investigate how, within this conceptual framework, the sulfur isotopic compositions of these reservoirs evolve over time, as well as assess the model sensitivity to different parameters (Section 3 and Section S5 in Supporting Information S1).We explore the appropriateness of Rayleigh distillation involving gravitational stratification, escape to space, and efficient mixing as an approximation for Io's sulfur isotope cycle; what processes enable efficient mixing; how large the sulfur reservoir is that interacts with the atmosphere; and whether we can constrain the average mantle melting rate (Section 4).Section 5 summarizes our key conclusions.

Model
The sulfur cycle on Io is complex: sulfur moves between different reservoirs by way of multiple processes, several of which involve the atmosphere (e.g., de Pater et al., 2021de Pater et al., , 2023; Figure 1), and most of these processes are at least potentially capable of fractionating sulfur isotopes.Based on this conceptual model, we build a numerical isotope-enabled box model to track the four stable isotopes of sulfur over time (Figure 2; Hughes, 2024).Although we do not explicitly track the chemical species of sulfur in the model, the isotopic fractionation factors and what reservoirs sulfur is transferred between are informed by the speciation of sulfur.We focus on 32 S and 34 S due to the recent measurements in the lower atmosphere (de Kleer et al., 2024) but include 33 S and 36 S because measurements of these isotopes may be made in the future.
First, we constrain the initial amount of sulfur in the mantle and its isotopic composition (Section 2.1).From this initial condition, our numerical model calculates the amount of sulfur and its isotopic composition for the mantle, crustal frosts, crustal silicates/sulfates, and space over time (see Section S3 in Supporting Information S1 for details).This is based on sulfur fluxes and isotopic fractionation factors ( 3n α, Equation 2 in Table 1) detailed in Section 2.2.Unless otherwise stated (i.e., for gravitational stratification, photo-dissociation, and electro-ionization), 3n α for n = 3 and 6 are calculated assuming canonical mass-dependent fractionation laws (Table 1) using Equation 3 in Table 1.The space reservoir is the cumulative material lost to space from photo-ionization and plasma interactions.It is tracked as a model reservoir for mass conservation but is not a meaningful physical reservoir because this "reservoir" does not return to the atmosphere and thus has no subsequent exchange with the other model reservoirs.The atmosphere reservoir in our model is treated differently to other reservoirs because the atmosphere represents a tiny proportion of the sulfur on Io (Table 2, Section S1.4 in Supporting Information S1), which would make the size of the time-step required for modeling unfeasible otherwise.At each time-step, sulfur is transferred from the mantle (and crustal frosts) due to outgassing (+sublimation) to the atmosphere, and then removed by photo-dissociation + space loss (+deposition) to crustal frosts and space.As the atmosphere size is constant over time (e.g., Tsang et al., 2012), these fluxes balance such that there is no sulfur in the atmosphere reservoir between time-steps.We feel this is a justified simplification given that the atmosphere size (1.8-3.2 × 10 9 mol S, Table 2) is many orders of magnitude smaller than the flux of sulfur into or out of the atmosphere over a time-step (for the 0.1 Myr time-step used, the flux is ∼10 18 mol S).However, the isotopic composition of the atmosphere at each time-step can be calculated given the assumption that it is in equilibrium with the crustal frosts (Section S1.4 in Supporting Information S1).

Amount of Sulfur and Its Isotopic Composition in the Initial Mantle
We assume that the initial bulk composition of Io is represented by either L/LL ordinary chondrite meteorites (Dreibus et al., 1995;Gao & Thiemens, 1993;Kuskov & Kronrod, 2001) or solar system proportions (Lodders, 2021;McKinnon, 2007) to give lower and upper bounds, respectively, on sulfur concentration (Table 2; Section S1.1 in Supporting Information S1).Currently, there are no constraints on the potential extent of sulfur evaporation from an early magma ocean for Io.Thus, we do not consider this process, although it would have little effect on the bulk sulfur isotopic composition (Wang et al., 2021).We constrain the sulfur content of the core and the initial mantle by combining constraints on Io's mantle/core density and core radius from Io's mean density and moment of inertia (Sohl et al., 2002) with experimental metal-silicate sulfur partition coefficients at Io's coremantle boundary (Keszthelyi & McEwen, 1997a;Keszthelyi & Suer, 2023;Keszthelyi et al., 2004Keszthelyi et al., , 2007;;Kuskov & Kronrod, 2000;Moore, 2001;Suer et al., 2017;Section S1.2 in Supporting Information S1).While the density and moment of inertia constraints alone permit 0%-100% of Io's initial sulfur to reside in the core, incorporating the partition coefficient constraints and assuming equilibrium results in 80%-97% of Io's bulk sulfur partitioning into the core (5.1-14.7 × 10 22 mol S, where the minimum assumes L/LL ordinary chondrite for Io's bulk composition, whilst the maximum assumes solar system proportions; Table 2).Given the small sulfur isotopic fractionation factor between metal and silicate ( 34 α = 0.9998; Labidi et al., 2016), there is minimal isotopic fractionation when the core differentiates from the mantle (Section S1.2 in Supporting Information S1).
Once sulfur in Io is partitioned between the core and the mantle, we assume that there is negligible interaction between these reservoirs.Hence, the core is not included in our numerical model.The amount of sulfur and its isotopic composition in the mantle reservoir are initial inputs into our numerical model.We assume the crust is 40 km thick as estimates range from 20 to 50 km based on the maximum height and total volume of tectonic mountains (Carr et al., 1998;Jaeger et al., 2003;Keszthelyi & McEwen, 1997b;Kirchoff & McKinnon, 2009;Ross & Schubert, 1985;Schenk et al., 2001).

Sulfur Fluxes and Isotopic Fractionation Factors
When the mantle partially melts (Section S2.1 in Supporting Information S1), sulfur initially held in the mantle dissolves into melt as sulfide, up to the limit of sulfide-saturation at ∼1,000 ppm (Battaglia et al., 2014).If the mantle is sulfide-undersaturated (i.e., due to previous sulfur extraction), we assume the melt sulfur concentration equals that of the mantle as the degree of mantle melting on Io is uncertain.There is negligible sulfur isotopic fractionation between the melt and the residual for the case where all sulfur is present as sulfide in both materials, based on the isotope composition of coexisting silicate and sulfide in natural samples ( 34 α = 1.0000 ± 0.0003; Labidi & Cartigny, 2016;Labidi et al., 2014;Mandeville et al., 2009).Some studies suggest that the mantle may be split into an upper and lower mantle (Spencer, Katz, Hewitt, May, & Keszthelyi, 2020), which we include as an option with an inter-mantle flux of sulfur between them.Our model assumes that mantle melts are either intruded into the crust as plutons or erupted at the surface and degassed; based on Spencer, Katz, and Hewitt (2020), we assume 80% is intruded (Section S2.2 in Supporting Information S1).Current resurfacing rates on Io are ∼1 cm/yr (e.g., Johnson et al., 1979), which are the erupted melts in our model from which we calculate the mantle melting and pluton emplacement rate (Sections S2.1 and S2.2 in Supporting Information S1).
Degassing is assumed to be complete due to the low pressure at Io's surface.Hence, all magmatic sulfur is transferred to the vapor for erupted melts with no net isotopic fractionation, whilst plutons do not degas (Section S2.2 in Supporting Information S1).
the surface, SO 2 sublimates/deposits with isotopic fractionation controlled by the vapor-pressure isotope effect.This is likely to be small relative to other fractionations considered based on measurements of CO 2 and SF 6 ice (e.g., 34 α ∼ 0.9969-1.0004;Eiler et al., 2000Eiler et al., , 2013: Section S2.9 in Supporting Information S1).We assume that any sulfur compounds not lost to space (next paragraph) or broken down by photo-dissociation are deposited as crustal frosts (or if there is insufficient sulfur, sublimated from the crustal frosts) such that the atmosphere size is constant over time (e.g., Tsang et al., 2012).
Various plasma and photon interactions result in sulfur in the atmosphere being lost to space, principally through reactions that generate ionic products, which are then accelerated away from Io by Jupiter's magnetosphere (Bagenal & Dols, 2020;Delamere et al., 2004; Section S2.12 in Supporting Information S1).The space loss rate is bounded by the measured supply rate to the torus and the estimated production rate (Bagenal & Dols, 2020;Delamere et al., 2004).Photo-ionization of SO 2 is not isotopically fractionating and is a relatively small proportion of the sulfur lost to space (e.g., Bagenal & Dols, 2020;Croteau et al., 2011;Saur et al., 1999;Section S2.11 in Supporting Information S1).The upper atmosphere will be gravitationally stratified (e.g., gradients in isotopic composition can form in <1.5 hr), which produces large isotopic fractionations between the homopause and the exobase ( 33 α = 0.9574, 34 α = 0.9168, and 36 α = 0.8404: de Kleer et al., 2024; Section S2.8 in Supporting Information S1: calculations based on Giunta et al., 2017;Seltzer et al., 2017).Gravitational stratification scales with the absolute difference between isotope masses rather than the relative difference and therefore can contribute to mass-anomalous fractionations as θ 33 g = 0.501 and θ 36 g = 2.0 (calculated using Equation 3 of Table 1) rather than the canonical values of 0.515 and 1.9, respectively (see also Dauphas & Schauble, 2016).Thermal stratification (also called the Soret effect) between the homopause and the exobase is predicted to be relatively small ( 34 α = 0.9998: based on experiments by Wullkopf (1956); Section S2.8 in Supporting Information S1).We assume that the plasma-neutral processes occurring (and their rates; Table 2) are those described in Figure 8 of Bagenal and Dols (2020).Electro-ionization has a negligible isotope effect (e.g., Basner et al., 1995;Section S2.13 in Supporting Information S1).We estimate the sulfur isotope fractionation factor for electro-dissociation ( 33 α = 0.996, 34 α = 0.994, and 36 α = 0.989) using the experiments of Ustinov and Grinenko (1971) (Section S2.13 in Supporting Information S1).Asymmetric charge exchange between S ++ or O + and SO 2 is assumed to be nonfractionating, but this simplifying assumption is unlikely to lead to systematic errors in our model as it is an insignificant loss process (Section S2.14 in Supporting Information S1).Resonant charge exchange between SO 2 + and SO 2 is unlikely to be isotopically fractionating as the collisional energies are much larger than the threshold energies required for charge exchange (e.g., Bodo et al., 2008;Hodges & Breig, 1993;Zhang et al., 2011;Section S2.14 in Supporting Information S1).Elastic collision between SO 2 and SO 2 + is not isotopically fractionating based on the high energy of the incoming SO 2 + from the torus compared to the escape energy (e.g., Chassefière & Leblanc, 2004;Johnson et al., 2000;Section S2.14 in Supporting Information S1).Pick-up ion escape is not isotopically fractionating because the energy of the pick-up ion is much larger than the energy required to escape Io's atmosphere (Chassefière & Leblanc, 2004;Dols et al., 2008;Section S2.15 in Supporting Information S1).Overall, the isotopic composition of material lost to space is dominated by the effects of gravitational separation, as the energies involved in the plasma interactions are high enough that negligible isotopic fractionation occurs.
All sulfur on the surface and in the crust is slowly buried as further resurfacing occurs (Section S2.16 in Supporting Information S1).The crustal thickness stays constant because it is controlled by temperature and melting.Hence, surface-deposited and crustal sulfur compounds are either returned to the mantle or recycled to the surface.We infer that sulfur that remains in the form of buried frost and fluids will be retained in the crust or transferred to the atmosphere, whereas sulfur as sulfates and sulfide accessory minerals in igneous rocks can subside to the crust-mantle boundary and return to the mantle.The amount of sulfur returned to the mantle depends on the amount of sulfur in the crustal silicates/sulfates and the amount of silicate added to the crust.A summary of the rates and isotopic fractionation factors is given in Table 3.

Amount of Sulfur and Its Isotopic Composition on Io Today
There are a few constraints on the amount of sulfur and its isotopic composition in various reservoirs on Io today as well as the sulfur flux from some processes that we can compare to outputs from our numerical model.SO 2 frost on the surface may be >1 km thick in some locations (Keszthelyi et al., 2004;Schenk & Bulmer, 1998  The 33 S/ 34 S Io's surface frosts is 0.13 ± 0.07 from infrared spectroscopy ( 32 S was saturated), which is within error of the solar system average (Howell et al., 1989; Table 2).de Kleer et al. ( 2024) measured the δ 34 S VCDT of SO 2 in the lower atmosphere as +347 ± 86‰ (0.0595 ± 0.0038 34 S/ 32 S; Table 2).The flux of SO 2 entering the atmosphere via volcanism (i.e., equivalent to outgassing) can be estimated from the observed atmospheric abundance and an estimated residence time (e.g., Jessup et al., 2004;Lellouch et al., 2003; Table 3; Section S2.4 in Supporting Information S1).We compare our model outputs to these constraints in Section 3.

Results
We use our numerical model of Io's sulfur cycle to constrain the combinations of conditions that can result in the δ 34 S VCDT value measured in Io's atmosphere today (de Kleer et al., 2024) and to explore model sensitivities.

Base Scenario: If Io Today Is Representative of Past Io
The "today" scenario assumes that exchange rates relevant for Io today have been constant throughout Io's history (Table 4).The amount of sulfur in the mantle is set such that ∼97% of the mantle sulfur has been lost after 4.57 billion years (Gyr) based on the Rayleigh distillation calculations by de Kleer et al. ( 2024).Here we describe the different reservoirs and how they evolve with time using our model under this scenario."Space-loss" is the sulfur lost to space in a given time-step, whilst "outgassing" is the SO 2 that is transferred to the atmosphere during a given time-step.We emphasize these model variables because they are potentially measurable by future remote and/or in situ techniques.
The behavior of the different reservoirs in the "today" scenario is shown as a function of time over Io's 4.57 Gyr history in Figure 3.As the model evolves from time zero, the amount of sulfur in the mantle reservoir decreases as sulfur is transferred to the crust and space, whilst the size of the space reservoir increases over time (Figure 3a).The amount of sulfur in the crust as frosts and sulfates/silicates reaches steady-state quickly (<30 Myr), with more sulfur in the form of sulfates/silicates than frosts (Figure 3a).The steady-state reservoir sizes of the crustal frosts and silicates are of a similar order of magnitude to estimates derived from the literature (∼10 20 mol S: Leone et al., 2011; Table 2).After ∼3.2 Gyr, the sizes of both crustal reservoirs begin to decrease because the sulfur flux from mantle melting decreases as the sulfur concentration of the mantle decreases below sulfide-saturation (the rate of sulfur lost-to-space does not change because it is constant in the model).Hence, the sulfur concentration of the mantle melt is no longer at sulfide-saturation but is equal to that of the mantle.As the mantle melt is no longer sulfide-saturated, the rate of mantle outgassing of sulfur decreases over time, even for a constant rate of mantle melting.The rate of outgassing in the model falls between the minimum and maximum estimates based on the estimates of residence time of SO 2 in the atmosphere (Table 3).
The δ 34 S VCDT for all reservoirs becomes more positive over time, that is, 34 S becomes progressively enriched relative to 32 S (Figure 3b).Sulfur isotope ratios of mantle, crust, and outgassed material are closely coupled at any given time, with crustal sulfates/silicates slightly isotopically heavier compared to the mantle, and crustal frosts further enriched in 34 S. The space-loss material is isotopically lighter (lower 34 S/ 32 S) compared to the mantle and Note.See Section S2 in Supporting Information S1 for details on the constraints on the rate and isotopic fractionation factor of each process.§ Rates listed here are model inputs and constraints from the literature that can be compared to model outputs: rates that are only calculated as part of a model run without an independent constraint are not listed here.‡ All rates are given in kg/s SO 2 -eq (i.e., if all sulfur was present as SO 2 ), where n/a indicates not applicable.# indicates 33 α and 36 α are calculated independently of 34 α: otherwise, 33 α and 36 α are calculated from 34 α assuming canonical mass-dependent law coefficients (Table 1) and Equation 3in Table 1."S" in reactions denotes unspecified speciation of sulfur.Process also affected by † gravitational and thermal stratification and/or *pick-up.crust because this material has undergone gravitational stratification in the atmosphere, which enriches the light isotope in the upper atmosphere that is then lost to space.This means that the space reservoir is also isotopically light.The mantle, crust, and outgassed material become progressively isotopically heavier because the light isotope is irreversibly lost to space and the isotopically heavy residual material is mixed back into the crust and mantle.Over time, the space reservoir and space-loss material become isotopically heavier with respect to their earlier values because the reservoirs providing the material being lost to space (i.e., mantle and crust) become heavier over time.Due to mass balance, the space reservoir is always isotopically lighter than the initial δ 34 S VCDT and tends toward the initial δ 34 S VCDT value as more sulfur is transferred to it.

Journal of Geophysical
Similar behavior to δ 34 S VCDT is observed for δ 33 S VCDT and δ 36 S VCDT in all reservoirs and fluxes (Figures 3c  and 3d).The magnitude of the isotopic enrichment compared to δ 34 S VCDT is approximately half for δ 33 S VCDT and double for δ 36 S VCDT , reflecting the mass differences between the different isotopes.The 33 S/ 34 S ratio of the crustal frost after 4.57 Gyr is 0.16, which is in the range measured spectroscopically (Howell et al., 1989; Table 2).We choose not to calculate and compare Δ 33 S and Δ 36 S values (e.g., Farquhar & Wing, 2003;Ono, 2017) because: (a) ∆-values mostly reflect the contrast between the mass law of the gravitational isotope effect and the somewhat arbitrary mass law of the canonical reference frame (Table 1); and (b) at the large δ-values predicted for all sulfur isotopes by our numerical model (Figures 3b-3d) ∆-values behave non-intuitively due to their definitions (discussed further in Section S4.1 in Supporting Information S1; see also Kaiser et al., 2004;Miller, 2002).

Alternate Scenarios: If Past Io Differs From Io Today
We next explore dropping the assumption that exchange rates were the same throughout Io's history as they are today and vary the amount of sulfur in the initial mantle.For each of these different scenarios, we alter an individual parameter of interest to track its effect on Io's isotopic signatures and keep all other parameters the same as the "today" scenario; the scenarios we consider are described in Table 4. Equivalent figures to Figure 3 for the "today" scenario are shown in Figures S10-S21 in Supporting Information S1 for the other scenarios.Below, we focus on the isotope values of the mantle, crustal frosts, and space-loss material after 4.57 Gyr (Figure 4).In our a Uses values representative of today and a mantle size such that ∼97% of the sulfur has been lost after 4.57 Gyr.b Mantle is segregated into an upper and lower mantle with a small flux of lower mantle into the upper mantle.c Scenarios using variable mantle melting rate use the mantle size of today (e.g., Hmm) or that used for BM (e.g., BM-Hmm).d Resurfacing rate changes every 100 Myr between 1 and 9 (where the last 100 Myr is 5 cm/yr), which averages to 5 cm/yr overall (equivalent to Mmm).
model, the lower atmosphere will have the same isotopic composition as the crustal frosts, but the upper atmosphere will vary due to gravitational separation (see Figure S9 in Supporting Information S1 for more details).
In "High space loss" (Hsl: blue square), we set the mantle size to 8.5 × 10 21 mol S (i.e., larger mantle size compared to "today") and increase the space loss rate such that ∼97% of the sulfur has been lost after 4.57 Gyr (i.e., same overall loss fraction as "today").The compositions of the mantle, frost, and space-loss material after 4.57 Gyr are very similar for Hsl and "today."This highlights that the value of the isotopic composition of the atmosphere today is primarily sensitive to the fraction of sulfur lost (i.e., size of reservoirs and rate of loss trade off against each other).There is little effect on the results when varying the intrusive:extrusive ratio of igneous material (i.e., "Low plutons" Lpl: gray square, where the ratio was 20:80 rather than 80:20 in "today") or turning off photo-dissociation (0pd: gray circle).
In "Big Mantle" (BM: red diamond), we increase the mantle size relative to "today" (same value as Hsl), which decreases the loss fraction after 4.57 Gyr.For BM, δ 34 S VCDT of the mantle, frost, and space-loss material and the δ 33 S VCDT of the frost are all less positive compared to "today."In 2M ("Segregated Mantle": blue circle), Io's mantle is split into a lower and upper mantle that do not mix except for a small flux from the lower mantle into the upper mantle.The amount of sulfur in the initial upper mantle plus that added to the upper mantle over 4.57 Gyr from the lower mantle is approximately equal to the amount of sulfur in the initial mantle for "today" to give similar loss fractions.In practice, this means that the total amount of mantle sulfur (i.e., across the lower and upper mantle) is more than double the mantle sulfur in the "today" scenario.This scenario results in a significantly less positive δ 34 S VCDT of the mantle, frost, and space-loss material (by ∼60‰) and δ 33 S VCDT of the frost (by ∼30‰) compared to the "today" scenario.Thus, it is challenging to achieve the high positive δ 34 S VCDT observed by de Kleer et al. (2024) if the upper mantle is continuously replenished from a reservoir with 0‰ δ 34 S VCDT .
Figure 3. Reservoir sizes and isotopic compositions in the "today" scenario, which applies exchange rates thought to be relevant for Io today as having been constant through time and assumes an initial amount of sulfur in the mantle reservoir such that ∼97% of the sulfur has been lost after 4.57 Gyr.Evolution over time of different reservoirs and fluxes for (a) total moles of sulfur and (b-d) δ 3n S VCDT using the left-hand axis and 3n S/ 32 S using the right-hand axis, where n = (b) 4, (c) 3, and (d) 6.
Results from our numerical model are in colored curves (solid: mantle = red, crustal frosts = yellow, crustal sulfates/silicates = blue, and space = gray; dash: outgassing = red and space loss = gray), while the black curves show the Rayleigh distillation results for the scenario described in the text (reservoir = solid, total loss = dash, and instantaneous loss = dot).The black horizontal line in the gray shaded region in (b) is the atmosphere value with error observed today from de Kleer et al. ( 2024) [dK+].
We turn off sulfur sequestration via gas-rock reactions in "no sequestration" (0sq: gray up-triangle), whilst in 0fr (gray down-triangle) we turn off frost remobilization.In the former case, sulfur is not returned to the mantle; while in the latter, frost is not released from the crust by heating from magmatism.In both cases, this causes all the sulfur to be transferred from the mantle to the crust very quickly, resulting in the amount of sulfur in the crust as frosts being two orders of magnitude larger than the estimates derived from Leone et al. (2011) (Table 2; Figures S14a and S15a in Supporting Information S1) for Io today.At the maximum crustal frost size during the model-run, this scenario would have >10% of the crust by mass being sulfur.For 0sq, the δ 34 S VCDT of the frost and space-loss material is slightly more positive (by ∼9 and ∼2‰, respectively) compared to "today" and frost δ 33 S VCDT is more positive (by ∼2‰) (note there is essentially no mantle sulfur by 4.57 Gyr).However, the outgassing rate for much Io's history is one-to-two orders of magnitude greater than the maximum in Table 3.For 0fr, the δ 34 S VCDT of the frost and space-loss material (again, there is essentially no mantle sulfur at 4.57 Gyr) is much more negative (by ∼103 and ∼72‰, respectively) and frost δ 33 S VCDT is more negative (by ∼42‰) compared to "today."Material is supplied to the atmosphere via sublimation of the crustal frosts after the mantle runs out of sulfur for this scenario.
For the mantle size used in "today" and BM, we vary the resurfacing rate: 1 (white and red diamonds), 5 (white and red circles: Mmm), and 10 (white and red up-triangles: Hmm) cm/yr.This has a small effect on the δ 34 S VCDT for frost (by <2‰) and space-loss material (<4‰), and frost δ 33 S VCDT (<2‰).The variation of the mantle δ 34 S VCDT is large for "today" (57‰ between 1 and 10 cm/yr) but small for BM (∼3‰).For the same average mantle melting rate (5 cm/yr), there is a slight difference in the isotopic composition of the different reservoirs between constant (Mmm; white circle) and oscillatory (Omm; gray square) rates (∼10‰ for δ 34 S VCDT ), but this effect is smaller than when changing the mantle melting rate.
From this sensitivity analysis, we find that the intrusive:extrusive ratio of igneous rocks, amount of sulfur sequestration via gas-rock reactions, and mantle melting rate (whether static or oscillatory) have little impact on the δ 34 S VCDT of the atmosphere after 4.57 Gyr.Therefore, these variables are not tightly constrained by our modeling of the observed atmospheric δ 34 S VCDT today.

How Big Is the Accessible Sulfur Reservoir on Io?
The size of Io's accessible sulfur reservoir (i.e., sulfur that can be moved between Io's interior and atmosphere) has a significant impact on how a given isotopic signature is interpreted in terms of the total mass of sulfur lost.Given the δ 34 S VCDT measurement from de Kleer et al. (2024), it is unlikely that an unmixed reservoir (i.e., a lower mantle) is a dominant input into the system (2M), because this would continuously add less fractionated material into the system and reduce how fractionated the isotopes are in the observed atmospheric reservoir (Figure 4a).Assuming the reservoir is the whole mantle, and that the rate of space loss has been the same for the lifetime of Io as it is today, requires the initial reservoir to contain 5.3-5.6 × 10 21 mol S, 12%-14% between our minimum and maximum estimates of the initial amount of sulfur in the mantle (Section 2.1).Therefore, the core contains 91%-92% of Io's initial sulfur if the composition is akin to L/LL ordinary chondrites or 96%-97% if in solar proportions.This would result in a core composition of 13.9 or 27.7 wt% S, radius of 784 or 888 km, and density of ∼6,688 or 5,738 kg/m 3 for L/LL ordinary chondrite or solar proportions, respectively.
However, models incorporating a higher space loss rate and a larger sulfur reservoir in the mantle (Hsl) also produce results consistent with the isotopic composition of Io's atmosphere today (de Kleer et al., 2024; Figure 4a).These scenarios would require there to be less sulfur in the core.Hence, improving our constraints on the size of Io's core-and hence the amount of sulfur available in the mantle-or the cumulative amount of material lost to space is required to allow us to resolve the balancing effects of space-loss rate and mantle reservoir size in our model.Additionally, future measurements by telescopes or spacecraft could improve our knowledge of the temperature structure of Io's atmosphere, which would improve our estimate of the isotopic fractionation due to gravitational stratification.

Rayleigh Distillation Approximates Sulfur Loss From Io: What Is the Mechanism for Efficient Mixing?
We compare our model results for "today" to the Rayleigh distillation (Rd) calculation presented by de Kleer et al. (2024; black curves in Figure 3).For Rd, given an initial isotope value for the whole system and a fractionation factor for the loss process, all that can be calculated is the fraction of material lost required to achieve a certain isotope ratio of the residue.To attribute a certain loss fraction over a given time requires an initial reservoir size and loss rate.To compare Rd to the "today" scenario, we assumed that the loss rate was the space-loss rate of today and the initial reservoir was the same initial mantle size.The space-loss material has very similar isotopic compositions for "today" and Rd (Figures 3b-3d).However, there are differences in the evolution of the frost and mantle δ 34 S VCDT , δ 33 S VCDT , and δ 36 S VCDT (Figures 3b-3d).The frost δ 34 S VCDT for "today" is lower than Rd for the same loss fraction (e.g., after 4.57 Gyr their values are +324‰ and +346‰, respectively; Figure 4).This is partly because the isotopic fractionation factor between material that is lost to space versus buried and recycled is slightly closer to 1 than pure gravitational stratification because of the other processes occurring.Additionally, mixing is not instantaneous as assumed for Rd, which results in less isotopic fractionation occurring over time.Overall, Rayleigh distillation is a good approximation and shows that gravitational stratification followed by loss of material to space from the upper atmosphere is the main fractionating mechanism, and that recycling is highly efficient.
Our model included two recycling processes (Figures 1 and 2).The shallow cycle involves sulfur being added to the crustal frosts via gas equilibria, photo-dissociation, and deposition, and released via frost remobilization driven by the heat from magmatism.The deep cycle involves sulfur being added to crustal sulfates via sequestration, where gaseous SO 2 reacts with silicate rock to form predominantly sulfates.The sulfates are then buried into the mantle and later released via mantle melting and outgassing.Using realistic rates for these various processes (Table 1) gives results close to the Rd calculation, showing that highly efficient recycling is occurring.Given that if there is no sequestration (0sq) or frost remobilization (0fr) the amount of sulfur in the crustal frosts is unrealistically high (i.e., >10% of the crust is frost), both deep and shallow recycling processes are likely occurring.Additionally, if there is no frost remobilization, the atmosphere is only supported by sublimation because the mantle runs out of sulfur early in Io's history, which does not agree with observations of SO 2 being emitted from volcanoes (e.g., McGrath et al., 2000).Hence, the shallow cycle is important for maintaining volcanic outgassing.
Deep recycling of sulfur provides a mechanism to oxidize Io's mantle over time, as we assumed sulfur is extracted from the mantle as sulfide (e.g., Battaglia et al., 2014) and returned to the mantle as sulfide in plutons but also as sulfate due to sequestration in the crust (e.g., Renggli et al., 2019).The overall redox balance of the mantle would depend on the cycles of other redoxsensitive elements such as iron, but if the mantle was oxidizing over time, this would cause the mantle melts to contain significant quantities of dissolved sulfate.At sulfide-saturation, for the same melt composition and temperature, the total sulfur content of a melt containing both dissolved sulfide and sulfate is higher than that of S 2-CSS (sulfide content at sulfide saturation: e.g., Hughes et al., 2023;Wieser & Gleeson, 2023); hence, mantle melts would carry more sulfur.This might help to explain the high SO 2 emissions from volcanoes on Io, which greatly exceed the sulfur budget expected from sulfide-saturated melts (e.g., up to 30% by mass, Cataldo et al., 2002).Additionally, mantle oxidation over time may explain the high oxygen fugacities inferred for Io's volcanic gases (between the Ni-NiO and hematite-magnetite buffers, although volcanic degassing is likely have increased their oxygen fugacity compared to the source mantle; Hughes et al., 2023), in addition to early hydrogen loss (Zolotov & Fegley, 1999).

Constraints on the History of Tidal Heating on Io
Using our bounds on the amount of sulfur in the mantle (Table 2) combined with the bounds on space-loss rates from observations today (Table 3) gives sulfur loss fractions independent of the isotope measurement that range from 0.10 to 2.80 (average being 0.41).These are compatible with the loss-fractions from de Kleer et al. (2024) based on sulfur isotopes (0.94-0.99), but the range is large, highlighting the utility of isotopes for this problem.Assuming Rayleigh fractionation, the smallest fraction of available sulfur lost on Io that is consistent with the smallest δ 34 S VCDT within the measurement uncertainties (+261‰; Table 2) is 0.94 ( f L , de Kleer et al., 2024).Combining this loss-fraction with the smallest amount of sulfur in the mantle within the constraints described above (3 × 10 21 mol S, M i ; Table 2), and the maximum current space loss rates consistent with models and measurements (3,671 kg/s SO 2 -eq, r sl ; Table 3) gives a minimum time interval (t) over which recycling has been occurring of 1.6 Gyr (t = f L M i /r sl ).Higher space-loss rates would shorten this time; if the space-loss rate was 10 times the maximum current space loss rate, recycling must still have been ongoing for at least 160 Myr.Better characterization of the present-day space-loss rate could be achieved using data from current and upcoming missions.
For a given loss-fraction, the difference between the δ 34 S VCDT of the frost and mantle is dependent on the average resurfacing rate and does not seem to be affected by other processes (Figure 5; Section 3.2).For higher resurfacing rates, the δ 34 S VCDT of the mantle and frosts are more similar due to more efficient recycling, whilst the opposite is true for lower resurfacing rates (Figure 4d).Hence, temporal and spatial variability (e.g., bulk atmosphere vs. plume measurements, or dayside vs. nightside measurements) of atmospheric δ 34 S VCDT (as sulfur will be sourced from potentially pure mantle outgassing to pure frost remobilization) could be used to infer the average resurfacing (and therefore mantle melting) rate over Io's history.The most positive δ 34 S VCDT value could be interpreted as representing the frost end-member and the least positive value the mantle end-member.From our modeling, a resurfacing rate of 1 cm/yr and an initial mantle sulfur amount of 5.4 × 10 21 mol S (the values for "today") produces a ∼65‰ δ 34 S VCDT difference between frost and mantle.Hence, if the actual ‰ difference is smaller than this, it suggests that today's rate is anomalously low, whilst a larger value suggests that the rate today is anomalously high.Such differences would be resolvable using an instrument capable of measuring δ 34 S VCDT with approximately ±10‰ precision.This assumes that the initial amount of sulfur in the mantle is 5.4 × 10 21 mol S; if this is not the case, the relationship between mantle melting rate and mantle/frost difference is different (e.g., for BM in Figure 5).This calculation becomes less sensitive the smaller the loss fraction.However, this highlights how isotopes on Io could be used to constrain the long-term average mantle melting rate.

Conclusions
We explore the parameter space able to reproduce Io's current atmospheric sulfur isotope ratio using a numerical model for the sulfur isotope cycle on Io.Based on Io's mean density and moment of inertia alone, the range of sulfur that could be in the core ranges from 0% to 100% of Io's initial sulfur.However, assuming metal-silicate equilibrium during planetary differentiation and using experimental constraints on sulfur partitioning between silicate and metal narrows the range to 80%-97% of Io's original sulfur that is locked away in the core, restricting the modeled parameter space for initial sulfur.The main factors controlling the isotopic composition of the atmosphere are mantle size, loss rate, and the isotopic fractionation associated with gravitational stratification.We confirm that Rayleigh distillation is a good approximation for the net effect of sulfur loss on the sulfur isotope ratio of atmospheric SO 2 on Io.However, our model also predicts the isotope ratios of the different reservoirs on Io over time, variables that are not considered in the Rayleigh approximation and that could be targets for future measurements.The efficient recycling on Io that is required occurs for plausible rates of sulfur sequestration into the crust and burial into the mantle combined with frost remobilization driven by the heat from magmatism.This highlights the important role that sulfate formation via SO 2 -silicate reactions plays in recycling sulfur on Io, which may have caused oxidation of Io's mantle over time.Under the most conservative assumptions (i.e., requiring the least mass-loss and hence the smallest initial mantle), efficient recycling on Io, which is driven by tidal heating, must have occurred for at least 1.6 Gyr if the current space-loss rates are indicative of past space-loss rates.We find that Io's atmospheric sulfur is too isotopically fractionated at the current day to be accounted for at Io's current mass loss if the core contains <91% of Io's initial sulfur.The difference in δ 34 S VCDT between the mantle and crustal frosts/atmosphere is predicted to depend on the average mantle melting rate, with decreasing average mantle melting rate increasing the difference between the reservoirs.Hence, future measurements of the temporal and spatial variability of the δ 34 S VCDT of the atmosphere and/or surface materials may reveal the history of Io's tidal heating and whether the rate of mantle melting today is anomalous or typical.

Figure 1 .
Figure 1.Conceptual model for the sulfur cycle on Io.

Figure 2 .
Figure 2. Outline of our numerical model of Io's sulfur isotope cycle based on the conceptual model in Figure 1, indicating how sulfur fluxes between reservoirs are calculated.
;Turtle et al., 2001) and we derive an estimate for the mass of sulfur present in the crustal silicates and frosts/fluids

Figure 5 .
Figure 5. Difference between the δ 34 S VCDT for the frost and mantle against resurfacing rate.Symbols and colors as described in Figure 4.

Table 1
Values, Formulas, and Definitions of Isotope Values and Fractionation FactorsNote.VCDT = Vienna Canyon Diablo Troilite, the standard for sulfur isotopes.n = 3, 4, or 6 depending on the sulfur isotope of interest.Note that the sulfur isotope ratios for VCDT are the inverse (i.e., 32 S/ 3n S rather than 3n S/ 32 S) of the ratios used throughout the text to preserve the values as stated in the reference.

Table 2
Sizes and Isotopic Compositions of Sulfur Reservoirs on Io

Table 4
Details of Scenarios Explored