Carbon Fluxes and Primary Magma CO₂ Contents Along the Global Mid-Ocean Ridge System

1 Introduction

The 56,000-km length of the global mid-ocean ridge system represents the dominant locus of magma production and heat loss on Earth. Decades of deep-water oceanographic research have revealed that new seafloor is produced over a limited range of crustal thicknesses (6 ± 1 km, Behn & Grove, 2015; Klein & Langmuir, 1987; Van Avendonk et al., 2017; White et al., 1992) and at variable spreading rates that are well-described globally (e.g., Gale et al., 2013). As a result, magma production rates at mid-ocean ridges are well understood and the fluxes of nonvolatile elements emerging from the mantle at ridges are well constrained. These same research cruises have retrieved tens of thousands of individual submarine mid-ocean ridge basalt (MORB) samples, and the study of these samples has revealed that chemical variations in MORB magmas are produced by polybaric melting of mantle sources at variable temperatures, whereas local deviations from global trends are primarily due to variations in mantle composition (Gale et al., 2014; Kelley et al., 2013; Kinzler & Grove, 1992; Klein & Langmuir, 1987; Langmuir et al., 1992).

Carbon is present in MORBs in trace amounts as CO₂-rich vesicles and as carbonate ions dissolved in quenched seafloor glass (Cartigny et al., 2008; Dixon & Stolper, 1995; Javoy & Pineau, 1991; Marty & Tolstikhin, 1998; Michael, 1995). The solubility of CO₂ in basaltic magma is low at hydrostatic pressures corresponding to seafloor water depths (e.g., Dixon et al., 1995; Jendrzejewski et al., 1997; Shishkina et al., 2010), and thus nearly all magmas erupted at mid-ocean ridges have degassed the majority of their carbon. This flux of CO₂ into the oceans contributes to the global surface carbon cycle, and changes in basaltic magma production and associated degassing have been proposed as one of several important forcing mechanisms that have influenced past global climate variations (Burley & Katz, 2015; Robock, 2000; Shindell et al., 2003; Tolstoy, 2015; Zhong et al., 2010). At the same time, it is now widely recognized that most of the Earth's carbon is present in the Earth's deep interior, rather than at Earth's surface (Marty et al., 2013; Marty & Tolstikhin, 1998). Importantly, basaltic magmatism and subduction of oceanic crust represent the two largest outgassing and ingassing fluxes at the interface between the surface and deep-Earth carbon cycles.

These observations provide the motivation of the present study, in which we apply the rapid, precise, and accurate measurement of volatiles by secondary ion mass spectrometry (SIMS; Hauri et al., 2002) to characterize the extensive worldwide sample collection of MORBs to examine in detail the flux of CO₂ from the global ridge system. In this study we report new SIMS data on 753 geographically distributed MORB samples, and combine these data with a quality-controlled database of published data from an additional 2,446 samples (Figure 1). This combined data set permits a description of CO₂ fluxes, and estimated primary magma CO₂ contents, for 387 of the 458 global ridge 89 segments with major element data, in addition to examination of the variability of primary MORB CO₂ contents, the role of carbon in the melting processes that produce MORB, and the role of the oceanic crust in the context of the Earth's deep carbon cycle.

2 Methods

Here we briefly describe the methods used in this study for sample preparation, SIMS analysis of volatile elements, major elements, and trace elements.

3 Published Data Quality Control

Data from PetDB (www.earthchem.org/petdb) for MORB glasses and melt inclusions, consisting of samples containing data for at least one volatile (CO₂, H₂O, F, S, Cl) were downloaded in late 2017; the data are compiled, along with the newly reported data from this study, in supporting information Table S1, and the references for the original data are included in the supporting information Table S5. All corresponding data (major elements, trace elements, isotopes) were included in this download. Concentrations were recalculated as weight percent for H₂O, and parts per million by weight for CO₂, F, S, and Cl for data reported in other units or as oxides (e.g., C data were converted to CO₂ ppm; SO₂ and SO₃ data were converted to S ppm), and the database was checked for obvious typographic errors by comparison with the original references. Samples with CO₂ and/or H₂O were retained in the database if the samples were determined by SIMS, Fourier-transform infrared spectroscopy, or post-1995 vacuum extraction of the volatile fraction dissolved in the glass; vesicle-only analyses were discarded due to the inability to use the data to calculate a glass concentration. Loss on ignition data for H₂O were also discarded, as were four analyses of melt inclusions from Helo et al. (2011) with >3,900 ppm CO₂ that were judged to have been contaminated by carbon coat during SIMS analysis. Samples with F data were retained only if they were determined by SIMS or ion-selective electrode methods (Bach et al., 1994; Schilling et al., 1983); electron probe microanalysis (EPMA) for F is subject to interference from overlapping Fe peaks in basaltic glasses, and these data were thus rejected from the database. Samples with S or Cl data were retained if they were determined by SIMS or EPMA.

We also compiled previously published Fe³⁺/ΣFe ratios of MORB glasses for this study. For this compilation, we report Fe-speciation data collected using micro X-ray absorption near-edge structure (μXANES) spectroscopy, with all data collected at beamline X26A of the National Synchrotron Light Source, Brookhaven National Laboratory, following methods described by Cottrell et al. (2009) and Zhang et al. (2018). The Zhang et al. (2018) study assessed the effects of recoil-free fraction on the room temperature Mössbauer spectra, against which the XANES centroid energy was calibrated by Cottrell et al. (2009). Room temperature Mössbauer spectra overestimate the Fe³⁺/∑Fe ratios of basaltic glasses, and the updated calibration of Zhang et al. (2018) accounts for this effect. The Fe³⁺/ΣFe ratios compiled here include those reported by Cottrell and Kelley (2011) and Le Voyer et al. (2015), both updated using the new calibration of Zhang et al. (2018).

Interlaboratory bias between the DTM SIMS method and published data sets is essentially nonexistent, for several reasons. SIMS calibrations for H₂O and CO₂ rely on the use of an extensive suite of MORB and back-arc basin glasses previously analyzed by vacuum manometry and/or Fourier-transform infrared spectroscopy (references in Hauri et al., 2002) originating from laboratories that provided the vast majority of previously published data archived in PetDB. The SIMS method for S and Cl is more precise and has lower detection limits than EPMA. The precision of bias factors between data sets is limited by the reproducibility of the EPMA data and is thus indeterminate. For F the SIMS calibration, which relies on only four glasses (Hauri et al., 2002), was checked against the data sets of Rowe and Schilling (1979), Schilling et al. (1983), and Bach et al. (1994) by SIMS analysis of samples in common with these studies, which used a halogen-selective electrode wet chemical method. We found excellent agreement (±8% 2σ) with the data set of Bach et al. (1994) but much poorer agreement with the other, significantly older, data sets that used the same method for F analysis. Thus, the comparison for F found no interlaboratory bias but does show clear examples of variable extents of reproducibility of other methods compared with SIMS, which likely permeate much of the published data for fluorine and, to a lesser extent, sulfur and chlorine. In general, interlaboratory comparisons agree to better than 10% for CO₂, H₂O, and F, and by a somewhat worse and variable degree for S and Cl.

The supporting information contains the data used in this paper and is available through the EarthChem Library (https://doi.org/10.1594/IEDA/111195). Table S1 contains the results of new analyses of MORB glasses, together with associated data on major and trace elements, and isotopic compositions. In Table S2, these new data are compiled along with the quality-controlled published data in a single table. Table S3 is a compilation of data (volatiles, major elements, trace elements, isotope compositions of host glass) for MORB melt inclusions, for comparison to the submarine glass data. Table S4 is a modification of the segment-averaged tables from Gale et al. (2013, 2014) that includes data for segment metrics (length, spreading rate, depth, etc.) as well as averaged data for major elements, trace elements, volatiles, and isotopic compositions; this table also includes estimated primary magma CO₂ concentrations, and calculated CO₂ fluxes on a per-segment basis, with estimated uncertainties shown in yellow highlight. Table S5 contains an alphabetic list of all the publications that contained the original data that were used in the data compilations of all of the data tables.

4 Results

The samples in this compiled volatile data set were assigned to individual ridge segments based on a comparison of their sampling coordinates with the global spreading center segment catalog of Gale et al. (2013; Figure 1). In constructing geochemical averages of global spreading segments, Gale et al. (2013) included in their compilation samples within 10 km of the ridge axis. At the slowest spreading rates, a distance of 10 km off-axis corresponds to an age of approximately 1.5 million years, and so in our study individual samples were categorized as on-axis if their apparent age was less than or equal to 1.5 million years. Samples older or farther than this were categorized as either off-axis or as transform samples depending on their location relative to the nearest ridge axis or transform fault. Segment average volatile concentrations were then calculated for all segments that contained at least one on-axis sample analyzed for H₂O, CO₂, F, S, and/or Cl. The PetDB data download also contained a small amount of new major element, trace element, and isotope data published after the compilations of Gale et al. (2013) and Class and Lehnert (2012), and this permitted the calculation of several new segment averages for major elements, trace elements, and isotopes to augment the database. For consistency, these updated segment averages were calculated using the same methods described in Gale et al. (2013); however, it should be mentioned that only samples with volatile data were used to calculate these additional new segment averages. Finally, we recalculated the segment-averaged major and trace element data for segment EPRR31 to account for two additional samples from this segment analyzed by Cottrell and Kelley (2011); samples from this segment have a bimodal distribution of depleted (n = 3) and enriched (n = 3) MORB. By doing this, the segment averages for all data (major elements, trace elements, volatiles) all come from the same six samples.

The MORB glasses measured in this study range from 2 to 420 ppm CO₂, 0.04 to 1.37 wt% H₂O, 41 to 1,131 ppm F, 121 to 2,380 ppm S, and 1.7 to 5,540 ppm Cl. These new volatile data extend the range of MORB H₂O, F, and Cl abundances, while the data for CO₂ and S cover most of the range defined by previously published data for MORB (up to 600 ppm CO₂ and 3,300 ppm S). H₂O and F are broadly correlated in the global data (Figure 2a), with an average H₂O/F ratio of 14.3 but with a large range displayed by individual samples (0.56 to 216). Sulfur is correlated with total Fe (FeO*) (Figure 2b), reflecting the dependence of sulfide solubility on iron content in mafic melts (Keith et al., 2017; O'Neill & Mavrogenes, 2002). Chlorine is generally not well correlated with H₂O, but a lower bound to the global data is apparent with a Cl/H₂O ratio in the range of 0.01 to 0.05. Most of the samples with higher Cl/H₂O ratios have likely undergone assimilation of Cl-rich seawater-derived components (Le Roux et al., 2006; Le Voyer et al., 2015; Michael & Schilling, 1989). There does not appear to be any relationship of the volatiles with the oxidation state of Fe in MORBs examined in this study.

Segment averages for the volatile data were corrected to MgO = 8 wt% and to compositions in equilibrium with Fo₉₀ olivine, in order to compare with the fractionation-normalized database of Gale et al. (2014). This correction for crystal fractionation cannot account for the irreversible loss of CO₂ by degassing; however, we will estimate primary magma CO₂ contents using Ba and Rb as proxies, in a subsequent section (below). The scarcity of available volatile data at any given ridge segment, coupled with potential CO₂, H₂O, or S loss by magmatic outgassing, resulted in poorly defined liquid lines of descent (LLDs) for the volatile data. Thus the 8% MgO corrections for the volatile segment averages were done by calculating a fractionation factor from the ratio of the segment average P₂O₅(8) divided by the segment average P₂O₅ (uncorrected) using the data given in Gale et al. (2013, 2014). P₂O₅ was chosen because it typically forms better-correlated and more linear LLDs than K₂O; nevertheless, for segments lacking P₂O₅ data the K₂O(8) and K₂O data were used. An identical procedure was used to further correct the volatile segment averages to Fo₉₀ values, by determining a fractionation factor from the ratio of K₂O(90) to K₂O(8) for each segment reported by Gale et al. (2014).

Uncertainties from numerous sources were quantified throughout the examination of the data. Analytical uncertainties for CO₂ and other elements are typically minor (a few percent relative) compared to the compositional heterogeneity within individual ridge segments; the standard deviations of segment averages were quantified in Gale et al. (2013) while standard deviations for compositions corrected to 8 wt% MgO were quantified in Gale et al. (2014). For the purposes of calculating segment magma fluxes, we have assigned a ±4 km uncertainty on ridge length, ±3 mm/year uncertainty in spreading rates, and a variable uncertainty on crustal thickness (ranging from 3% to 50%) determined by the uncertainty of a regression of estimated errors in crustal thickness versus absolute thickness determined seismically (Behn & Grove, 2015; Van Avendonk et al., 2017).

The uncertainty that is most difficult to quantify is our assumption of single average values for the ratios CO₂/Ba (81.3 ± 28%) and CO₂/Rb (991 ± 13%) for MORB; these values are the averages of the population averages of six groups of CO₂-undersaturated MORB glasses and melt inclusions (Hauri et al., 2017). To the extent that normal MORB (N-MORB) can be characterized as a mixture of depleted MORB (D-MORB) and enriched MORB (E-MORB; e.g., Donnelly et al., 2004), this assumption implies that D-MORB and E-MORB sources have the same CO₂/Ba and CO₂/Rb ratios. Two studies support this assumption. Shimizu et al. (2016) estimated CO₂ contents in Pacific E-MORB from their Cl contents and the constant CO₂/Cl ratio in Siqueiros Fracture Zone D-MORB melt inclusions, yielding CO₂/Ba ratios of 105 ± 61 for Pacific D-MORB and 134 ± 23 for Pacific EMORB. Similarly, in the Atlantic the D-MORB melt inclusion suite of Le Voyer et al. (2017) has CO₂/Ba of 97 ± 10, very close to the CO₂/Ba of 106 ± 34 of the highly enriched E-MORB popping rocks from 14°N on the Mid-Atlantic ridge (MAR; Cartigny et al., 2008; Hauri et al., 2017). Despite the possibility that this assumption may eventually be contradicted, it presently provides a data-supported basis from which to estimate the primitive, pre-eruptive CO₂ content for degassed MORB. We emphasize here that a key underpinning of our estimations of undegassed global MORB CO₂ contents is the assumed constancy of CO₂/Ba and/or CO₂/Rb across the full compositional spectrum of MORBs; uncertainties on the CO₂ concentrations estimated here will increase if these ratios are later shown to vary. We calculate undegassed MORB CO₂ concentrations by multiplying these ratios by the Ba and/or Rb concentrations measured in the basalts. In our calculations of primary CO₂ content, segment-scale CO₂ fluxes, and overall ridge flux we have combined all of these sources of error in quadrature, resulting in the stated uncertainties listed with the data in the yellow-highlighted columns of Table S4.

5 Degassing of CO₂ From MORB

MORBs are erupted on the seafloor under the hydrostatic pressure of the overlying water column. The solubility of CO₂ in MORB is strongly dependent on pressure (Dixon & Stolper, 1995), which would predict that MORB CO₂ concentrations should correlate with the pressure of eruption on the seafloor. Yet MORBs often have dissolved CO₂ contents that exceed the equilibrium solubility predicted at the eruption pressure. In most cases, this oversaturation is due to rapid melt ascent during eruption, during which diffusion of carbon in the melt is too slow to maintain equilibrium between the dissolved CO₃²⁻ in the melt and the CO₂ vapor within the vesicles (Chavrit et al., 2012; Javoy & Pineau, 1991).

The vapor solubility model of Dixon et al. (1995) has long been used to constrain vapor saturation pressure (P_sat) in MORBs because these authors determined CO₂ and H₂O solubility specifically on MORB melt compositions. When we use measured CO₂ and H₂O contents in MORBs to calculate P_sat using the Dixon et al. (1995) model, and compare these with the pressures of MORB sample collection (P_sample), the distribution of the P_sat/P_sample ratio is strongly peaked near unity (median = 1.37) but is clearly skewed toward high P_sat/P_sample ratios (Figure 3). This skew is likely due to vapor oversaturation in MORBs driven by high magma ascent rates, as discussed above. At the same time, there is a sharp cut-off below unity, with few samples displaying P_sat/P_sample significantly <1.

A similar comparison of P_sample with P_sat predicted by other experimental solubility studies provides a rough approximation of the accuracy of these solubility models. The model of Papale et al. (2006) produces P_sat/P_sample ratios that are 60% lower on average compared with Dixon et al. (1995), and suggests that most MORB are vapor undersaturated. Vesicle abundances observed in MORB, however, contradict this prediction (Chavrit et al., 2012), indicating that the Papale et al. (2006) model likely underestimates vapor saturation pressures in MORB. The Papale et al. (2006) model was calibrated over a very wide range in melt compositions, whereas the Dixon et al. (1995) model was calibrated specifically on MORB compositions at pressures <1 kbar, so disagreement between these models is perhaps not surprising. The solubility model of Iacono-Marziano et al. (2012) provides values of P_sat for MORB that are within 2% of those calculated from the Dixon et al. (1995) model. More recently, the CO₂-H₂O solubility model of Ghiorso and Gualda (2015) results in P_sat/P_sample ratios that are 20% lower on average compared with Dixon et al. (1995), with a peak of the P_sat/P_sample ratio distribution closer to unity (median 1.10) and a lower limit of the ratio of ~0.6. The Ghiorso and Gualda (2015) model thus also predicts more samples than Dixon et al. (1995) to be vapor undersaturated. Both models provide reasonable estimates of vapor saturation pressures, but we prefer the model of Dixon et al. (1995) because it is specific to MORB and it predicts that fewer MORBs are vapor undersaturated, which is consistent with the presence of vesicles in the vast majority of MORBs.

The ratios of volatile to nonvolatile elements have been useful to examine processes of loss or gain of volatile elements, such as degassing, dehydration, and shallow-level assimilation that do not affect nonvolatile elements that otherwise share similar mineral-melt partitioning behavior during partial melting and fractional crystallization. Further constraints on the extent of degassing of individual MORB samples can thus be obtained by examination of the ratios of CO₂/Rb and CO₂/Ba, which exhibit a narrow range in genuine vapor-undersaturated MORB glasses and melt inclusions. Rosenthal et al. (2015), Michael and Graham (2015), and Hauri et al. (2017) recently reviewed the data for vapor-undersaturated MORB; in this compilation we use average ratios CO₂/Rb (991 ± 26%) and CO₂/Ba (81.3 ± 56%) to represent undegassed MORBs, which are calculated as the average of the population averages for six groups of vapor-undersaturated samples given in Table 1 of Hauri et al. (2017). CO₂/Rb and CO₂/Ba ratios are much more uniform than CO₂/Nb or CO₂/Th ratios in undegassed MORBs. The ratio of the measured CO₂/Ba ratio in degassed MORBs to the average CO₂/Ba ratio for undegassed MORB gives a measure of the amount of CO₂ lost by degassing. The extent of degassing of typical MORB samples estimated in this way demonstrates that a majority of samples degassed >70% of their initial inventory of CO₂ (Figure 4). Approximately 3% of MORB samples, however, exhibit CO₂/Rb ratios close to the undegassed ratio, and most of these samples have been identified previously (Hauri et al., 2017; Le Voyer et al., 2017; Michael & Graham, 2015; Saal et al., 2002; Shimizu et al., 2016). Based on this analysis, emplacement of MORB on the seafloor typically results in eruptive degassing into the oceans of ≥70% of the primary carbon carried by MORB magmas.

The CO₂/Ba ratios of melt inclusions from mid-ocean ridges also bear importantly on the carbon budget of MORB magmas. The olivine-hosted melt inclusion data of Saal et al. (2002; Siqueiros FZ, East Pacific Rise), Le Voyer et al. (2017; equatorial mid-Atlantic Ridge), and Hauri et al. (2017; northern Iceland), for example, all display excellent correlations of CO₂ with Ba (and with Rb, where measured), demonstrating the undegassed nature of these melt inclusion populations. An additional 113 melt inclusions have been measured for CO₂ and Ba from the Gakkel Ridge (Shaw et al., 2010; Wanless et al., 2014), the Juan de Fuca Ridge (Helo et al., 2011; Wanless & Shaw, 2012), and from the Lucky Strike segment of the mid-Atlantic Ridge (Wanless et al., 2015). Of these 113 melt inclusions, 22 exhibit CO₂/Ba ratios within 20% of the average undegassed ratio of 81.3 (Figure 4b). The majority have likely experienced CO₂ degassing either prior to melt inclusion entrapment or into secondary vapor bubbles that often form within melt inclusions during post-entrapment cooling (e.g., Moore et al., 2015) but are difficult to measure quantitatively. More data could be readily obtained from MORB melt inclusions through concentrated efforts to quantify the effect of vapor bubbles within melt inclusions (e.g., Aster et al., 2016) and by targeting ridge segments that erupt samples on the seafloor that have experienced little to no degassing (e.g., Le Voyer et al., 2017).

Matthews et al. (2017) present an alternative interpretation of certain suites of MORBs with CO₂ that is well correlated with highly incompatible trace elements. Such suites are often interpreted as undegassed with respect to CO₂, because the degassing process fractionates CO₂ from other highly incompatible elements. Matthews et al. (2017) instead envisages a process by which magmas emerging from the melting regime stall at a depth where enriched magma batches become oversaturated in CO₂, leading to degassing. This degassing is followed by limited mixing at shallower levels, which can produce correlated CO₂/Ba ratios in MORBs that are lower than the original mantle source CO₂/Ba ratio. Yet this process is limited in its ability to explain CO₂/Ba ratios in highly depleted MORB, many of which arrive on the seafloor with low CO₂ contents but also vapor undersaturated. Such magmas would need to degas CO₂ at pressures significantly lower than their seafloor collection pressures (e.g., Siquerios data of Saal et al. (2002); Quebrada-Discovery-Gofar Transform data of Shimizu et al. (2016); the vapor-undersaturated data of Michael and Graham (2015)), or at pressures within the eruptive Layer 2 of the oceanic crust (equatorial MAR data of Le Voyer et al. (2017)). For the more intermediate-enriched data sets (Borgarhraun, Iceland, and Top 10% MORB data sets of Hauri et al., 2017; 14°N MAR, Cartigny et al., 2008) the highest CO₂ contents vary widely, and each sample suite would need to have stalled and degassed at a very specific depth (e.g., not just the Moho depth) in order to produce such similar CO₂/Ba ratios as reported in Hauri et al. (2017). We thus conclude that, although the processes in Matthews et al. (2017) are plausible, they are unlikely to have operated in the samples that currently constrain the CO₂/Ba ratio of MORB sources.

6 Variations in Primary MORB CO₂ Contents

The restricted ranges in CO₂/Rb and CO₂/Ba ratios of undegassed MORB permit an estimation of the initial (undegassed) CO₂ concentration in MORB primary magmas, given an estimate of their primary Rb and/or Ba abundances. The studies of Gale et al. (2013, 2014) provide such estimates based on the compositions of MORB magmas averaged by ridge segment, and corrected for shallow-level crystal fractionation to MgO = 8 wt% and also further corrected to be in equilibrium with Fo₉₀ olivine. The data compilations of Gale et al. (2014) and this study provide primary MORB magma estimates for 387 segments that contain data for Ba and/or Rb, out of the 711 global spreading ridge segments identified by Gale et al. (2013). These segment-averaged concentrations were then multiplied by the undegassed-average ratios of CO₂/Rb (991) and CO₂/Ba (81.3) to calculate segment-averaged primary CO₂ concentrations. MORB primary CO₂ thus calculated varies over a factor of 180, from 104 ppm to 1.90 wt%. These calculated primary CO₂ contents further enable calculation of CO₂/Nb and CO₂/Th ratios for all segments that have segment-averaged Nb and Th data. These latter ratios are used to compare with prior studies that used CO₂/Nb ratios (e.g., Cartigny et al., 2008; Saal et al., 2002), while Th is another highly incompatible element that has a mantle melt partition coefficient similar to CO₂ (Rosenthal et al., 2015), and there exists Th data of sufficient quality in MORB and MORB melt inclusions to make a useful comparison (Hauri et al., 2017).

Primary MORB CO₂ abundances can be compared with other indices of depletion and enrichment in MORB such as the K₂O/TiO₂ ratio, which is low in depleted MORB and high in enriched MORB. The lowest calculated primary MORB CO₂ abundance is 104 ppm CO₂ at a depleted segment of the Galapagos Spreading Center (segment GALA16 with K₂O/TiO₂ of 0.05), and the highest is 1.90 wt% CO₂ at a highly enriched segment of the Southwest Indian Ridge (segment SWIR9 with K₂O/TiO₂ = 0.66). The distribution of MORB primary CO₂ contents is strongly skewed with a median of 1107 ppm, a mean of 2139 ppm, and a mode, determined using a kernel density analysis with an optimal bandwidth of 318 (Bowman & Azzalini, 1997), of 621 ppm CO₂ (see Figure S1 in the Supporting Information). In general, the segment-averaged primary CO₂ contents, as derived from reconstructed primary Ba and Rb concentrations, correlate well with segment-averaged primary Th and K₂O contents (Figure 5). Such correlations are most reasonably interpreted as signatures of variation in mantle source composition, rather than variable mantle temperature or spreading rate. In this view, variations in the primary MORB CO₂ contents determined here would also reflect mantle source variations.

Segment-averaged CO₂/Nb ratios calculated above vary by more than a factor of 10, from a low of 152 at the very depleted segment GALA16 on the Galapagos Spreading Center, to a high of 1,760 at enriched segment SWIR22 on the Southwest Indian Ridge. This range in CO₂/Nb contrasts with the much more limited range in CO₂/Nb measured in undegassed MORB glasses or melt inclusion populations (230–730, Cartigny et al., 2008; Hauri et al., 2017; Saal et al., 2002). We emphasize again that the CO₂/Nb ratios constrained here effectively reflect variations in Ba/Nb or Rb/Nb in MORB because of the assumptions adopted in the calculation of CO₂ concentrations. Yet the calculated segment-averaged CO₂/Nb ratios display rough trends with segment-averaged Sr and Nd isotopes (Figure 6) that are consistent with the relationships observed by Hauri et al. (2017) and further support the inference that primary CO₂ contents of MORB, inasmuch as these are inseparable from Ba and Rb, are a function of mantle composition. It remains to be seen if this large range in estimated CO₂/Nb ratios is reflected in actual sample measurements, but the data in this study provide a first-order prediction and help identify specific ridge segments to focus on to test the accuracy of these estimates. The variations are consistent with the dichotomy of CO₂/Nb and CO₂/Th ratios between Atlantic and Pacific Ocean basins, which could be a result of a greater amount of recycled subduction zone components in the Atlantic Ocean basin compared with the relatively older Pacific Ocean basin (Dixon et al., 2017; Hauri et al., 2017).

7 The Flux of CO₂ From the Global Ridge System

Crustal thickness at mid-ocean ridges is known to vary from a few kilometers to over 20 km at northern Iceland (Menke, 1999), controlled by a combination of the local magma production rate and the spreading rate. The study of Behn and Grove (2015) summarized seismically determined measurements of crustal thickness at ridge axes, where it varies from 2 km at the slow-spreading Gakkel Ridge to 30 km at central Iceland; additional data for crust <20 Ma old was taken from Van Avendonk et al. (2017). Seismically determined crustal thickness correlates well with mean ridge depth (Figure 7), reflecting isostatic compensation of MORB crust at the ridge. These constraints permit an estimation of the crustal thickness at all of the 711 segments defined by Gale et al. (2013) based on the axial depth. We then calculate a segment-averaged magma production rate by multiplying the crustal thickness by the mean segment spreading rate, which we then use to calculate the fluxes of CO₂ (and other elements) on a segment-by-segment basis. There is one exception to the relationship between ridge depth and crustal thickness at the northernmost three segments of the Alarcon Ridge, where the ridge crest shallows as it begins to intersect the incipient continental rift in the Gulf of California; here we have chosen to assign an average crustal thickness of 6.4 km to these three ridge segments (EPRR1, EPRR2, and EPRR3).

Carbon dioxide fluxes, normalized by ridge length, vary by a factor of >100 because they are the products of local primary magma CO₂ contents and local magma production rates, which are uncorrelated. Fluxes of CO₂ vary from a low of 1.52 × 10⁶ mol/year/km at segment SWIR17 (which is both depleted and slow-spreading) to a high of 4.74 × 10⁸ mol/year/km at segment JUAN1 (which has high primary CO₂ and moderate spreading rate). High CO₂ fluxes are independent of spreading rate, occurring at both slow- and fast-spreading ridge segments. High CO₂ fluxes also occur on segments both near to and far from hotspots, and there is a clear tendency for estimated CO₂ fluxes to increase toward hotspots (Figure 8), in concert with enrichment in highly incompatible elements (i.e., Ba, Rb). Provided that the mechanism of mantle enrichment for hotspots does not fractionate CO₂ from Ba or Rb, then this observation supports our perspective that variations in the MORB primary magma CO₂ content are mostly driven by variations in mantle composition rather than magma production dynamics or mantle temperature.

When summed over the length of the global mid-ocean ridge system, we calculate the global MORB magma production rate to be 16.5 km³/year, significantly lower than the widely cited 21 km³/year estimated by Crisp (1984); this difference in magma production rate alone accounts for 27% of our total flux calculation. The total flux of CO₂, summed over the widely varying fluxes at individual ridge segments (Figure 9), is  × 10¹² mol/year. Importantly, for individual ridge segments with no Ba or Rb data from which to constrain CO2, we assume an average magmatic CO2 content calculated from the log-normal mean of those segments with Ba or Rb constraints. This estimate lies within the range of CO₂ fluxes estimated by Resing et al. (2004) from CO₂ and ³He fluxes into the water column above ridges, and by Barry et al. (2014) in their study of Iceland CO₂ fluxes. Our estimate overlaps the lower average flux found by Saal et al. (2002; 0.93 ± 0.28 × 10¹² mol/year), who based their estimate on the CO₂ content of the depleted Siqueiros Fracture Zone melt inclusions. Our estimate is lower than the estimates of Marty and Tolstikhin (1998), Cartigny et al. (2008), Michael and Graham (2015), Le Voyer et al. (2017), and Hauri et al. (2017) based on a variety of approaches, but these prior studies utilize much more limited data sets than we have employed here. Our estimate is substantially higher than the study of Chavrit et al. (2014) who considered vesicle abundances in MORB; their underestimation is likely biased by CO₂ contained in vesicles that were lost from samples by degassing. In general there is a good correlation of segment-scale CO₂ flux with the segment-averaged primary CO₂ content (Figure 9), while the correlation of CO₂ flux with magma flux is much weaker.

The MORB flux of CO₂ is nearly quantitatively delivered to the overlying water column, during eruptive degassing and later by hydrothermal circulation, where it becomes part of the marine bicarbonate cycle. Short-term global effects of large variations in the ridge CO₂ flux are, however, muted by the ~100,000-year residence time of carbon dissolved in seawater (Broecker & Peng, 1986) and thus variations in the ridge CO₂ flux are likely to influence short-term surface climate conditions only if large fluxes bypass the water column in areas where the ridge is shallow or exposed above sea level. On the longer timescale, the ridge CO₂ flux is dwarfed by the deposition of carbonate sediments and organic carbon on the seafloor and formation of carbonate alteration veins within the crustal and mantle sections of the oceanic lithosphere, making the seafloor a net sink for carbon by a large factor (Alt & Shanks, 2003; Alt & Teagle, 1999; Hayes & Waldbauer, 2006; Kelemen & Manning, 2015). Much (but not all) of this carbon is present in forms, and at locations, where it can be returned to the deep crust and/or mantle in subduction zones.

8 Influence of Carbon on Mantle Melting

Several studies have demonstrated that the presence of carbonate in the mantle can strongly influence the depth of melt initiation beneath ridges (Dasgupta & Hirschmann, 2006, 2007; Dasgupta et al., 2013; Keller et al., 2017; Stagno et al., 2013), as incipient melting of carbonate may represent the first step in a polybaric process that ultimately results in eruption of MORB. To first order, the 2-order-of-magnitude range in primary MORB CO₂ content must reflect a large variation in the carbon content of the MORB source mantle. Le Voyer et al. (2017) hypothesized how changes in the carbon content of the mantle could result in geophysical observations illuminating the depth beneath ridges where low, but possibly interconnected, melt fractions (near 0.05%; Minarik & Watson, 1995) could first be detected at constant potential temperature. Specifically, they showed that changing the CO₂ content of relatively oxidized and anhydrous mantle from 20 ppm to 140 ppm increased the depth where interconnected melt (F~ 0.05%) would exist from 75 km to 275 km. If a combination of seafloor seismic and electromagnetic methods could elucidate clearly such low melt fractions (Eilon & Abers, 2017; Gaillard et al., 2008; Sifré et al., 2014), then a test of this hypothesis might be possible by comparing the geophysical signals at ridge segments characterized by widely different primary magma CO₂ contents as identified here.

This test is complicated, however, if carbonate is not the stable phase of carbon in the mantle at the depth of incipient melting. While direct analyses of MORBs and their peridotitic residues suggest that mantle oxygen fugacity (fO₂) is near the quartz-fayalite-magnetite buffer, where carbonate would be the stable phase (Bezos & Humler, 2005; Cottrell & Kelley, 2011; Birner et al., 2018), oxygen fugacity (fO₂) may decrease with depth in the convecting oceanic mantle, as observed in the cratonic xenolith record and supported by theory and experiment (Eggler & Baker, 1982; O'Neill et al., 1993; Ballhaus, 1995; Frost & McCammon, 2008; Rohrbach et al., 2007; Stagno and Frost, 2010; Stagno et al., 2013). It has been suggested that this decrease in mantle fO₂ would reduce carbonate to graphite, diamond, or metal alloy, and this would prevent carbon from fluxing the mantle. In this scenario, the existence of an interconnected melt network would not occur until upwelling mantle became oxidized enough to stabilize carbonate (~ QFM-2) or reached the onset of silicate melting, (e.g., Stagno et al., 2013). However, even in the limit where the mantle is too reduced to stabilize carbonate, carbon would dissolve in silicate melts and become exhausted within the first 1–2% of melting for typical carbon concentrations that characterize the MORB mantle (20–200 ppm C; Le Voyer et al., 2017; Rosenthal et al., 2015; Saal et al., 2002). Thus, whether the mantle is reduced or oxidized, and regardless of whether carbon actively fluxes mantle melting or not, we expect carbon to behave like a highly incompatible element and produce the correlations of C with Ba, Rb, Th, and Nb that characterize suites of undegassed melt inclusions for typical mantle carbon concentrations (cf. Hauri et al., 2017). That we observe carbon's incompatible behavior up to much higher carbon concentrations supports other inferences, based on MORB and peridotite analyses, that graphite is not stable at shallow depths beneath the ridge.

9 Conclusions

We have determined the concentration of CO₂ and other volatiles in 753 MORB glasses by SIMS from all the major ocean basins and, when combined with data from prior studies, these data provide global coverage of carbon fluxes from mid-ocean ridges. Assuming constant ratios of CO₂/Ba and CO₂/Rb, we estimate primary magma CO₂ concentrations that vary from 104 ppm to 1.90 wt%, MORB CO₂ fluxes that range from 3.3 × 10⁷ to 4.0 × 10¹⁰ mol/year at the scale of individual ridge segments, and a global integrated CO₂ flux of  × 10¹² mol/year, overlapping the global CO₂ flux estimated from depleted MORB (0.93 ± 0.28 × 10¹² mol/year, Saal et al., 2002) but considerably lower than the flux estimated from E-MORB popping rock at 14°N on the mid-Atlantic Ridge (15 × 10¹² mol/year, Javoy & Pineau, 1991). The accuracy of the CO₂ fluxes determined here is due to consideration of the heterogeneity of trace element abundances in MORB along the global ridge system, with the assumption that our global CO₂/Ba and CO₂/Rb ratios are both representative and homogeneous in the MORB mantle.

The primary magma abundances, and segment-scale fluxes, of mid-ocean ridge CO₂ are heterogeneously distributed both along the ridge system and between major ocean basins. Areas of both high flux and high primary CO₂ are predicted at hotspots (Iceland, Azores, Bouvet, Afar, Reunion, Easter) and regionally along the Gakkel Ridge (Arctic). On the East Pacific Rise, these fluxes are correlated with deep areas of elevated mantle temperatures (Dalton et al., 2014) and high magma production rates. High primary CO₂ contents are, however, mainly predicted for areas of enriched mantle composition, which may be associated with long-term enrichment of the mantle as measured by radiogenic isotopes.

Large variations in primary MORB CO₂ imply similar large variations in mantle carbon content beneath ridges. If the upper mantle is sufficiently oxidized, these variations would result in large changes in the depth at which interconnected melt (F ~ 0.05%) would exist (from 75 to 275 km). Such large depth differences might be resolvable using geophysical methods by comparing ridge segments characterized by vastly different CO₂ contents as we have described here.