Magmatic Controls on Volcanic Sulfur Emissions at the Iceland Hotspot

Outgassing of sulfur (as SO2) is one of the principal hazards posed by volcanic eruptions. However, S emission potentials of most volcanoes globally are poorly constrained due to a short observational record and an incomplete understanding of the magmatic processes that influence pre‐eruptive S concentrations. Here, we use a compilation of published and new data from melt inclusions (MIs)—which can preserve magmatic S concentrations prior to eruptive degassing—from the Iceland hotspot to evaluate the effects of mantle melting and crustal magmatic processes on the S budgets of Icelandic melts. We use MI data to estimate S emission potentials (ΔSmax, in ppm S) for 73 eruptions from 22 of Iceland's presently active ∼33 volcanic systems. We show that the S systematics of Icelandic melts are strongly regulated by the sulfide solubility limit. Sulfide‐saturated conditions during lower‐degree mantle melting, prevalent at off‐rift zones, likely explains observed decoupling between S and Cl. During magmatic differentiation, a local maximum in modeled sulfide solubility occurs in evolved basalts (4–6 wt.% MgO), coinciding with highest MI S concentrations. Highest ΔSmax values (2,100–2,600 ppm) are found in the Hekla 1913 CE, Eldgjá 939 CE, and Surtsey 1963–1967 CE eruptions in the South Iceland Volcanic Zone. Our results extend the record of volcanic sulfur emissions back in time and can be used to assess volcanic gas hazards at Icelandic volcanoes where no direct measurements are available. Broadly, the results underline the governing role of sulfide saturation during melting and magma differentiation in controlling the eruptible S contents of Icelandic magmas.


Introduction
Volcanic eruptions release vast quantities of volatiles (H 2 O, CO 2 , S, Cl, and F) into the atmosphere.On geological timescales, volcanic outgassing regulates the composition of Earth's atmosphere (Gaillard et al., 2021).On short timescales, volcanic gases pose one of the principal hazards of eruptions.When gas plumes of large eruptions reach the stratosphere, volcanic sulfur (mainly released as SO 2 ) forms aerosols that block sunlight and cool the surface, with the potential to cause short-term hemispheric or even global climate perturbations (Robock, 2000).In the troposphere, SO 2 reacts with water and oxidizes to form sulfuric acid and sulfate (SO 4 2 ) aerosols, with both the gas and aerosol contributing to volcanic pollution that is one of the main environmental and health hazards for eruptions of all sizes (Carlsen et al., 2021;Gíslason et al., 2015;Ilyinskaya et al., 2017;Oppenheimer et al., 2011;Schmidt et al., 2011;Stefánsson et al., 2017;Stewart et al., 2022).Assessing the sulfur release potential of future eruptions based on knowledge from past eruptions is a critical component of both long and short term hazard assessments.
Volatiles in volcanic gas plumes can be measured by satellite instruments (e.g., Carboni et al., 2019;Carn et al., 2017;Esse et al., 2023;Schmidt et al., 2015) or by ground-based measurements (e.g., Donovan et al., 2023;Ilyinskaya et al., 2017;Kern et al., 2020;Pfeffer et al., 2018Pfeffer et al., , 2024;;Scott et al., 2023;Vignelles et al., 2016).However, direct measurements are only available for a relatively small selection of post-1970s eruptions globally and in Iceland (Oppenheimer et al., 2011).Volatile emissions from past eruptions can be gauged retroactively by indirect means.The "petrological method" (Anderson, 1974;Devine et al., 1984) utilizes the capacity of crystalhosted silicate glass melt inclusions (MIs) to retain high volatile concentrations, which in the right circumstances indicate the undegassed, pre-eruptive state of their host magmas (Lowenstern, 1995).The key advantage of the petrological method is that volatile emissions from any past eruption where glassy MIs are preserved, and the erupted mass of magma is known, can be estimated.
However, the volatile emission potentials of most of the 33 active volcanoes in Iceland (www.icelandicvolcanos.is)remain undetermined, leaving a large gap in our understanding of regional volatile systematics and limiting the ability to prepare for future volcanic degassing hazards.Furthermore, effects of petrological processes that determine pre-eruptive S contents, that is, mantle melting dynamics, mantle source heterogeneity and S concentration, and timing of sulfide saturation during crustal magma evolution are incompletely understood for Icelandic magmas.For example, chalcophile and platinum group element systematics of Icelandic lavas (Momme et al., 2003), along with partial melting models (Ding & Dasgupta, 2018) suggest that mantle melting at rift zones is likely sulfide-undersaturated, but that rift melts reach sulfide saturation during crustal evolution (Halldórsson et al., 2018;Momme et al., 2003;Ranta et al., 2022).However, these results are based on limited data sets and are ambiguous with regards to the timing of sulfide saturation and the topology of the sulfide saturation surface ("sulfide content at sulfide saturation" = S 2-CSS; e.g., Ding et al., 2023;Hughes et al., 2023;Nash et al., 2019;O'Neill, 2021;O'Neill & Mavrogenes, 2002;Shima & Naldrett, 1975;Smythe et al., 2017) in the P-T-X space.The S 2-CSS, in turn, essentially controls the theoretical maximum eruptible S concentrations in melts where dissolved S species are dominated by S 2 over S 6+ (Ding et al., 2023;Hughes et al., 2023) such as those in Iceland (Ranta et al., 2022).
Here, we present a database of MI and matrix glass (MG) volatile concentrations from 84 Icelandic eruptions, compiled from published data as well as new data from eight eruptions.We use the database to review the pre-eruptive sulfur systematics of Icelandic volcanoes, from mantle melting to crustal sulfide melt immiscibility.We estimate eruption-specific maximum S emission potentials (in ppm) for 73 eruptions where MI S data are available, and calculate total S emission estimates (in Tg) for a subset of 45 eruptions with known volume, using the petrological method.Our results can be used to inform forecasts of hazards for future eruptions.For example, an estimate of the S emission potential for a given volcano/volcanic scenario combined with an estimate of effusion rates and injection heights can be used as input data for gas dispersal models in pre-eruptive or early eruptive phases (Barsotti, 2020), most importantly at volcanoes where no prior direct measurements are available.

Geological Setting
Iceland is an atypically volcanically active subaerial segment of the slow-spreading Mid-Atlantic Ridge, an ocean island hotspot located at the divergent plate boundary between North America and Eurasia (Figure 1).Excess magmatism is caused by interaction between the spreading center and a hot mantle upwelling, the Iceland mantle plume (e.g., Ito et al., 1996;Schilling, 1973;Shorttle et al., 2014;Wolfe et al., 1997).At present, 33 volcanic systems in Iceland are considered active, with a collective average eruption interval of 4-5 years (Thordarson & Larsen, 2007;www.icelandicvolcanos.is).The activity is concentrated in several rift segments, here divided into the Eastern, Northern and Western Rift Zones (ERZ, NRZ and WRZ) and the Reykjanes Peninsula (RP), the Mid-Iceland Belt, the propagating rift South Iceland Volcanic Zone (SIVZ) and two off-rift zones, Snaefellsnes Volcanic Zone (SNVZ) and Öraefajökull Volcanic Belt (ÖVB) (Figure 1).Holocene eruption types vary from common effusive basaltic fissure eruptions (which may have explosive phases) to less frequent explosive silicic eruptions (Thordarson & Larsen, 2007).Many volcanoes are located entirely or partly beneath glaciers; subglacial eruptions are often associated with a major explosive component (Larsen, 2002).Silicic volcanics constitute about 10% of the exposed bedrock (Jónasson, 2007).

Samples and Electron Probe Microanalysis
New major element and S and Cl concentration data for glassy MIs (hosted in olivine, plagioclase and clinopyroxene) and matrix glasses are presented here from eight eruptions from five volcanoes (yellow triangles in Figure 1).Seven of these, from the SIVZ, SNVZ, and ÖVB were targeted specifically to patch a shortage of S data from these volcanic zones.(a) Fjallsendahraun (∼1362 CE; also known as Frambruni) belongs to the Bárðarbunga volcanic system (Sigmarsson & Halldórsson, 2015;Svavarsdóttir et al., 2017;Þórðardóttir, 2020) and is one of the largest Holocene fissure eruptions in Iceland.It is a tholeiitic lava flow with an estimated volume of ∼4 km 3 (Thordarson & Larsen, 2007).Two scoria samples (FJ19-5a, FJ19-6a) from two different craters from the eruptive fissure were used in this study.(b) An unnamed Holocene trachybasalt lava at Djúpalónssandur in the SNVZ from the Snaefellsjökull volcano (Hardarson, 1993).The sample (DJUP1) is scoriaceous lava crust.(c) Trachybasalt/ basaltic trachyandesite Eldfell eruption that started 23 January 1973 in the island of Heimaey in the SIVZ (Furman et al., 1991;Jakobsson et al., 1973).Two tephra samples were used; ELD-2 is from the first day of eruption (January 23) and ELD-1 was collected from a "vagabond" crater (collapsed piece of a crater transported by a lava flow) formed during a later stage of the eruption.(d) A scoria sample (BFR-1) was collected from one of the craters of the basaltic Grábrókarhraun lava near the village of Bifröst, which is the easternmost Holocene eruption of the SNVZ and thought to originate from the Ljósufjöll volcanic system (Hardarson, 1993).Sulfur and Cl data are also presented for MIs from tephra of the ∼4 ka Berserkjahraun eruption of Ljósufjöll, for which major element data were previously reported by Kahl et al. (2021).(e) Two basaltic (ÖRA-1 and ÖRA-3) and one basaltic trachyandesite lava (ÖRA-2) of unknown age from Öraefajökull were collected from separate stratigraphic units by the Svínafellsjökull glacier (J.Helgason, 2000).The lavas have a thick (>1 cm) glassy crust and likely erupted subglacially.Major element, S and Cl analyses were performed on hand-picked glasses and crystals by electron probe microanalysis (EPMA) with the JEOL JXA-8230 SuperProbe at the Institute of Earth Sciences, University of Iceland, which is equipped with five wavelength-dispersive spectrometers.The analytical settings were identical to those described in Caracciolo et al. (2020) and Ranta et al. (2022).Precision and accuracy of the EPMA measurements and instrumental drift were monitored by analyzing the basaltic glass standard VG-A99 and a Lipari obsidian at the beginning and end of each session.The VG-A99 values for S (137 ± 27 ppm) and Cl (196 ± 22 ppm; both 1σ, n = 12) agree within mutual uncertainty at 1σ level with published values of S = 96 ± 63 ppm (Thornber et al., 2002), S = 160 ± 40 ppm (Jégo & Dasgupta, 2013) and Cl = 240 ± 20 ppm (Dixon & Clague, 2001) (Table S7 in Supporting Information S2).

Compilation of Published Data
The Iceland Melt Inclusion Catalog (IMIC; Where available, associated MG analyses of the host lavas are also included in the database.

Data Type, Quality, and Filtering
Volatile concentrations in volcanic glasses are typically determined using microanalytical techniques that enable the analysis of MIs as small as ∼5 μm in diameter.The most common analytical techniques are secondary ion mass spectrometry (SIMS; H 2 O, CO 2 , F, S, Cl), electron probe microanalysis (EPMA; F, S, Cl), Fourier transform infrared spectroscopy (FTIR; H 2 O, CO 2 ) and Raman spectroscopy (H 2 O, CO 2 ).Data from all these techniques are represented in IMIC, and the methods used by each study are specified in Table 1.Volatile abundances determined by the "volatiles by difference" method are not included due to typically low accuracy (Hughes et al., 2019).To retain the data integrity, no filtering was done, and published data are reported in IMIC as in the original publications.However, volatile concentrations reported for MIs that have been rehomogenized under atmospheric pressure are considered unreliable due to possible volatile loss during reheating (Caracciolo et al., 2020(Caracciolo et al., , 2022;;Venugopal et al., 2020) and are thus not considered further.

Post-Entrapment Process (PEP) Correction
Major and trace element and volatile concentrations in MIs can be modified by post-entrapment crystallization and/or diffusive re-equilibration (e.g., Danyushevsky et al., 2000).Post-Entrapment Process (PEP) can significantly affect the concentrations of compatible elements in the residual melt (e.g., MgO; Figure S1a in Supporting Information S1) but the effect is minor (<20%) for incompatible elements, including S (Figure S1b in Supporting Information S1).Some, but not all, studies perform PEP corrections to reconstruct original MI compositions.However, the correction schemes vary between studies and depend on the host-crystal.When possible, both the measured concentration data (i.e., prior to PEP corrections) as well as the reported PEP-corrected values are reported in IMIC.Thus, the catalog user has access to unprocessed data and can opt to make PEP corrections using the method of their choice.
MIs generally preserve high volatile concentrations that are used as an approximation of the pre-eruptive state (C 0 ) of a magma.In turn, matrix glasses -which most commonly either quenched rapidly after leaving the vent (tephra) or at the surface of a lava flow (lava crust)-provide an estimate of the post-eruptive volatile contents (C 1 ).Because the exact crystallinity of lavas is not reported for the majority of the eruptions in IMIC, we opt to use X cryst = 0 for all ΔC calculations, leading to a bias toward high calculated ΔC values.The X cryst value ranges from <0.05 for common, relatively aphyric basaltic eruptions (e.g., Holuhraun 2014Holuhraun -2015;;Halldórsson et al., 2018) to up to >0.3 in plagioclase ultraphyric basaltic eruptions (Hansen & Grönvold, 2000).For sulfur, S 0 is chosen differently for melt compositions above and below MgO = 4 wt.% because of the pronounced effect of sulfide saturation on S concentrations in melts below that approximate threshold (Ranta et al., 2022;Section 5.2).For melt compositions with MgO > 4 wt.%, S 0 is taken as the highest measured MI concentration (after discarding outliers, always outside 75% confidence interval) following empirical evidence that this approach leads to close correspondence with observed emissions (Bali et al., 2018;Pfeffer et al., 2018Pfeffer et al., , 2024)).For melts with MgO < 4 wt.%(low-MgO basalts to rhyolites), where maximum S contents tend to be associated with higher-MgO MIs hosted in more primitive crystal cargo relative to the carrier melt, S 0 is instead chosen as the mean MI S content.For all melt compositions, S 1 is chosen as the minimum S concentration measured in matrix glasses.The minimum concentration is used due to an analytical bias toward higher concentrations than bulk erupted material, because analyses are generally performed on the fastest quenched (and hence least degassed) glass.Because S concentrations in tephra glasses from a single eruption may vary depending on how fast the sample was quenched (Liu et al., 2018;Sigmarsson et al., 2013), it is useful to differentiate between the apparent sulfur release at the vent and the maximum sulfur emission potential which reflects degassing of all pre-eruptive sulfur (i.e., S 1 = 0 in Equation 2).
If the total dense rock estimate volume, V DRE (m 3 ) (equivalent of vesicle-free solidified lava), and density, ρ (kg/ m 3 ), of the erupted products are known, the total mass of sulfur release, m S (kg), can be calculated as The densities (ρ) of the solidified lava and the melt are assumed to be 2,710 kg/m 3 for basalts, 2,600 kg/m 3 for intermediate rocks and 2,470 kg/m 3 for silicic rocks (Hartley & Maclennan, 2018;Sharma et al., 2008).For tephras with no reported V DRE , we make the approximation  S7 of Ranta et al., 2022).All MI and matrix glass data are plotted raw and unnormalized.
Geochemistry, Geophysics, Geosystems where V is the estimated volume of tephra deposits and/or the lava flow volume and b is a scaling factor that accounts for vesicularity (Table S2).

Caveats and Sources of Error of the Petrological Method
Several caveats apply when using the petrological method: 1. MIs do not always accurately record the volatile composition of the pre-eruptive melt.Volatile concentrations in MIs can strictly represent the pre-eruptive state of the melt only if they were trapped from the carrier melt shortly before eruption and were quenched rapidly enough to prevent diffusive modifications or postentrapment crystallization.More commonly, MIs formed during an earlier stage of the magma's history.
Then, decoupling of MI and host melt after entrapment partly shields the MI from subsequent magmatic processes (e.g., fractional crystallization, magma recharge, gas fluxing, vapor or sulfide saturation) that may have caused the host magma to either lose or gain volatiles.The common observation that MIs and their host minerals tend to be more primitive and chemically heterogeneous than the carrier melts (e.g., Maclennan, 2008a;Thomson & Maclennan, 2013) reflects such decoupling.This problem is partly overcome in our ΔS calculations by the choice of maximum MI S concentration, which is typically found in MIs with MgO contents that closely matches the carrier melt.2. Different volatile species need to be treated separately because of their differing solubilities and diffusion rates in melts (Baker et al., 2005), and additional factors that affect their concentrations in melts and their preservation in MIs, listed below in detail.3. Accurate CO 2 measurements in MIs are complicated by the necessity to account for vapor bubbles and mineral precipitates on the bubble walls, which can host more than 90% of the CO 2 contained in the MI (Hartley et al., 2014;Moore et al., 2015;Rasmussen et al., 2020;Schiavi et al., 2020), and decrepitation, which leads to loss of CO 2 (Maclennan, 2017).Bubbles in MIs may also host S as gas or precipitates (Venugopal et al., 2020), but the effect on total MI S contents may be insignificant (Korneeva et al., 2023;Rasmussen et al., 2020).Bubbles are present in 17%-65% of investigated MIs in eight studies in IMIC where they are reported (Hartley et al., 2014;Hauri et al., 2018;Matthews et al., 2021;Miller et al., 2019;Moune et al., 2007Moune et al., , 2012;;Neave et al., 2014Neave et al., , 2017)).4. Water, due to the high diffusivity of H + , is especially prone to re-equilibration which can lead to both increased and decreased H 2 O contents in the MI (Barth & Plank, 2021;Bucholz et al., 2013;Gaetani et al., 2012;Hauri, 2002).This means that MIs are more likely to have equilibrated with their carrier melts with respect to H 2 O, thus representing true pre-eruptive contents, although the same effect can generate artificially low MI H 2 O concentrations via post-eruptive re-equilibration with a degassed host melt. 5. Concentrations of F, Cl and S-all slow-diffusing elements relative to H (Baker et al., 2005;Freda et al., 2005)-are here assumed to be only affected by post-entrapment crystallization, which leads to an increase of their concentrations in the residual melt.This distillation effect leads to only a small increase in the measured concentrations (<5%) for the typically modest (<15%) reported amounts of post-entrapment crystallization and is always less than 20% for published data (Figure S1b in Supporting Information S1). 6.A caveat of the previous point is that boundary layer effects during MI entrapment can lead to diffusive fractionation that generate higher concentrations of slow-diffusing elements in the MI relative to the original melt (Baker, 2008).In experiments, MIs with S and Cl enrichments of up to 50% and 40%, respectively, have been observed (Baker, 2008).7. Sulfide globules in MIs may buffer the glass S concentrations.Sulfide globules in Icelandic MIs occur (Ranta et al., 2022) but are relatively rare (observed in less than approximately 2% of basaltic MIs according to the collective experience of the authors).These may either signal trapping of melts that were already sulfidesaturated prior to MI entrapment but could also form post-entrapment if the trapped melt reaches S 2-CSS, for example, by post-entrapment crystallization (Audétat et al., 2018) or diffusive Fe-loss (Lerner et al., 2021).In the latter case, the measured MI glass S concentration by microanalytical methods would lead to underestimation of the original melt S content.This effect could be overcome by rehomogenization of MIs (Korneeva et al., 2023), or by bulk LA-ICP-MS measurements of MIs (Rottier & Audétat, 2019) in setups where S can be quantitatively determined.8.The petrological method only takes into account the volatile emission capacity of the erupted melt itself.
Actual volcanic volatile emissions may additionally involve input from magmatic S-bearing volatile phases or external sources of volatiles (see Section 5.2).
Geochemistry, Geophysics, Geosystems 10.1029/2024GC011443 9. The petrological method also neglects post-eruptive outgassing from cooling lava, which is an additional source of volcanic pollution (Donovan et al., 2023;Óskarsson et al., 1984;Pfeffer et al., 2018;Sigmarsson et al., 2020;Simmons et al., 2017;Thordarson et al., 1996).This process disproportionately affects F and Cl, which have initial concentrations far below their melt saturation values.However, as their concentrations in residual melts are driven up by microlite crystallization in cooling lava fields, they partition in a greater proportion to a developing vapor phase, and may even reach saturation after extreme crystal fractionation.Thus, petrological F and Cl emission estimates-often negligible (Figure S4 in Supporting Information S1) -should be taken as minimum values.10.Slow diffusivity of S in silicate melts may cause disequilibrium degassing (Lerner et al., 2021), which may explain commonly observed MG S contents above the vapor solubility limit of melts at surface pressures.The residual sulfur is progressively degassed from a cooling lava field, as seen in progressively waning S concentrations in matrix glasses with increasing distance to the vent and measurable SO 2 degassing of lava fields observed by field measurements (Donovan et al., 2023;Lerner et al., 2021).For this reason, we consider ΔS max to be a more reliable indicator of total SO 2 release.
Minding these caveats, the petrological method is broadly considered to be a reliable way to measure pre-eruptive S concentrations (e.g., Devine et al., 1984;Lerner et al., 2021;Wallace & Edmonds, 2011) and has proven to be an accurate means for estimating SO 2 emissions of Icelandic basaltic eruptions (Bali et al., 2018;Pfeffer et al., 2018Pfeffer et al., , 2024)).
Notably, the Eldfell MIs lie outside of the compositional coverage of previously published MI data, being more alkaline and extending from trachybasaltic to basaltic trachyandesitic compositions (the latter exemplified by the ELD-2 tephra glass from 23 January 1973, the first day of the eruption; Figure 2c).These compositional trends echo the bulk rock compositions of the eruption (Furman et al., 1991;Jakobsson et al., 1973;Sigvaldason & Óskarsson, 1976).

Sulfur Concentrations of Icelandic Melt Inclusions
Sulfur concentrations in Icelandic MIs in IMIC vary between 3 and 2,630 ppm (Figure 3).The S concentrations are highest in basaltic MIs with a mean of 1,010 ± 380 ppm, decreasing with increasing SiO 2 to 650 ± 300 ppm in intermediate to 80 ± 60 ppm in silicic MIs (Figure 3a).Most variability is seen for basaltic MIs, where S is lower (∼300-1,300 ppm) in primitive high-MgO (>10 wt.%) basalts and increases with decreasing MgO toward peak S contents (∼1,000-2,630 ppm) at about 5 wt.%MgO (Figure 3b).3a).Subaerially erupted tephra matrix glasses have generally low S concentrations (from below detection to few 100s of ppm) relative to MIs from the same eruption.

Effects of Crustal Magma Evolution on S and S 2-CSS
Fractional crystallization during crustal magma evolution leads to increasing concentrations of dissolved volatiles, which tend to behave as incompatible elements in silicate melts, but only until volatile-bearing fluid, immiscible melt or solid phases are saturated (e.g., Hughes et al., 2023;O'Neill & Mavrogenes, 2002;Shima & Naldrett, 1975).The effects of S degassing and sulfide melt immiscibility on pre-eruptive S concentrations of Icelandic magmas are discussed below.
With respect to vapor saturation, S has an intermediate solubility in basaltic melts, being more soluble than CO 2 but less soluble than H 2 O at a given pressure (Wallace & Edmonds, 2011).Basaltic melts in Iceland are H 2 O-poor (<1.5 wt.%, Figure S2a in Supporting Information S1) but contain enough CO 2 to form a CO 2 -dominated vapor phase already at deep crustal, or even sub-Moho-levels (Hartley et al., 2014;Matthews et al., 2021;Miller et al., 2019).Modeling of decompression degassing of Icelandic basalts using the Sulfur_X model (Ding et al., 2023) shows that although minor amount of S partitions to the CO 2 -rich vapor phase at depth, >99% of S remain dissolved in the melt at typical magma storage pressures (∼1-10 kbar; Figure 4; Baxter et al., 2023;Ding et al., 2023;Neave et al., 2019).Bulk of the dissolved sulfur in Icelandic basaltic melts is degassed at very low pressures of <0.2 kbar (Figure 4; see also Ranta et al., 2023;Sigmundsson et al., 2020).From a volcano monitoring perspective, it is important to note that more H 2 O-rich basaltic andesite melts (such as the Hekla 2000 eruption) can start degassing sulfur at greater pressures of ∼0.8 kbar (Figure 4).Thus, the MI S record, at least for basaltic MIs, should reflect near-undegassed melt contents at the time of entrapment.
However, formation of an immiscible sulfide liquid at the sulfide solubility limit has for long been recognized as an important factor that buffers the S (and chalcophile element) budgets of both arc-related (Jenner et al., 2010;Lee et al., 2012) and mid-ocean ridge basalts (MORBs) (Mathez, 1976), as well as many ocean plateau and hotspot lavas (Labidi et al., 2015;Reekie et al., 2019;Wieser et al., 2020) including Iceland (Momme et al., 2003;Moune et al., 2009;Ranta et al., 2022).In Icelandic melts, a smoking gun evidence for sulfide saturation is provided by decreasing concentrations of Cu-which is incompatible in silicate minerals but strongly compatible in the forming monosulfide liquid (Kiseeva & Wood, 2015)-with decreasing MgO below approximately 6 wt.% MgO in subglacial pillow rim glasses (Ranta et al., 2022).Whether higher-MgO melts reach sulfide saturation is not entirely clear from previous data (Ranta et al., 2022).
The choice of S 2-CSS model does not radically affect the following discussion: for example, the O'Neill (2021) and Li and Zhang (2022) models generally yield similar-shaped S 2-CSS curves at slightly higher and lower S 2-CSS (Figures 5a and 5b), respectively, but within the estimated 1σ uncertainty (∼27%) of the Smythe  et al., 2020; Óskarsson et al., 1994), P = 2 kbar, T calculated with the major element composition-based liquid thermometer of Putirka (2008; Equation 14) and a sulfide composition fixed at Fe/(Fe + Ni + Cu) sulf = 0.65 (Table S8 in Supporting Information S2).
The modeling highlights the dominant effects of FeO and T on the modeled S 2-CSS and the effect of Fe 3+ /Fe T on S T CSS.Decreasing temperature lowers the S 2-CSS, but also causes crystal fractionation and decreasing MgO  Halldórsson et al., 2018) and Hekla 2000 (3.9 kbar, calculated from an estimated magma storage depth of >14 km and mean crustal density of 2,860 kg/m 3 ; Höskuldsson et al., 2007).To illustrate the effect of higher H 2 O-concentrations on S degassing, we modeled hypothetical degassing paths for the Holuhraun melt with H 2 O concentrations of 1.5 and 3 wt.%(two dashed gray curves).Higher H 2 O-concentrations allows the onset of S degassing at greater pressures (cf.Ding et al., 2023).However, water contents above 1.5 wt.%, required for significant degassing to occur above the approximate lower limit of 1 kbar of basaltic magma reservoirs in Iceland (Neave & Putirka, 2017) (Moune et al., 2007), 500 ppm CO 2 (arbitrarily chosen due to lack of data) and 990 ppm S for Hekla 2000 (Moune et al., 2007).Melt temperatures were chosen as 1,170°C for Holuhraun (Halldórsson et al., 2018) and 1,100°C for Hekla (Höskuldsson et al., 2007).Major element compositions used were averages of matrix glasses (Table S2).
(Figure 5a).The opposing effects of FeO and T on the S 2-CSS lead to zig-zagging liquid lines of descent in a FeO versus S 2-CSS diagram for a typical tholeiitic fractionation path (Figure 5b).At early stages of tholeiitic basaltic melt evolution characteristic for Icelandic rift zones, fractional crystallization of olivine ± clinopyroxene ± plagioclase causes melt FeO to increase with decreasing T and decreasing MgO.Following the onset of magnetite crystallization at ∼5 wt.% MgO, both FeO and T decrease causing S 2-CSS to drop (Figures 5a and 5b).2017) model.S T CSS was calculated for fixed Fe 3+ /Fe T of 0.10 and 0.15, representing the approximate variability found in Icelandic subglacial and submarine glasses (Novella et al., 2020;Óskarsson et al., 1994).The S T CSS model takes into account S 6+ /S T calculated according to O'Neill and Mavrogenes (2022); Figure S7 in Supporting Information S1, and the S 6+ CAS (sulfate concentration at anhydrite saturation, not shown) is significantly higher than estimated S 6+ concentrations (Zajacz & Tsay, 2019; see Table S8 in Supporting Information S2 for model parameters).Because the S T CSS and S 2-CSS overlap for Fe 3+ /Fe T = 0.10 (panel a), this S T CSS calculation is omitted in the other panels for clarity.Similarly, only S 2-CSS calculated for Fe 3+ /Fe T = 0.10 is shown, as the difference to S 2-CSS calculated for Fe 3+ /Fe T = 0.15 is small (<5%; Figure S8 in Supporting Information S1).An example S 2-CSS along a fractional crystallization path is shown for a PEPcorrected MI composition (NAL-356-O5-MI1; Ranta, 2022), produced by isobaric (3 kbar) closed-system fractional crystallization in the Petrolog3 software (Danyushevsky & Plechov, 2011).Dashed, solid and dotted lines show the S 2-CSS calculated after Smythe et al. (2017), O'Neill (2021), and Li and Zhang (2022), respectively.The fractional crystallization path illustrates the main features seen in empirical data, but reaches magnetite saturation at slightly higher MgO.Negative effect of P on S 2-CSS means that the S 2-CSS for MIs trapped at greater depths than the modeled 2 kbar may be overestimated, and vice versa.For example, the S 2-CSS calculated at 6 kbar are between 6% and 10% lower than at 2 kbar (Figure S8 in Supporting Information S1).All S 2-CSS and S T CSS models and the indicative 1σ error of the Smythe et al. (2017) model of approximately ±27% (±405 ppm at 1,500 ppm, shown in panel a) were calculated using the open source PySulfSat code (Wieser & Gleeson, 2023).
The S 6+ /S T and S T CSS are strongly dependent on redox conditions: For the lower Fe 3+ /Fe T of 0.10, modeled S 6+ / S T are low (<0.10 at all compositions; Figure S7 in Supporting Information S1) and S 2-CSS is almost equal to S T CSS (Figure 5a).At more oxidizing conditions (Fe 3+ /Fe T = 0.15), modeled S 6+ /S T are high in high-MgO melts (up to 0.40 at MgO = 10 wt.%) but drop with decreasing MgO to less than 0.20 for most MI below 8 wt.% MgO.These values are compatible with low measured S 6+ /S T (0.03 ± 0.11 to 0.09 ± 0.20) in three Icelandic basaltic subglacial glasses (MgO between 5.1 and 8.2 wt.%) with Fe 3+ /Fe T = 0.13-15 (Taracsák et al., 2021).
Importantly, magma evolution leads to a local maximum in both calculated S 2-CSS and S T CSS and observed S concentrations (1,200-2,600 ppm) plotted against MgO, seen in evolved basaltic MIs at about 4-6 wt.% MgO (Figure 5a), coinciding with maximum melt FeO contents (>14 wt.%, Figure S3a in Supporting Information S1).Subsequent melt evolution takes place under sulfide-saturated conditions, leading to systematically decreasing S due to decreasing S  5b).These effects explain why moderate S concentrations are observed in basaltic andesites, and lowest S are systematically found in dacites and rhyolites (Figures 3a and 8a).
That the calculated S 2-CSS of the MIs recreates some of the main features traced by the upper limit of measured MI S concentrations at a given MgO (Figure 8a) suggests that at least the maximum melt S content is limited by S 2-CSS at all MgO contents.However, in general, of MIs with MgO >6 wt.%, most have S concentrations below S 2-CSS and all have S concentrations below S T CSS calculated for Fe 3+ /Fe T = 0.15 (Figure 5d).The local S 2-CSS and S T CSS maximum at about 4-6 wt.% MgO is compatible with a similar MgO threshold for sulfide saturation suggested previously for Icelandic melts (Ranta et al., 2022) and is similar to that observed in many arc basalts globally based on Cu systematics (Zhang et al., 2021;Zhao et al., 2022).The range of possible fractional crystallization paths likely allows both sulfide-saturated and undersaturated high-MgO melts with lower Fe 3+ /Fe T (≈0.10), whereas higher Fe 3+ /Fe T (≈0.15) high-MgO are more likely to be sulfide-undersaturated.Most MIs below 4 wt.%MgO appear to be sulfide-saturated regardless of their model Fe 3+ /Fe T (Figure 5d).

Mantle Source, Partial Melting and Regional Variability
Sulfur contents of primary basaltic melts are essentially controlled by the S content and the extent of melting of their mantle source, and its oxygen fugacity.The mantle beneath Iceland is thought to have an oxygen fugacity close to, or slightly below the FMQ (Novella et al., 2020;Óskarsson et al., 1994;van der Meer et al., 2021).While more oxidized melts have higher S T CSS and are capable of dissolving more S at high MgO (Figure 5a), the measured MI S contents show good agreement with S 2-CSS calculated assuming Fe 3+ /Fe T of 0.10.Thus, Icelandic primitive melts appear to have a have low S 6+ /S T .A number of studies have shown that the mantle beneath Iceland is heterogeneous with respect to the volatile origin (and probably, volatile concentrations) on a regional scale (e.g., Füri et al., 2010;Halldórsson, Barnes, et al., 2016, Halldórsson, Hilton, et al., 2016;Harðardóttir et al., 2018;Hilton et al., 2000;Kurz et al., 1985;Marshall et al., 2022;Matthews et al., 2021;Miller et al., 2019;Nichols et al., 2002;Poreda et al., 1986;Ranta et al., 2022).Given this, and that the degree of mantle melting is largely controlled by the tectonic setting and lithospheric thickness-melting degree being higher near the rift zones, and lower in the off-rift zones (Harðardóttir et al., 2022)-regional variations are expected to be expressed in MI volatile concentrations.For example, 2% partial melting of a hypothetical mantle component A with 200 ppm H 2 O (assumed to be perfectly incompatible) will yield a primary basaltic melt with 1 wt.%H 2 O, whereas 20% partial melting produces a melt with 0.1 wt.% H 2 O-a factor of 10 less (assuming accumulated fractional melting).On the other hand, if mantle component B has 1,000 ppm H 2 O (i.e., 5× mantle component A), any equal degrees of melting of A and B will lead to a 5× difference in the H 2 O contents of their respective primary melts.These differences are well within the suggested variability in both the degree of melting across different volcanic regions of Iceland (e.g., Koornneef et al., 2012) and the range of volatile concentrations inferred for the enriched and depleted mantle sources of mid-ocean ridge basalts (MORBs; Shimizu et al., 2016).Harðardóttir et al., 2022;Koornneef et al., 2012;Schilling, 1973).However, regional systematics do not appear to govern magmatic S contents, which are similar at both rift and off-rift regions at similar MgO or SiO 2 contents (Figure 3) and do not correlate with La/Yb (Figure 6b; R 2 = 0.08).To investigate which additional process controls melt S concentrations, we inspect the MI S/Cl ratio, which should stay approximately constant during partial melting and fractional crystallization due to the similar incompatibility of both elements in relevant silicate minerals (Baker et al., 2022(Baker et al., , 2023;;Callegaro et al., 2020).Off-rift MIs generally have low S/Cl ratios between 0.5 and 5, lower than either the depleted (S/Cl ≈ 250) or enriched (S/Cl ≈ 7.5) mantle endmembers inferred for the depleted MORB source mantle (D-DMM and E-DMM, respectively; Shimizu et al., 2016;Figure 6c).By contrast, highly depleted MIs in rift basalts (products of high-degree melting) have high S/Cl ratios (>100) approaching the S/Cl signature of D-DMM (Figure 6c).The decoupling of S and Cl systematics suggests that S concentrations in primary melts are buffered by sulfide saturation during mantle melting, in accord with previous work (Ding & Dasgupta, 2018;Momme et al., 2003).
We modeled the effect of partial mantle melting on S and Cl concentrations and S/Cl ratios of primary melts under the assumption that sulfur in the melting region is hosted exclusively in sulfides (Ding & Dasgupta, 2018).Two mantle endmember compositions were used to simulate depleted and enriched mantle components under Iceland (e.g., Gurenko & Chaussidon, 1995;Hanan et al., 2000;Koornneef et al., 2012;Shorttle & Maclennan, 2011) which were assumed to have S, Cl, La, and Yb concentrations of the D-DMM and E-DMM (Workman & Hart, 2005;Shimizu et al., 2016; Figures 6c and 7).The melting model of Lee et al. (2012) implemented in the open source Python3 tool PySulfSat (Wieser & Gleeson, 2023) was used for the calculations.Full model parameters are given in Table S8 in Supporting Information S2.
Similar to the results of Ding and Dasgupta (2018), our modeling indicates that sulfide saturation persists under low-to moderate-degree mantle melting (F < 0.08-0.16).Thus, off-rift melting (typically F < 0.1; Harðardóttir et al., 2022) likely occurs under sulfide-saturated conditions, whereas higher melt degrees under rift zones (typically F > 0.16, Harðardóttir et al., 2022) are generally sufficient to completely exhaust sulfides from the mantle source (Figures 7a and 7b).During sulfide-saturated melting (off-rift), the S   (Shimizu et al., 2016) of the depleted mantle source of mid-ocean ridge basalts (DMM) under initially sulfide-saturated conditions (see Figure 7).The resulting fractionation of the S/Cl ratio suggests that source heterogeneity and variable degree of melting could explain the bulk of the observed melt inclusion S/Cl variability.Model parameters are given in Table S8 in Supporting Information S2.Error bars indicate average reproducibility for S (2σ = ±118 ppm; see Figure 3) and Cl (2σ = ±42 ppm).The latter is based on average reproducibility reported for VG-2 (USNM 111240-52) and VG-A99 (USNM 113498-1) reference glasses (similar Cl concentrations of 200-300 ppm) by Dixon and Clague (2001), Ranta et al. (2022) and this study.
concentrations of accumulated fractional melts are approximately constant (Figure 7a).With increasing meltdegree, Cl concentrations gradually decrease, causing the S/Cl ratio to increase until sulfide originally present in the mantle source is exhausted completely.At this point, the S/Cl ratio of primitive melts approaches that of the mantle source and remains approximately constant during further melting (Figure 7b).Thus, the variable S/Cl ratios of high-MgO MIs reflect both the degree of melting and the S and Cl composition of the source mantle.
The S/Cl variability of Icelandic MIs can be explained by melting of D-DMM and E-DMM with S and Cl concentrations of 100 and 0.4 ppm (D-DMM) and 165 and 22 ppm (E-DMM), respectively (Shimizu et al., 2016).These mantle concentrations are within previously estimated peridotite S concentrations of 100-150 ppm for Iceland (Ding & Dasgupta, 2018) and S of 50-200 ppm and Cl of 3-15 ppm for the mantle source(s) of the Holuhraun 2014-2015 eruption (Bali et al., 2018).For the used model parameters, the lowest S/Cl ratios are produced by low-degree melting (e.g., S/Cl = 1 at F = 2%) of E-DMM.The majority of S/Cl and La/Yb ratios of the most enriched off-rift lavas can be matched by 2%-7% melting of E-DMM (Figure 7c).By contrast, high S/Cl of up to 250 is produced by high-degree melting (F > 10%) of the D-DMM component at sulfide-undersaturated conditions, which aligns well with highly depleted MIs found in rift lavas (e.g., Bali et al., 2018).MIs that fall below the partial melting curves in the S/Cl versus La/Yb space could be explained by binary mixing of enriched and depleted melts (mixing curve for 2% E-DMM melt and 20% D-DMM is shown for illustrative purposes in Figure 7c).Mixing of heterogeneous melts derived from geochemically enriched and depleted mantle components-taking place at multiple levels in the magma plumbing system concurrent with fractional crystallization-is a common and well-documented process in Iceland (e.g., Halldórsson et al., 2022;Maclennan, 2008aMaclennan, , 2008b;;Matthews et al., 2017;Shorttle & Maclennan, 2011).In MIs, these processes are expressed as a diversity of trace element and radiogenic isotope signatures found even in single eruptions (e.g., Maclennan, 2008b).

Magmatic Volatile Phases and Excess Degassing
An inherent assumption of the petrological method is that volatiles that degas during an eruption are solely derived by exsolution of dissolved S from the erupted lava volume.This assumption does not account for volatiles that may coexist with magmas in a separate magmatic volatile phase (MVP), which may be an accumulative product of deep degassing and/or fluids derived from a larger volume of intrusive magma that has reached volatile saturation in the magma storage region (Edmonds & Woods, 2018;Hartley et al., 2014;Ranta et al., 2021;Shinohara, 2008).Thus, the petrological estimate (only accounting for the dissolved volatiles) may underestimate the total volatile emissions (dissolved volatiles + MVP).One of the main factors controlling the volatile saturation state of a magma, and thereby, the development of a MVP, is the pre-eruptive magma storage pressure (Figure 4).
Underestimation of S emissions by the petrological method has been frequently observed at arc volcanoes, something known as the "excess degassing problem" (e.g., Andres et al., 1991;Gerlach & McGee, 1994;Sharma et al., 2004;Shinohara, 2008;Su et al., 2016).To account for existence of a separated MVP, a modified petrological method was proposed by Scaillet and Pichavant (2003) and Scaillet et al. (2003), where 1-5 wt.% of free gas phase at depth was added to the petrological estimate.However, such correction is not needed for basaltic eruptions in Iceland, as the amount of MVP in these magma reservoirs is likely to be significantly lower due to their low H 2 O contents (0.05-1.5 wt.%; Nichols et al., 2002;Figure S2a in Supporting Information S1) relative to arc magmas (>2 wt.%; Plank et al., 2013).At typical Icelandic magma reservoir pressures of 1-10 kbar (Baxter et al., 2023;Neave & Putirka, 2017), virtually all S (and H 2 O, Cl, and F) in the relatively dry basaltic melts remains dissolved in the melt at the observed concentration range (e.g., Wallace & Edmonds, 2011; Figure 4).Thus, the petrological method is suitable for calculating S emissions of basaltic Icelandic volcanoes.
By contrast, it is likely that relatively long-lived silicic magma reservoirs in Iceland develop MVPs (Ranta et al., 2021).For example, the silicic MIs from Hekla have very high H 2 O contents (up to 6.5 wt.%; Portnyagin et al., 2012).Thus, silicic magmas in Iceland could develop S-bearing MVPs in analog with similarly H 2 O-rich arc eruptions (Scaillet et al., 2003), potentially providing significant sources of excess degassing that are not accounted for by the petrological method.However, as no silicic eruptions have occurred during recent history in Iceland, S degassing (and that of other volatiles) accompanying silicic volcanism in Iceland remains poorly constrained.An interesting observation in that regard is that matrix glasses in the silicic Askja 1875 CE eruption have higher S content than the corresponding MIs (Figure 3a), which could indicate an enrichment in S caused by influx of deeper-derived gases from mafic intrusions or recharge melts (Ö.Helgason et al., 1992;Watkins et al., 2017).An important task for future research is to better understand both volcanic S and halogen (particularly F) pollution hazards of silicic eruptions in Iceland.
Petrological estimates of SO 2 have further been shown to agree reasonably well with satellite-based measurements using the Total Ozone Mapping Instrument (TOMS; Sharma et al., 2004) for the Krafla 1984 (0.64 ± 0.19 vs. 0.4 + 0.04-0.15Mt SO 2 ) and Hekla 2000 (0.48 ± 0.14 vs. 0.10 ± 0.05 Mt SO 2 ) eruptions (Sharma et al., 2004).By contrast, similar comparisons at arc volcanoes sometimes show orders of magnitude higher TOMS SO 2 emissions relative to petrological estimates (Sharma et al., 2004).The close match between direct  2012) was implemented in PySulfSat (Wieser & Gleeson, 2023) to emulate fractional melting of a depleted (D-DMM) and mildly enriched (E-DMM) source lithology previously suggested for mid-ocean ridge basalts (Shimizu et al., 2016).The teal and amber model curves in panels (b) and (c) track the melt concentrations of accumulated fractional melts.The dotted pink curve in (c) shows binary mixing of a 2% E-DMM melt and a 20% D-DMM melt.The good match between observed and modeled S/Cl versus La/Yb trends suggests that melt inclusion (MI) S/Cl ratios reflect various degrees of melting of the D-DMM and E-DMM mantle components and subsequent melt mixing.Model parameters are reported in Table S8 in Supporting Information S2.Only MIs with MgO >5 wt.% are plotted in panel (c) to minimize the effect of sulfide melt immiscibility during crustal magma evolution.
Geochemistry, Geophysics, Geosystems 10.1029/2024GC011443 SO 2 measurements with the petrological method for hitherto observed Icelandic fissure eruptions provides an empirical validation of the petrological method as an accurate means to estimate volcanic SO 2 emissions of basaltic fissure eruptions in Iceland.

Sulfur Emission Potentials of Icelandic Eruptions
The sulfur emission potentials (ΔS max , ΔS vent ) for 73 Icelandic eruptions, calculated with the petrological method, are summarized in Table 2 and Figures 8 and 9.An extensive summary including calculation details, mean major element compositions and petrological estimates for the other volatiles is available in Table S2 and Figure S10 in Supporting Information S1.In the discussion and the figures, we opt to use ΔS max , which includes both vent degassing as well as degassing from the cooling lava field, and has been empirically shown to be close to total S emissions for Icelandic eruptions (Bali et al., 2018;Pfeffer et al., 2024).

Largest S Emitters in Iceland
The largest pre-eruptive sulfur concentrations in Iceland (ΔS max > 1,200 ppm) are systematically found in evolved basalts (MgO between 4 and 8 wt.%; Figure 8a) irrespective of volcanic zone (Figure 9b).This is also the most common compositional range of Icelandic eruptions (Figures 2a and 8a), including the most frequently erupting Holocene volcanoes Grímsvötn, Bárðarbunga and Katla and the largest Holocene eruptions by volume (Þjórsárhraun, Eldgjá 939 CE and Laki 1783 CE; Table 2).Basaltic eruptions at the higher end of this range (6-8 wt.% MgO), characterizing, for example, Bárðarbunga and RP eruptions, have somewhat lower ΔS max relative to the more evolved (4-6 wt.% MgO) eruptions typical of Katla and Grímsvötn (Figure 8a).
Exceptionally high sulfur emission potentials (ΔS max > 1,800 ppm) are seen in five eruptions, all within a narrow range in MgO between 4.5 and 6 wt.% (Figure 8a).This coincides with and even exceeds the apparent sulfide solubility peak in Icelandic melts (Figure 5a; Ranta et al., 2022).Four of these eruptions, Eldgjá 939 CE, Hekla 1913CE, Surtsey 1963-1967CE and Fimmvörðuháls 2010 CE, are located in the SIVZ and have a mildly alkaline composition (SiO 2 < 48 wt.%, TAS > 3 wt.%; Figure 9a), although alkalinity itself does not correlate with S 2-CSS (Figure S6 in Supporting Information S1).Thus, up to 50% higher S emissions may be expected for basaltic eruptions in the SIVZ relative to average basaltic rift zone eruptions of similar volume.
Compared to silicic arc volcanoes (Oppenheimer et al., 2011;Scaillet et al., 2003), Icelandic eruptions tend to have higher ΔS (Figure 10), and may thus emit more sulfur per unit mass of erupted material.This is largely a function of the dominant eruption composition in Iceland (relatively evolved basalts) coinciding with a maximum in dissolved S contents.However, basaltic arc magmas can have even higher S contents (several thousands of ppm) due to their generally more oxidized nature (Muth & Wallace, 2022;Rasmussen et al., 2022;Wallace & Edmonds, 2011).Notably, the total S emissions and emission rates of basaltic eruptions are principally controlled by the masses and mass fluxes of erupted material-which both vary by 3-5 orders of magnitude-rather than ΔS, which only varies by a factor of ∼2-3× (Figure 10).
The SO 2 hazard potential of volcanoes that consistently erupt lavas with similar compositions-like Bárðarbunga, Grímsvötn, Krýsuvík, and Brennisteinsföll (Figures 8a and 9b)-can be forecasted with less uncertainty based on our ΔS estimates, with uncertainty remaining due to factors including lava effusion rate and weather conditions.However, for volcanoes like Hekla that produce silicic, intermediate and evolved basaltic eruptions spanning the whole range of Icelandic S emission potentials, the pre-eruptive S content remains a further source of uncertainty.This uncertainty can be addressed using multiple scenarios before an eruption composition is known.

Knowledge Gaps and Suggestions for Future Research
Some large gaps remain in the regional coverage of available volatile MI data from Iceland.Of the large active volcanic systems in Iceland, there is a notable lack of data from Langjökull and Hofsjökull.Katla, one of the most productive volcanoes during the Holocene (Thordarson & Larsen, 2007), has only 16 entries in the database, and the apparent exceptionally high MI S contents of the Eldgjá 939 CE eruption of Katla (Thordarson et al., 2001) warrant more research.There is also a lack of S data from Torfajökull, the largest producer of silicic magmas among active volcanoes in Iceland (Gunnarsson et al., 1998).These would all be desirable targets for future MI studies.That the ΔS max of some SIVZ melts exceed the S 2-CSS by up to 50% is notable.This may suggest that SIVZ magmas are more oxidized than remaining Icelandic magmas and have higher proportion of dissolved S 6+ relative to S 2 , which would increase the total sulfur solubility (e.g., Hughes et al., 2023;Nash et al., 2019;Wallace & Edmonds, 2011).Spinel-olivine oxybarometry of two 500-700 ka ankaramite eruptions of the Eyjafjallajökull volcano (∆FMQ = +0-0.5;Nikkola et al., 2019), and Fe 3+ /Fe T of up to 0.155 in tephra glass of the Surtsey 1963-1967 CE eruption (Schipper & Moussallam, 2017) suggest that SIVZ melts are slightly more oxidized than both  S2). e Likely source craters of the ∼8.6 ka Þjórsárhraun lava (Caracciolo et al., 2020).
rift basalts and the off-rift SNVZ magmas (van der Meer et al., 2021).However, based on S T CSS calculations, even higher melt Fe 3+ /Fe T of 0.17-0.19are required for the SIVZ melts to attain maximum observed ΔS max values (Figure S9 in Supporting Information S1).To better resolve the issue, more MI data from basaltic eruptions in the SIVZ, along with S 6+ /S T and/or Fe 3+ /Fe T determinations of Icelandic glasses are needed, accompanied by directly measured SO 2 emissions of future eruptions.

Conclusions
This study presents an overview of volcanic S emissions of past Icelandic eruptions and the magmatic processes that control pre-eruptive S concentrations in Icelandic magmas.The results can be used for building S emission scenarios of future eruptions.The main conclusions are: 1. Compiled MI S concentrations of Icelandic eruptions show values between 10 and 2,630 ppm and vary systematically with melt MgO and SiO 2 contents.Primitive basalts (MgO > 8 wt.%) have moderate S of 400- 1,500 ppm.Sulfur contents increase with decreasing MgO to a peak in evolved basalts (900-2,600 ppm at MgO = 4-6 wt.%; Figure 8a).Moderate S concentrations and ΔS are observed in basaltic andesites, and lowest S and ΔS values are found in dacites and rhyolites.2. The S/Cl ratios of Icelandic MIs are highly variable (0.5-300) but anticorrelate with indices of melt degree and source enrichment (e.g., La/ Yb).A partial melting model accounting for sulfide-saturated mantle melting successfully recreates the observed S/Cl versus La/Yb trend.The data are compatible with the depleted and enriched components of the Iceland mantle having S and Cl concentrations similar to previously suggested endmembers of the MORB-source mantle (D-DMM and E-DMM).3. Modeled S 2-CSS and measured MI S data show similar patterns against indices of magma differentiation (FeO, MgO, SiO 2 ), demonstrating that pre-eruptive sulfur concentrations and hence, eruptible S content of evolving Icelandic magmas in the crust are limited by the S 2-CSS.4. Due to the buffering effect of mantle sulfides during melting and the limiting effect of S 2-CSS during crustal evolution, mantle S heterogeneity and melt degree only have a minor effect on the S contents of primary melts, in contrast to other volatiles.5. Sulfur emission potentials (pre-eruptive S contents) are highly predictable for eruptions of volcanoes that commonly erupt melts of similar major element compositions, like Grímsvötn, Bárðarbunga and the RP volcanic systems.For volcanoes with a compositionally variable repertoire, like Hekla, the S emissions potentials are variable and multiple gas hazard scenarios must be constructed for each possible eruption type.6.Total S emissions and S emission rates of basalts-important parameters for volcanic hazard assessmentsdepend more on eruptive volume and mass flux (which vary by ∼3-5 orders of magnitude) than on preeruptive S contents (which only vary by a factor of ∼2-3×).7. Volcanic S emissions of silicic eruptions (which have not been yet directly measured in Iceland) are likely to be vastly underestimated by the petrological method, because it does not account for excess S-bearing magmatic volatile phases.8. Three of the four highest sulfur emission potentials were measured in evolved basaltic eruptions from the SIVZ volcanic systems Hekla, Katla and Surtsey.This indicates that SIVZ melts may have up to 50% higher S emission potential than melts of similar compositions in other volcanic zones.However, due to scarce data for these locations, more work is needed to confirm this finding and to understand what causes the SIVZ MIs to attain S contents above apparent S 2-CSS limits.9. Icelandic eruptions tend to have high S emission potentials because Icelandic melts tend to erupt at compositions where S 2-CSS is at its maximum.(Oppenheimer et al., 2011;Scaillet et al., 2003 and references therein;Kern et al., 2020).The comparison demonstrates that Icelandic eruptions typically have high sulfur emission potentials (ΔS max , dashed gray lines).

Figure 2 .
Figure 2. Overview of the Iceland melt inclusion catalog.(a) SiO 2 versus MgO.(b) Total alkali (Na 2 O + K 2 O) versus SiO 2 (TAS).The A and B panels show published MI data (matrix glasses, embayments, and pillow glasses are omitted).There is a broad overlap of MIs with whole-rock data from Iceland from the compilation of Harðardóttir et al. (2022; gray dots in (b)).Histograms in (a) and (b) highlight that the distribution of the data is heavily skewed toward basaltic compositions.(c) A cropped TAS diagram with a detailed overview of new data.MI = melt inclusion.MG = matrix glass.All MI and MG oxide data are plotted raw and unnormalized.

Figure 4 .
Figure 4. Sulfur degassing.Modeled closed-system sulfur degassing paths are shown for the basaltic Holuhraun 2014-2015CE eruption (black solid line), which was a typical rift zone fissure eruption, and the basaltic andesite Hekla 2000 CE eruption (red solid line), a common composition for the Hekla volcano.Assuming pre-eruptive melt H 2 O and S contents inferred from melt inclusions (MIs), the model predicts insignificant S degassing at the estimated magma storage pressures for both Holuhraun (2.3 ± 1.4 kbar;Halldórsson et al., 2018) and Hekla 2000 (3.9 kbar, calculated from an estimated magma storage depth of >14 km and mean crustal density of 2,860 kg/m 3 ;Höskuldsson et al., 2007).To illustrate the effect of higher H 2 O-concentrations on S degassing, we modeled hypothetical degassing paths for the Holuhraun melt with H 2 O concentrations of 1.5 and 3 wt.%(two dashed gray curves).Higher H 2 O-concentrations allows the onset of S degassing at greater pressures (cf.Ding et al., 2023).However, water contents above 1.5 wt.%, required for significant degassing to occur above the approximate lower limit of 1 kbar of basaltic magma reservoirs in Iceland(Neave & Putirka, 2017), are not seen in Icelandic basaltic MIs (FigureS2ain Supporting Information S1).This implies that S concentrations of basaltic MIs in Iceland have not been affected by degassing.Degassing paths were calculated with the open source COHS-degassing model Sulfur_X(Ding et al., 2023) using the COH model ofNewman and Lowenstern (2002), S speciation model of O'Neill and Mavrogenes (2022), oxygen fugacity buffered at ∆FMQ = 0 and crystallization turned off.Pre-eruptive volatile concentrations were chosen as 0.61 wt.% H 2 O, 1,700 ppm CO 2 and 1,700 ppm S for Holuhraun(Bali et al., 2018)  and 2.4 wt.% H 2 O(Moune et al., 2007), 500 ppm CO 2 (arbitrarily chosen due to lack of data) and 990 ppm S for Hekla 2000(Moune et al., 2007).Melt temperatures were chosen as 1,170°C for Holuhraun(Halldórsson et al., 2018) and 1,100°C for Hekla(Höskuldsson et al., 2007).Major element compositions used were averages of matrix glasses (TableS2).

Figure 5 .
Figure 5. Modeled S 2-CSS (sulfide content at sulfide saturation) and S T CSS (total sulfur content at sulfide saturation) of Icelandic melt inclusions (MIs) versus (a) MgO, (b) FeO T , (c) SiO 2 , and (d) measured S concentrations.Circles, colored after the FeO (a) and (c) or MgO (b) and (d) content, show S 2-CSS for raw MI data, calculated using the Smythe et al. (2017) model.S T CSS was calculated for fixed Fe 3+ /Fe T of 0.10 and 0.15, representing the approximate variability found in Icelandic subglacial and submarine glasses(Novella et al., 2020;Óskarsson et al., 1994).The S T CSS model takes into account S 6+ /S T calculated according to O'Neill and Mavrogenes (2022); FigureS7in Supporting Information S1, and the S 6+ CAS (sulfate concentration at anhydrite saturation, not shown) is significantly higher than estimated S 6+ concentrations(Zajacz & Tsay, 2019; see TableS8in Supporting Information S2 for model parameters).Because the S T CSS and S 2-CSS overlap for Fe 3+ /Fe T = 0.10 (panel a), this S T CSS calculation is omitted in the other panels for clarity.Similarly, only S 2-CSS calculated for Fe 3+ /Fe T = 0.10 is shown, as the difference to S 2-CSS calculated for Fe 3+ /Fe T = 0.15 is small (<5%; FigureS8in Supporting Information S1).An example S 2-CSS along a fractional crystallization path is shown for a PEPcorrected MI composition (NAL-356-O5-MI1;Ranta, 2022), produced by isobaric (3 kbar) closed-system fractional crystallization in the Petrolog3 software(Danyushevsky & Plechov, 2011).Dashed, solid and dotted lines show the S 2-CSS calculated afterSmythe et al. (2017), O'Neill (2021), andLi and Zhang (2022), respectively.The fractional crystallization path illustrates the main features seen in empirical data, but reaches magnetite saturation at slightly higher MgO.Negative effect of P on S 2-CSS means that the S 2-CSS for MIs trapped at greater depths than the modeled 2 kbar may be overestimated, and vice versa.For example, the S 2-CSS calculated at 6 kbar are between 6% and 10% lower than at 2 kbar (FigureS8in Supporting Information S1).All S 2-CSS and S T CSS models and the indicative 1σ error of theSmythe et al. (2017) model of approximately ±27% (±405 ppm at 1,500 ppm, shown in panel a) were calculated using the open source PySulfSat code(Wieser & Gleeson, 2023).

Figure 6 .
Figure 6.Effect of mantle melting and source heterogeneity on (a) Cl versus La/Yb, (b) S versus La/Yb, and (c) S versus Cl systematics.Good correlation of Cl, and poor correlation of S with La/Yb indicates decoupling of S from Cl during the mantle melting process.The teal and amber model curves in panel (c) track fractional mantle melting of depleted (D-DMM) and enriched (E-DMM; amber) endmembers(Shimizu et al., 2016) of the depleted mantle source of mid-ocean ridge basalts (DMM) under initially sulfide-saturated conditions (see Figure7).The resulting fractionation of the S/Cl ratio suggests that source heterogeneity and variable degree of melting could explain the bulk of the observed melt inclusion S/Cl variability.Model parameters are given in TableS8in Supporting Information S2.Error bars indicate average reproducibility for S (2σ = ±118 ppm; see Figure3) and Cl (2σ = ±42 ppm).The latter is based on average reproducibility reported for VG-2 (USNM 111240-52) and VG-A99 (USNM 113498-1) reference glasses (similar Cl concentrations of 200-300 ppm) byDixon and Clague (2001),Ranta et al. (2022) and this study.

Figure 7 .
Figure 7. Effect of sulfide-saturated partial mantle melting on (a) S and Cl of aggregate fractional melts versus F (melt fraction), (b) S/Cl versus melting fraction (f), and (c) S/Cl versus La/Yb.Melting model ofLee et al. (2012) was implemented in PySulfSat(Wieser & Gleeson, 2023) to emulate fractional melting of a depleted (D-DMM) and mildly enriched (E-DMM) source lithology previously suggested for mid-ocean ridge basalts(Shimizu et al., 2016).The teal and amber model curves in panels (b) and (c) track the melt concentrations of accumulated fractional melts.The dotted pink curve in (c) shows binary mixing of a 2% E-DMM melt and a 20% D-DMM melt.The good match between observed and modeled S/Cl versus La/Yb trends suggests that melt inclusion (MI) S/Cl ratios reflect various degrees of melting of the D-DMM and E-DMM mantle components and subsequent melt mixing.Model parameters are reported in TableS8in Supporting Information S2.Only MIs with MgO >5 wt.% are plotted in panel (c) to minimize the effect of sulfide melt immiscibility during crustal magma evolution.

Figure 8 .
Figure 8. Sulfur emission potentials.(a) ΔS max versus MgO.Pale blue dots show S 2-CSS for raw melt inclusion data calculated after Smythe et al. (2017) in PySulfSat(Wieser & Gleeson, 2023).(b) Total sulfur emission potential ΔS max compared to vent sulfur emissions (ΔS vent ).There is no obvious correlation between eruption type (explosive, effusive, phreatomagmatic) and ΔS max or ∆S vent .Plotted major element compositions of eruptions in panel (a) and Figure9aare based on average matrix glass compositions (TableS2).Error bars are based on the reproducibility of VG-2 and VG-A99 reference glasses (see Figure3caption).

Figure 9 .
Figure 9. Summary figure of the sulfur emission potentials (ΔS max ) of Icelandic volcanoes.(a) TAS diagram.Circles represent individual eruptions and are shaded from blue to red, and increase in size, with increasing ΔS max .(b) ΔS max sorted by volcano.Volcanoes are grouped by volcanic zone and increasing ΔS max .

Figure 10 .
Figure10.Sulfur emissions of Icelandic eruptions based on the petrological method.Estimated S emissions from notable arc and Hawaiian eruptions with eruption mass of >1 Tg are shown for comparison(Oppenheimer et al., 2011; Scaillet et al., 2003 and references therein;Kern et al., 2020).The comparison demonstrates that Icelandic eruptions typically have high sulfur emission potentials (ΔS max , dashed gray lines).
Table S1) was compiled from published geochemical MI data from Icelandic eruptions.The catalog includes, to the authors' best knowledge, all MI data sets published at the time of writing (28 March 2024; n = 38; references are provided in Table 1) that include measurements of one or more of the five major volatile elements in silicate melts: H (here given as H 2 O), C (given as CO 2 ), F, S, and Cl, or B.

Table 1
Overview of Data in the Iceland Melt Inclusion Catalog

Geochemistry, Geophysics, Geosystems 10.1029/2024GC011443 3.5. The Petrological Estimate of Eruptive Volatile Emissions 3.5.1. Calculation of Volatile Fluxes Based on the Petrological Method
At the low-MgO side of this inflexion point, S contents steadily decrease toward intermediate and silicic MIs, which have 0-240 ppm S. Highest MI S contents occur in the SIVZ eruptions of Hekla (1910 CE), Katla (Eldgjá 939 CE) and Surtsey (1963-1967 CE) (Figure For the modeling parameters, we used reported MI major element compositions, Fe 3+ /Fe tot of 0.10 and 0.15 (representing the range measured in Icelandic subglacial and submarine basaltic glasses, Novella Newman and Lowenstern (2002)basaltic MIs (FigureS2ain Supporting Information S1).This implies that S concentrations of basaltic MIs in Iceland have not been affected by degassing.Degassing paths were calculated with the open source COHS-degassing model Sulfur_X(Ding et al., 2023)using the COH model ofNewman and Lowenstern (2002), S speciation model of O'Neill and Mavrogenes (2022), oxygen fugacity buffered at ∆FMQ = 0 and crystallization turned off.Pre-eruptive volatile concentrations were chosen as 0.61 wt.% H 2 O, 1,700 ppm CO 2 and 1,700 ppm S for Holuhraun(Bali et al., 2018)and 2.4 wt.% H 2 O (O'Neill & Mavrogenes, 2002;Smythe et al., 2017;O'Neill & Mavrogenes, 2002;Smythe et al., 2017; Figure Yb ratios in Iceland indicate decreasing melt degree, deeper melting and higher proportion of melts derived from an enriched mantle component (e.g., regional controls appear to influence H 2 O and Cl, and to a lesser extent F concentrations of Icelandic basalts.Highest H 2 O (up to 1.6 wt.% at Ljósufjöll and 1.5% in Surtsey), Cl (up to 2,410 ppm at Fimmvörðuháls, 1,880 ppm in Snaefellsjökull, 1,260 ppm at Ljósufjöll) and F (2,850 ppm, Ljósufjöll) concentrations are all seen in MIs of off-rift basalts in the SIVZ and SNVZ, whereas rift basalts (NRZ, ERZ, WRZ, RP) typically have lower concentrations of H 2 O (<0.8 wt.%), Cl (<500 ppm) and F (<800 ppm).The Cl concentrations in basaltic MIs correlate with La/Yb (Figure6a, R 2 = 0.76).Higher La/

Table 2
Sulfur Emission Potentials Data references for each eruption are shown in the electronic version of Table2.bPublishedvolumeestimates(TableS2).Values in italics indicate DRE values recalculated from tephra volume (TableS2).cMeanMI S concentrations, instead of maximum, are used for calculating ΔS vent and ΔS max for silicic and intermediate eruptions.The uncertainty of electron probe microanalysis (EPMA) S measurements are not consistently reported in the source literature.Relative instrumental errors for EPMA S measurements in our study are on average ±3.7% for MIs (higher S) and 17% for matrix glasses (lower S; TableS3in Supporting Information S2).d SiO 2 , MgO, Na 2 O + K 2 O are average, unnormalized values of MG data (Table a ER thanks the national volcanic hazard assessment program for Iceland (Gosvá), funded by the National Avalanche and Landslide Fund for generous funding that made this work possible, as well as support from the Research Council of Finlandfunded FINMAG project (Grant 338714).SAH acknowledges support from the Icelandic Research Fund (Grant 196139-051).Daniel Ben-Yehoshua and Rob Askew are thanked for bringing home and sharing samples from their Svínafellsjökull expeditions.Sigrún Karlsdóttir is thanked for coordinating the project.The Liquidus petrology discussion group at Cambridge are thanked for helpful comments on a preprint version of the manuscript.Marie Edmonds is acknowledged for smooth editorial handling.The authors are thankful to Ery Hughes and Paul Wallace for thorough and constructive reviews, which significantly improved the quality of the manuscript.