Diapirism of carbonate platforms subducted into the upper mantle

The fate of subducted carbon is of great importance for understanding deep carbon budgets on geological time scales (Dasgupta and Hirschmann, 2010; Pall et al., 2018). Carbon enters subduction zones primarily as sedimentary carbonate minerals (calcite or dolomite) in various forms, including carbonate-rich pelite, abundant calcite fracture fill in greywackes, and various forms of organic carbonate precipitation, among others (Edmond and Huh, 2003; van der Ploeg et al., 2019). One of the potentially most important sinks of carbon, subduction of large carbonate platforms (Handford and Loucks, 1993), remains poorly understood. Some carbon is incorporated in fold-and-thrust belts into the mid- to lower crust of thickened collisional orogens (Selverstone and Gutzler, 1993). Subducted carbonate is thought to be released from the forearc (Stewart and Ague, 2020) and the subarc region (Plank and Manning, 2019) by decarbonation reactions, some involving melt (Duncan and Dasgupta, 2014). The majority of that released carbon, liberated as CO₂, makes its way back to the surface either via arc magmatism (Lee and Lackey, 2015) or through fault pathways as gas (Italiano et al., 2017). A small component of subducted carbon remains in the mantle lithosphere (Kelemen and Manning, 2015) or is transferred to the deep mantle, and, if deep enough, it can turn into diamonds and other high-pressure carbon-bearing phases (Shirey et al., 2013).

Direct evidence for carbonate in the mantle comes from peridotite xenoliths that contain primary carbonate (Ionov et al., 1995; Ducea et al., 2005; Chen et al., 2018); these were found in xenolith suites from numerous localities and ages globally. Their primary origin is argued based on textural evidence, and the subducted, initially sedimentary, origin is demonstrated based on stable isotopes (Ducea et al., 2005). Not all carbon coming from the mantle is recycled (Stachel et al., 2009).

Carbon dioxide budgets in subduction have been quantified in recent studies (e.g., Kelemen and Manning, 2015; Galvez and Pubellier, 2019; Plank and Manning, 2019). These studies acknowledge that Mediterranean subduction systems could significantly change these budgets via subduction of large amounts of carbonate. Various basins of the Tethys realm in Europe (van Hinsbergen et al., 2020) were relatively narrow and shallow basins (some were true oceanic basins; others, highly extended continental; Fig. 1) commonly covered by marine carbonate sequences of Mesozoic or Cenozoic age. The low latitude and mostly east-west orientation of these margins made them prone to carbonate formation, very different from most modern subduction systems that are found around the Pacific Ocean and that incorporate nontrivial, but much smaller, amounts of carbonate in their budget (Plank and Langmuir, 1988). Much of that carbonate mass is seen in Alpine-aged European fold-and-thrust belts as well as undeformed sequences on the continental margins that escaped subduction. One particularly important aspect of carbon pathways into the mantle that has only recently been investigated (Chen et al., 2021) is the possible mechanical transport by diapirism of large carbonate masses in the uppermost mantle above subduction zones. This process is similar to the relamination hypothesis for siliciclastic sediments (Hacker et al., 2011; Miller and Behn, 2012), which is based on earlier geodynamic forward models (Currie et al., 2007).

In our study, we used geodynamic models to investigate the pathways and effects of subduction of large carbonate sequences (hundreds of kilometers wide, and a minimum of 1 km thick). These models simulate carbonate sequences that could be involved in interaction with subduction environments: (1) shallow-water carbonate sequences (including here carbonate platforms, ramps, mud-mounds, and rimmed and non-rimmed shelves) and pelagic carbonate platforms developed along continental margins; and (2) pelagic carbonate seafloor sediments developed on oceanic crusts. Figure 1 shows two snapshot paleogeographic reconstructions of the extent of carbonate platforms in the Late Triassic and the Late Jurassic, respectively, in the Tethyan domain of today's Europe (from Cosentino et al., 2010). Only a small fraction of this mass of carbonate remains preserved in the fold-and-thrust belts of Alpine Europe or on foreland basin margins and other platformal regions that escaped subduction (e.g., Moesia [in the Balkans]). Geologic constraints for the size and fate of these carbonate volumes are primarily drawn from the European Alps and Carpathian Mountains, and are described in the Supplemental Material¹.

We used the SOPALE numerical code (http://geodynamics.oceanography.dal.ca/sopaledoc.html) to model the coupled thermal-mechanical evolution of the lithosphere–upper mantle system (Fig. 2; Figs. S1–S3 in the Supplemental Material). The two-dimensional model domain is 2000 km wide and 600 km deep. The initial model geometry is shown in Figure S1 and consists of an oceanic plate between two continental plates. Similar results would be obtained if a highly extended and thinned continental plate, such as those found in passive-margin continental shelves, were used instead of a true oceanic plate. The oceanic plate (7-km-thick oceanic crust and 64-km-thick mantle lithosphere) is overlain by a 1 km layer of carbonates, making it a 72-km-thick plate corresponding to a ca. 30 Ma oceanic plate (see the Supplemental Material for justification of the choice of thickness of the carbonate layer), and subducts below a 50-km-thick continental plate at a rate of 2.5 cm/yr. The 1 km thickness is constrained by geologic observations but is also the lower limit for computation here; thinner layers cannot be modeled yet due to computational limitations. Additional effects such as extensive serpentinization of the mantle in the downgoing plate can enhance buoyancy of the plate Other effects such as a thinner carbonate layer would decrease the effect; they are not considered here. Our focus is the carbonate layer atop the oceanic or highly attenuated continental plate, assumed to be predominantly carbonate mixed with 20% pelagic and terrigenous sediments. The density of this layer is 2700 kg/m³ at 200 °C, making it 100 kg/m³ less dense than upper and mid continental crust and ~650 kg/m³ less dense than mantle when materials are at the same temperature. Details on parametrization and constraints of the reference model and alternative models, as well as model initiation and run details, are presented in the Supplemental Material.

The extent and composition of syn–carbonate diapirism partial melts predicted to be found above the most volatile-saturated solidus in Figure 3 were modeled with the pMELTS algorithm (Ghiorso et al., 2002) between 1 GPa and 3 GPa over a temperature range of 900–1300 °C at nickel–nickel oxide (NNO) oxygen fugacity. Starting composition consisted of dry (Hirose and Kushiro, 1993) and wet (0.5 wt% H₂O) peridotite (average natural lherzolite with no carbonate) mixed with 2.5% to 15% CaCO₃ in increments of 2.5%. Each mixture was run at constant pressure, starting with 1 GPa in steps of 0.5 GPa with 3 °C increments from initial to final temperature. Carbonate in the dry mix produces a gradual increase in melt fraction to as much as 12% at 10% carbonate in the mix, especially in the low pressure–high temperature region. However, melt fraction does not exceed 4% between 1.5 and 2.5 GPa. The presence of 0.5% H₂O further increases melt fraction, although not significantly.

Within the model domain, material movement is driven by both the kinematically imposed convergence and internal buoyancy forces arising from thermal and compositional variations. Figure 2, and Animation S1 in the Supplemental Material, show the evolution of the model. The model starts at 0 m.y., following initiation of subduction (Fig. 2A). The incoming carbonates accumulate at the plate margin, initially forming a wide wedge, and with continued plate convergence, the carbonates are carried into the mantle by the subducting plate (Fig. 2B). At ~14 m.y., the carbonate layer in the mantle becomes unstable and begins to buoyantly detach from the oceanic plate as a diapir that initiates at 80–100 km depth (Fig. 2C). Similar results are obtained if the thickness of the carbonate is larger (we performed model runs with carbonate as much as 2 km thick [not pictured in Figure 2]). Diapirism also takes place at lower density contrasts between carbonate and the upper mantle, to as little as 250 kg/m³ contrast (not pictured). Thicker upper-plate lithosphere would lead to delayed and deeper detachment, but models show that diapirs always form under the general scenario envisioned in our study. The upwelling diapir rapidly ascends through the mantle wedge and penetrates the continental mantle lithosphere (Fig. 2D). Ongoing subduction creates continuous diapirism that results in (1) a ponding of the carbonates either at the lithosphere-asthenosphere boundary or at the continental Moho, and (2) cooling of the mantle wedge (Figs. 2E and 2F). Additional models show that carbonate diapirism and accumulation in the mantle wedge corner readily occur for a wide range of conditions (Fig. S3).

Our models predict that carbonate diapirs mix with the convective mantle and cool the mantle wedge during subduction before being added to the bottom of the continental lithosphere. A weak upper-plate lithosphere model allows for more thorough mixing of the existing mantle lithosphere (some of which is engaged in a delamination process from the backarc side) than a strong lithosphere model. For a strong lithosphere, the asthenospheric peridotite-carbonate mix does not penetrate the existing mantle lithosphere, but rather attaches to its bottom and adds to the total thickness of mantle lithosphere.

Seismic velocity calculations demonstrate that over a realistic range of pressures and temperatures for the uppermost mantle, peridotite-carbonate mixes produce a seismically slow mantle, with as much as 3% lower velocities in the lowermost lithosphere relative to asthenospheric mantle. It is unlikely that such a low-density material would be prone to subsequent convective removal, so this may add carbonate to the bottom of, or into, the mantle lithosphere.

Vertical thermal profiles above where the slab is located at 100 km depth at 14 and 23 m.y. after subduction initiation (Fig. 2) are shown in Figure 3 for the mantle wedge compared to the carbonated fertile peridotite solidus and the H₂O + CO₂–saturated peridotite solidus (Falloon and Green, 1989). This diagram shows that carbonate diapirism in the mantle wedge for relatively short-lived subduction systems retards melting because cold materials are being mixed into the mantle wedge. The H₂O + CO₂–saturated peridotite solidus is depressed relative to the water-only saturated peridotite, but the addition of cold materials from the downgoing slab counteracts that in terms of melt productivity. Given that the fluid-saturated fertile peridotite solidi shown in Figure 3 indicate the initiation of partial melting at the lowest temperatures possible and that most likely the actual peridotite in the wedge is less fertile (i.e., less clinopyroxene rich) and is likely not volatile saturated, and very likely not water saturated, the amount of melting in the mantle wedge of these systems is likely insignificant if carbonate material is introduced diapirically. The maximum predicted partial melt fraction under saturated conditions (H₂O and CO₂) is ~4%. For a near-pure limestone diapir, there is no water involved, but we assume that the presence of minor silicic sediments can add some water to the system, otherwise melting would not occur. Over time, after consumption of the carbonate platform, the wedge heats up and regains a hotter thermal regime by additional convection, and conductive heating of the carbonate diapir and its attachment to the bottom of the lithosphere take place (see below).

Forward melting modeling results show that minor amounts of partial melts are either calcio-carbonatites at 2.5 GPa or highly undersaturated feldspathoid-bearing leucitite and phonolites at 1.5 GPa. It is unlikely that the mantle wedge was ever less than ~45 km thick in a Tethyan scenario, given that the upper plates of such subduction systems were continental. Our models suggest that under a scenario of subduction of a sizable amount of carbonate, melting is retarded in most situations, which is consistent with the paucity of magmatism in places like the European Alps. Rather than producing a normal-size magmatic arc, such regions (the Carpathian Mountains, Alps, etc.) have either (1) no magmatism at all for long periods in times of convergence, or (2) very limited and alkali-rich magmatism that may be explained by the processes described here.

The pMELTS algorithm is not particularly accurate and at low melt fraction; a better path forward in understanding the extent and composition of melts is via experiments. Recently, Chen et al. (2021) provided some preliminary experimental data on precisely the process and physical conditions we envision here, and they obtained remarkably similar results to the pMELTS predictions (see the Supplemental Material).

We quantify how much CO₂ returns to the lithosphere in the upper plate by solid-state diapirism versus the amount of carbonate lost to the forearc and the fraction that is delivered to the deeper mantle along the subduction plane. A large fraction, ~35%, remains in the forearc; ~10% descends deeper into the mantle. The return via diapirism to the upper-plate lithosphere or the bottom of the upper plate is predominant; it appears to make up approximately 40%–50% of the subducted carbonate mass in the models presented here. This mechanism is important, needs to be included in the mass-balance budgets of subducted CO₂ (Plank and Manning, 2019), and predicts that the shallow supra-subduction mantle can be an important long-term sink for CO₂ under this scenario. For example, if one considers a 500-km-long carbonate (100% calcite) platform being subducted at a rate of 5 cm/yr and that ~40% of the subducted material returns to the upper plate through diapirism, as envisioned here, ~5 Mt of carbon per year are added to the upper plate from only one large platform.

Future tests need to quantify the amount of mantle-derived materials that may have isotopic signatures indicative of carbonate additions in regions like the Alpine-Carpathian chain; recent basaltic fields exist to test for that (Wilson and Downes, 2006). Ca isotope values (DePaolo, 2004) in sedimentary carbonates are typically distinct from mantle values (Amsellem et al., 2020). Mantle lithospheric xenoliths, also common regionally (Downs 2001), can be investigated for interaction with cryptic carbonate liquids using the techniques put forward for other regions (Chen et al., 2018). Frozen calcio-carbonatites in xenoliths from the Tethyan realm, like those described in detail from the Persani volcanic field (Romania; Chalot-Prat and Arnold, 1999), can be studied to gain information on the origin of carbonate, the magnitude of carbonate recycling, and the timing of percolation of the lithospheric mantle.

We acknowledge support from a Babes-Bolyai University (Romania) Fellowship (grant CNFIS-FDI-2021–0061) and Romanian UEFISCDI (Executive Unit for the Financing of Higher Education, Research, Development and Innovation) grant PN-III-P4-ID-PCCF-2016-0014. C. Currie acknowledges support and computing resources from WestGrid (Victoria, British Columbia, Canada) and Compute Canada (Toronto).