Local rapid exhumation during the long-lived Grenville orogeny

Duration is an essential dimension of geological processes. For orogenic metamorphism, in particular, the rates of heating and burial are closely linked to the kinematics of plate convergence. Convergence in continental collisions typically leads to regional metamorphism over time scales longer than 10 m.y. (Viete and Lister, 2017). The 10 m.y. lower limit is not arbitrary, as it is the time scale for heat conduction associated with a thermal anomaly for which the characteristic dimension is ~15 km in orogenic crust. In contrast, fast metamorphism (<1–10 m.y.) records local geologic events in response to lithospheric processes beyond orogenic thickening (England and Thompson, 1984), which can facilitate a holistic interpretation of tectonic evolution (Viete and Lister, 2017).

The emergence of modern-style plate tectonics characterized by one-sided subduction was a result of a change in the underlying geodynamic environment, which was also manifested in the style of metamorphism (Palin et al., 2020). The proposed ages for the onset of modern-style plate tectonics cover a broad spectrum from the Hadean to the Neoproterozoic, encompassing the entire Precambrian (e.g., Florence and Spear, 1995; Komiya et al., 1999; Stern, 2005; Cawood et al., 2006; Hopkins et al., 2008; Palin and White, 2016; Pauly et al., 2016; Stern et al., 2016; Tang et al., 2016). Metamorphism in Precambrian orogens has been shown to have proceeded more slowly than Phanerozoic metamorphic events (e.g., Scibiorski et al., 2015; Chowdhury et al., 2021; Zou et al., 2021, and summary therein). The appearance of rapid metamorphism and the shift to a bimodal distribution of metamorphic pressure-temperature (P-T) conditions (Holder et al., 2019) could reflect the transition to modern-style plate tectonics (Chowdhury et al., 2021). As such, quantification of metamorphic time scales, especially for Precambrian rocks, is helpful in further constraining the age of the onset of modern-style plate tectonics.

However, quantification of metamorphic rates is a long-standing problem. It is particularly challenging in Precambrian terranes, and few examples of Precambrian short-duration metamorphism have been documented (e.g., Guevara et al., 2017; Zou et al., 2020). The rarity of documented ancient rapid metamorphism could reflect the distinct tectonic regimes of Precambrian Earth, or it could be an artifact caused by the limited temporal resolution of current radiometric geochronology, or a combination of both (Viete and Lister, 2017; Zou et al., 2020). One potential approach to overcome this challenge is diffusion speedometry, the temporal resolution of which is independent of the age of the rock. Numerous studies have successfully constrained the time scales of high-grade metamorphism by applying garnet diffusion speedometry to Precambrian rocks up to Archean in age (e.g., Ague and Baxter, 2007; Raimbourg et al., 2007; Guevara et al., 2017; Chu et al., 2018; Chowdhury and Chakraborty, 2019; Zou et al., 2020).

The Grenville Province of North America (Fig. 1A) records a Mesoproterozoic continental collision during the assembly of the supercontinent Rodinia. Compared with orogenies in other supercontinental cycles, Grenvillian magmatism was distinctively enriched in Nb, Y, and Zr, indicative of a high-temperature intraplate signature (Moecher and Samson, 2006; Liu et al., 2017). Anorthosite massifs, peralkaline intrusions, dike swarms, extensive A-type granite intrusions, and massive granulite-facies terranes distinguish the Grenville Province from Phanerozoic orogens (e.g., Rivers et al., 2012). This abnormally hot orogeny might reflect prolonged heating of the subcontinental mantle beneath a long-lived supercontinental insulating lid (Cawood and Hawkesworth, 2014). Geochemical proxies of crustal thickness reveal cycles of crustal thickening and thinning during the Ottawan phase (Brudner et al., 2022). The thinning events could have been driven by gravitational collapse of weak overthickened lithosphere (Rivers, 2012), thermal erosion of the lower crust (Corrigan and Hanmer, 1997), or delamination of the lower crust (Lieu and Stern, 2019). Although the exact mechanism remains unclear, these repetitive thinning events resulted from rheological weakening of the crust by the hot underlying mantle.

The metamorphic histories of crustal high-pressure (HP) rocks serve as a field-based benchmark for the inferred kinematics of orogenic thickening and thinning during the Ottawan orogeny (Brudner et al., 2022). Pods and lenses of retrogressed eclogite with early Ottawan ages constitute the Western Grenville high-pressure belt in southern Ontario, Canada (Fig. 1B; Ketchum and Davidson, 2000; Rivers et al., 2002). These rocks are interpreted to have been uplifted from as deep as 50–60 km depth (Hynes and Rivers, 2010). One of these occurrences is the massive metabasite cropping out as discrete masses and ridges in the highly deformed, felsic to intermediate gneisses of the Mattawa domain (Figs. 1B and 1C). Previous studies have reported garnet-symplectite assemblages that have been interpreted as relicts of eclogite-facies metamorphism (Easton, 2006; Tom-Ying, 2015; Cao et al., 2021). The plagioclase + diopside symplectite is recognized as a common decompression texture in HP rocks and is inferred to be pseudomorphs after omphacite or omphacite + phengite (Anderson and Moecher, 2007; Chu et al., 2016). Cao et al. (2021) identified various melt-related microstructures and constructed a P-T path from peak pressure of 1.7 ± 0.1 GPa at 680 ± 50 °C to peak temperature of 920 ± 50 °C at 1.3 ± 0.1 GPa, further confirming peak metamorphism at eclogite facies. Despite the high temperatures that allowed partial melting, garnet in the metabasites still exhibits distinctive relicts of compositional growth zoning that is rarely preserved in granulite-facies rocks. This suggests a transient peak metamorphic pressure stage followed by rapid exhumation and cooling.

In this study, we examined the petrographic relationships of the least-retrogressed eclogite and conducted systematic analyses of mineral compositions. We integrated the P-T evolution derived from phase equilibria analysis and various thermometers with the time scales estimated based on diffusion modeling of zoned garnet. The metamorphic time scales were then used to constrain the rates of decompression, on the basis of which we further investigated the implications for the tectonic framework of the western Grenville Province.

The Grenville Province on the southeastern margin of the Canadian Shield hosts lithologic units of diverse ages and metamorphic grades (Fig. 1A). Metamorphic terranes of Grenvillian age extend through Northern Ireland and Scotland to the Sveconorwegian Province in Scandinavia (Davidson, 1998; Lorenz et al., 2012). Based on its metamorphic and structural architecture, the Grenville orogen has been considered to be a Mesoproterozoic prototype of the continental collision resulting in the Himalayan-Tibetan orogen (e.g., Hynes and Rivers, 2010; Indares, 2020).

The subdivision of the Grenville Province has been variably revised across provincial and geological borders. To be consistent with relevant publications and survey reports, we retain the old terminology in reference to western Grenville in Ontario and western Quebec (Wynne-Edwards, 1972) and the lithotectonic division of the Central Gneiss belt in the western Grenville Province by Ketchum and Davidson (2000). The Central Gneiss belt (Fig. 1B) is composed of amphibolite to granulite gneisses and hosts some of the most deeply seated rocks with Archean to Mesoproterozoic protolith ages (Easton, 1986). The SE-dipping Allochthon Boundary thrust, which is a deformation zone within the Central Gneiss belt, was interpreted as a major structural boundary dividing the Central Gneiss belt into (par-) autochthonous and allochthonous belts (Ketchum and Davidson, 2000). The parautochthonous terranes in the footwall of the Allochthon Boundary thrust were derived from Laurentian cratonic and supracrustal rocks (Easton, 1986). The allochthonous terranes comprise exotic thrust slices metamorphosed to greenschist, amphibolite, and granulite facies (Easton, 2000). They were transported laterally and tectonically overlie the parautochthonous belt. The structurally lowermost units of the allochthonous stack are collectively known as the Algonquin–Lac Dumoine terrane (Fig. 1B).

The Grenvillian collisional events were likely superimposed on preexisting active magmatic arcs and metamorphic terranes of the Pinwarian (ca. 1495–1445 Ma), Elzevirian (ca. 1245–1225 Ma), and Shawinigan (ca. 1140–1180 Ma) orogens. Multiple phases of deformation and metamorphism could account for the vast extent of high-grade terranes, and active arc and back-arc magmatism prior to the beginning of continental collision would have contributed to a significantly enhanced geothermal gradient in the lower crust (Rivers and Corrigan, 2000; Rivers et al., 2012). Within the broadly defined Grenville orogeny, two distinct pulses of continental collision, or orogenic phases, have been identified. The main collision Ottawan phase (1090–1020 Ma) was characterized by high-grade metamorphism, coeval magmatism, and a substantial amount of lateral transportation (allochthonous). The subsequent Rigolet phase (1005–980 Ma) was characterized by fold-and-thrust belts and plutonism following postconvergence gravitational spreading at the end of the Ottawan phase (Fig. 1A; Rivers, 1997; Rivers, 2009; McLelland et al., 2010; Indares, 2020). The latter Rigolet phase mainly affected parautochthonous foreland terranes in the footwall of the Allochthon Boundary thrust (Figs. 1A and 1B).

Occurrences of relict eclogite in the vicinity of the Allochthon Boundary thrust on the hanging-wall side, constituting the Western Grenville high-pressure belt (Fig. 1B; Rivers et al., 2002), have been documented in the Algonquin–Lac Dumoine terrane (Fig. 1B; Davidson, 1998; Rivers et al., 2002). Although these metabasites underwent extensive retrogression, omphacite in metabasite (Jamieson et al., 2003), plagioclase + clinopyroxene symplectite after omphacite (e.g., Davidson et al., 1982), and corundum + spinel + sapphirine + plagioclase intergrowth after kyanite (Davidson, 1990) indicate high-pressure metamorphic conditions. Indares and Dunning (1997) studied the retrogressed eclogite from Lac Dumoine terrane in western Quebec and determined the condition of granulite-facies overprinting (1.2–1.4 GPa at ~750 °C) on the basis of the interpretation of textures and mineral thermobarometry. They also interpreted the absence of retrograde zoning in metamorphic minerals due to “rapid” tectonic exhumation. The HP metabasites from the Shawanaga domain (Fig. 1B) record higher peak pressure of >1.5 GPa at ~750 °C, according to phase equilibrium modeling and Ti-in-zircon and Zr-in-rutile thermometers (sample 119 in Marsh and Culshaw, 2014). Similarly, pseudosection analyses applied to the HP metabasite from the Novar domain indicated peak pressure of ~1.5 GPa at ~750 °C and a clockwise P-T path (Marsh and Kelly, 2017). These metabasites yielded broadly Ottawan ages of eclogite-facies metamorphism (1.1–1.07 Ga; Indares and Dunning, 1997; Ketchum and Krogh, 1998; Rivers et al., 2002; Marsh and Culshaw, 2014; Marsh and Smye, 2017). In particular, five of six samples studied by Marsh and Culshaw (2014), from various domains of the western Grenville Central Gneiss belt, yielded zircon U-Pb ages between 1097 and 1085 Ma, so the likely time frame for HP metamorphism in the western Grenville Province is better constrained to the early Ottawan (ca. 1090 Ma; Figs. 1A and 1B).

The Mattawa domain is a minor component of the Algonquin terrane that structurally overlies the Tomiko terrane along the Allochthon Boundary thrust (Figs. 1B and 1C; Dickin and Guo, 2001; Dickin et al., 2008). This major structural boundary is supported by the protolith Nd model age difference between the Mattawa domain (1750 ± 55 Ma, 2σ) and the underlying parautochthonous units (2.6–1.9 Ga; Dickin and Guo, 2001; Easton, 2006). Geophysical evidence suggests that the Mattawa domain is a thin structural sliver emplaced over the Tomiko terrane, and the majority of it has been eroded (Easton, 2006). Major rock types of the Mattawa domain include quartz monzonite plutons, felsic to intermediate gneisses, and pods or layers of metagabbroic rocks associated with meta-anorthosite (Easton, 2006). The felsic gneisses are mainly fine-grained tonalitic to granodioritic with the mineral assemblage quartz + plagioclase + K-feldspar + biotite + garnet (Easton, 2006). The protolith signatures are volcanic-arc granitoids, differing from the A-type granitoids of the Tomiko terrane (Easton, 2006). The metabasites, resembling high-Fe tholeiite in composition (Easton, 2006), crop out as discrete masses and dikes within intermediate to felsic gneisses in the study area, suggesting that they were likely derived from mafic enclaves or dikes in the felsic country rock. The metagabbro consists of garnet porphyroblasts with polycrystalline plagioclase coronae hosted in a plagioclase + clinopyroxene symplectite matrix with various degrees of amphibolite-facies retrogression; this assemblage is interpreted to have been derived from an eclogite-facies precursor (Easton, 2006; Tom-Ying, 2015; Cao et al., 2021).

Quantitative wavelength-dispersive spectrometer (WDS) analyses, energy-dispersive spectrometer analyses (EDS), and backscattered electron (BSE) imaging were carried out using a JEOL JXA-6230 electron microprobe at the University of Toronto. The operating conditions were 15 kV (garnet) or 10 kV (other minerals) accelerating voltage and 10 nA (hydrous phases, plagioclase, and pyroxene) or 50 nA (garnet and rutile) beam current. Amphiboles and feldspars were analyzed using a defocused electron beam (5 μm) to minimize beam damage. The chemical maps were made at 15 kV, 100 nA, and 200 ms dwell time, with a 5 μm spatial resolution. Trace elements in rutile were counted for 100–200 s on peak to optimize counting statistics. Mineral standards used were well-defined natural and synthetic minerals. Representative mineral compositions are given in Tables 1–7.

Bulk-rock compositions were measured via X-ray fluorescence (XRF) on a Shimadzu XRF-1500 spectrometer at the Institute of Geology and Geophysics, Chinese Academy of Sciences (IGGCAS). Data were calibrated using the USGS AGV-1 standards, and the analytical precision was better than 5%. Prior to analyses, 0.5 g aliquots of sample powders were heated to ~1000 °C for 1.5 h in a furnace, and the loss on ignition (LOI) contents were calculated. Due to the oxidation of FeO to Fe₂O₃, calcination led to mass gains and negative LOIs. We recalculated volatile losses by assuming all iron oxides were FeO. These results are given in Table 8.

The studied samples were collected along Ontario Highway 553, north of Mattawa, Ontario, Canada (Fig. 1C). The outcrop is a road cut with an ~70 m intermittent exposure of metabasite. Felsic to intermediate gneisses crop out tens of meters to the north of the mafic outcrops, but the lithologic contact is not exposed. The metabasite is massive, with sparse felsic veinlets interpreted as leucosome (Fig. 2A) and dark nodules of biotite and amphibole aggregates. Samples XCG6A–XCG6D (e.g., Fig. 2B) were all collected from the same massive outcrop spanning ~10 m × 10 m. The samples were nearly identical with respect to mineral assemblages and textures. Cao et al. (2021) collected their sample 17-17 from the same outcrop.

The least-retrogressed metabasite consisted of garnet porphyroblasts with negligible polycrystalline plagioclase coronae (Fig. 2C) hosted in a matrix of lamellar clinopyroxene + plagioclase symplectite (Figs. 2D and 3). Sporadic dark amphibole + biotite patches were present (Fig. 2C); these large flakes of amphibole and biotite were randomly oriented. Smaller amphibole grains crosscut clinopyroxene and plagioclase grain boundaries (Figs. 2D and 3; Fig. S2¹). Accessory phases in the matrix included rutile, ilmenite, apatite, and minor epidote. Minor calcite and titanite were also present as inclusions in garnet (Cao et al., 2021).

Other metabasites showed more dominant amphibolite-facies overprinting in this area, with few relict eclogitic textures and minerals preserved. Plagioclase and amphibole replaced the dark-green symplectite matrix, resulting in a granoblastic texture. As a result, the matrix appeared to be much darker in more-retrogressed outcrops (Fig. 2E). Garnet was commonly resorbed and rimmed by plagioclase coronae. In more retrogressed samples, the plagioclase coronae were visibly thicker and showed a weak orientation of elongation (Fig. 2E). Garnet was almost completely resorbed in extreme cases (Fig. 2F).

The symplectite consisted of fine-grained lamellar intergrowths of clinopyroxene + plagioclase and late granular amphibole (Fig. 2D). The symplectite lamellae varied both in size and shape (Fig. S1). Well-preserved lamellae with uniform thickness and straight grain boundaries were 150–400 µm in length and 5–50 µm in thickness, whereas smaller grains appeared more vermiform (Fig. 3; Fig. S1). On average, the relative modes of plagioclase, clinopyroxene, and amphibole were 50 vol%, 38 vol%, and 12 vol%, respectively, estimated based on a series of BSE images using the image analysis software ImageJ (Fig. S1; Table 4A).

The compositions of the clinopyroxene in the symplectite intergrowth ranged from diopsidic to omphacitic with Na/(Na + Ca) ratios of 0.16–0.29 (Fig. 4A; Table 1). The Na content slightly decreased toward grain boundaries (Fig. 3D2; note the greener grain boundaries). The plagioclase in the symplectite was mostly oligoclase with an average An# of 16% (= Ca/[Ca + Na] × 100%; Fig. 4B; Table 2). Small euhedral grains of amphibole crosscut the clinopyroxene-plagioclase lamellae (e.g., Fig. 3D). With an average of 6.1 Si per formula unit (p.f.u.) and 0.6 Mg/(Mg + Fe²⁺) (Table 3), they plotted in the pargasite/tschermakite field according to the classification of Leake et al. (1997) (Fig. 4C).

Garnet grains were measured to be 1–5 mm in diameter. The garnet inner cores were irregular dark patches (Fig. 2C); larger garnet porphyroblasts grew over multiple dark cores. Chemical mapping was conducted on several representative garnets; relatively large garnet grains were selected to minimize the effect of retrogressive resorption. Three representative garnets from samples XCG6B and XCG6D are presented in Figure 5 and in Figures S2 and S3 (see ^{footnote 1}). The zoning pattern of the garnet in sample XCG6B was concentric, which was especially apparent in the partially preserved bell-shaped growth zoning of Mn (Fig. 5F). The X_SPS (= Mn/[Ca + Mg + Mn + Fe²⁺], where SPS indicates spessartine) of the core reached 0.018 and smoothly decreased to 0.008 toward the rim (Figs. 5H and 5I). The garnet core showed a patchy Ca pattern (Fig. 5D) and smoother Mg profiles (Fig. 5E). The Ti composition map showed a smaller poikilitic inner core with abundant rutile inclusions (Fig. 5G), matching the shape of the dark garnet core (Fig. 2C), the patchiness of the Ca map (Fig. 5D), and the area of high Mn concentration (Fig. 5F). The core was overgrown by a mantle characterized by high Ca, Mg, and Ti, where X_GRS (= Ca/[Ca + Mg + Mn + Fe²⁺], where GRS indicates grossular) reached as high as 0.2–0.24 (Figs. 5H and 5I; Table 5) at the euhedral boundary between the core and mantle (Fig. 5D). This euhedral trace of the X_GRS peak and outer rim (Fig. 5D) matched that of the high-Ti overgrowth (Fig. 5G). The Mg# (= Mg/[Mg + Fe²⁺]) profile increased from 0.25 to 0.35 rimward, with a minor downturn to 0.31 at the rim (Fig. 5H).

Garnet grains hosted variable amounts of rutile, monazite, ilmenite, and apatite inclusions. Polymineralic inclusions occurred near the garnet rims, and rarely in the cores. Garnet grains with well-preserved zoning patterns were largely free of polymineralic inclusions, and garnet grains with polymineralic inclusions were commonly strongly modified or absorbed. These polymineralic inclusions occurred in or close to the high-Ca mantle (fig. 6 in Cao et al., 2021). The polymineralic inclusions consisted of clinopyroxene, amphibole, quartz, titanite, and chlorite (Figs. 5B and 5C). Cao et al. (2021) reported polymineralic inclusions of plagioclase + K-feldspar ± biotite ± ilmenite ± amphibole ± biotite. Calcite was also found in polymineralic inclusions.

Amphibole was present as (1) garnet inclusions (Figs. 5B and 5C; Fig. S3D [see ^{footnote 1}]; Cao et al., 2021), (2) crosscutting phases in symplectites (Figs. 2D and 3), and (3) patchy prisms with biotite flakes (Figs. 2C and 2D). The latter two reflect variable degrees of retrogressive hydration of the rock. The small amphibole grains crosscutting the symplectite and the large amphibole grains in the matrix had overlapping compositions and plotted in the pargasite/tschermakite field (Fig. 4C; Table 3).

Biotite was present as a matrix phase, commonly associated with amphibole and ilmenite. The matrix biotite flakes were millimeter sized and showed good basal cleavage. Biotite had Mg# = 0.67–0.72 and was relatively rich in Ti (0.28–0.35 p.f.u.; Table 6). The highest Ti contents were found in the flakes where ilmenite was present. In addition, minor biotite was also found in the polymineralic inclusions in garnet porphyroblasts (Cao et al., 2021).

Plagioclase was associated with (1) garnet coronae (Fig. S2E [see ^{footnote 1}]), (2) symplectite intergrowths (Fig. 3), and (3) perthitic intergrowths with K-feldspar in polymineralic inclusions in garnet (Cao et al., 2021). The plagioclase lamellae in the symplectite were broadly oligoclase in composition (Fig. 4B) with an average An# of 16 (Table 2). The plagioclase in the coronae displayed a decreasing Na trend from the plagioclase interior toward the plagioclase-garnet grain boundary, representing the garnet resorption front (Fig. S2E).

Rutile was present in the assemblage both as a matrix phase and as inclusions in garnet. Rutile grains in the matrix were commonly rimmed by ilmenite as a result of decompression (e.g., Clarke et al., 1989). The inclusions were mostly subhedral and varied in size from ~10 to 100 µm (Fig. S2). Rutile was also found in polymineralic inclusions in garnet with apatite, quartz, and epidote. We filtered the electron probe microanalysis data from rutile inclusions in garnet using the method summarized in Ague and Eckert (2012). The Zr contents ranged from 688 to 1725 ppm in the rutile inclusions (Table 7), displaying two peaks in the histogram: 982 ± 25 ppm and 1539 ± 46 ppm (Fig. 6A). The Cr content was slightly higher than that of the host garnet. Aside from one outlier with an exceptionally high Nb content of 1265 ppm, there was no obvious pattern on both the Nb-Zr plot and Cr-Zr plot (Figs. 6B and 6C).

The evolution of mineral assemblages and partial melt in samples from the same outcrop is described in detail in Cao et al. (2021). This section focuses on the correlation between mineral compositions and the inferred metamorphic stages. The interpretations are then used to construct the P-T history.

The early prograde history is poorly preserved. Mg (Fig. 5E) and Fe contents of the cores are more homogenized than Ca, but the Mn map still shows a pronounced bell-shaped zoning pattern. Given the granulite-facies metamorphic condition and fast diffusion of Mn in garnet (e.g., Chu and Ague, 2015), such a Mn pattern suggests that the overall duration of the high-temperature stage must have been relatively brief (Caddick et al., 2010). The garnet core borders a Ca- and Ti-rich mantle with a euhedral boundary (contour 1 in Fig. 7A). Where the mantle overgrowth is thicker, X_GRS reaches 0.24 (Figs. 5I and 7A1). We speculate that the whole overgrowth was as rich in Ca, and the thinner mantle was then subject to diffusion modification. The Ca content decreases toward the rim; the Mg# peak (“2” in Fig. 7A2; up to 0.35) is closer to the rim than the Ca peak. Mg# of garnet typically increases as metamorphic grade increases. Thus, we regard the mantle as growth zoning. This growth zone is also characterized by higher Ti contents (Fig. 7A3), likely resulting from the decomposition of Ti-bearing phases (e.g., biotite) during prograde melting (Fig. 7C).

Minerals from the prograde-to-peak assemblage, including tschermakitic amphibole (e.g., Fig. S2; Cao et al., 2021) and rutile, are found in garnet inclusions (Fig. 7D). The polymineralic inclusions, characterized by euhedral negative shape or sharp cusps into garnet, are interpreted as small parcels of melts included in garnet at supersolidus conditions (Figs. 7C and 7D). The proximity of polymineralic inclusions to the high-Ca mantle, or the garnet rim where zonation is not preserved, suggests that the corresponding zones of garnet grew in supersolidus conditions. Except for minor plagioclase in solidified melt inclusions, plagioclase is exclusively found as coronae surrounding resorbed garnet or lamellae in symplectite, so plagioclase was not stable in the peak-pressure assemblage (Fig. 7D).

Symplectite in the matrix is inferred to have replaced Na-rich clinopyroxene precursors. The decomposition reaction is largely isochemical (Anderson and Moecher, 2007). In this study, we used the modal abundances and average compositions (Fig. 3) of clinopyroxene, plagioclase, and amphibole to calculate reintegrated compositions to approximate that of the original clinopyroxene. The water in amphibole was discarded as rehydration during symplectite formation or later amphibole replacement. Anderson and Moecher (2007) speculated that the source of hydroxyl in symplectitic amphibole is potentially derived from hydroxyl stored in M₂ vacancies of Ca-Eskola. Based on this speculation, the H₂O content of the clinopyroxene precursor would still be low (<0.15 wt%). The reintegrated compositions are close to the clinopyroxene formulae (Table 4; for another example, see Chu et al., 2016) and suggest that the precursors were omphacite with a considerable jadeite component (Na/[Na + Ca)] = 0.49 ± 0.02; Fig. 4A). If the amphibole component were ignored in the reintegration, the precursor SiO₂ content (~60 wt%) would be too high for a clinopyroxene group mineral. Because the amphibole grains have lower Na/(Na + Ca) ratios (0.28–0.31; Table 3), the compositional integration ignoring amphibole would lead to even higher jadeite fractions of the omphacite precursor. In summary, the peak assemblage consisted of garnet + omphacite + amphibole + rutile ± biotite with no plagioclase.

During decompression from eclogite-facies to granulite-facies conditions, garnet broke down partially, as evident by resorbed margins (Fig. 7A). In the matrix, homogeneous omphacite decomposed to more calcic diopside and more sodic plagioclase, forming symplectitic intergrowths in relatively anhydrous conditions (“3” in Fig. 7B). Retrograde amphibole crosscut clinopyroxene and plagioclase grain boundaries and partially replaced symplectite (“4” in Fig. 7B; Fig. 7D). The partial melt solidified as the P-T path crossed the solidus, and the solidifying melt released water to drive further retrogressive hydration. Where water was locally abundant in the assemblage, amphibole + biotite domains replaced the matrix omphacite or symplectite.

The degree of resorption varies among garnet grains. The garnet grains in another sample (XCG6D) are less euhedral in shape (Figs. S2 and S3) and show greater resorption. The garnet grains surrounded by thicker plagioclase coronae show more pronounced Mn upturns at the rims (Fig. S2C and S3C). In Garnet_1 in XCG6D (Fig. S2), an off-center core was revealed by the chemical maps of Mn (Fig. S2C). Although the peak Ca content (X_GRS = 0.176) is lower than in XCG6B (Fig. 5), the general zoning patterns are similar. Garnet_3 in sample XCG6D (Fig. S3) is more homogenized and shows a decrease in Mg and Ca and an increase in Fe and Mn toward the rim (Fig. S3F). A ring of high-Ca mantle of variable thickness (Fig. S3A) matches other garnets.

Rutile is a high-pressure phase, which breaks down to form ilmenite as pressure decreases (Fig. 7C). Both rutile and ilmenite were found in garnet inclusions as well as matrix phases. Most rutile inclusions in garnet appeared pristine; others coexisted with ilmenite. When closely inspected, some of the “pristine” rutile showed thin (<1 μm) lamellae that were brighter in BSE. These lamellae are likely nanoscale retrograde ilmenite or magnetite, although these rutile data do not necessarily show higher Fe contents. These rutile inclusions are more susceptible to Zr loss; on the other hand, ilmenite excludes Zr, so rutile near ilmenite replacement could be enriched in Zr (Whitney et al., 2015). Trace-element contents of the rutile inclusions showed no relationship with the distance to the garnet core. During prograde metamorphism, garnet grows and captures newly formed or recrystallized rutile grains, as well as detrital or inherited rutile grains. The presence of detrital or inherited rutile grains may disturb the correlation of rutile Zr contents with garnet growth zonation. During retrograde metamorphism, Zr in all the rutile grains tends to diffuse (e.g., Štípská et al., 2014; Usuki et al., 2017). Potential Zr loss has little correlation with the distance to the garnet core; it is mainly controlled by the density of garnet cracks and defects. Thus, the higher-Zr peak (1539 ppm; Fig. 6A) should be considered as a lower limit of the Zr content of rutile in the peak assemblage, and we speculate that the lower-Zr peak (982 ppm; Fig. 6A) represents a mixture of reset rutile and inherited grains that survived high-temperature metamorphism.

Based on the interpretations of the mineral assemblages and mineral zonation, we modeled the phase relations to assess the P-T evolution more quantitatively. The phase equilibria modeling employed THERMOCALC (version 3.45; Powell and Holland, 1988), the internally consistent thermodynamic data set ds62 (Holland and Powell, 2011), and compatible activity models. The activity-composition models used in this study were: tonalite-trondhjemite melt, sodic-calcic clinopyroxene (omphacite), and clino-amphibole (Green et al., 2016); garnet, biotite, orthopyroxene, and ilmenite (White et al., 2014a, 2014b); plagioclase and alkali feldspar (Holland and Powell, 2003); and spinel-magnetite group (White et al., 2002). The mineral abbreviations are after Holland and Powell (2011).

The pseudosection was calculated in the model chemical system of Na₂O-CaO-K₂O-FeO-MgO-MnO-Al₂O₃-SiO₂-H₂O-TiO₂-Fe₂O₃ (NCKFMMnASHTO), for which the bulk-rock composition was estimated from the results of X-ray fluorescence (Table 8). The majority of P was hosted in apatite, and the corresponding Ca was removed from the bulk composition. Excess “O,” standing for Fe³⁺, was estimated based on the Fe³⁺ content of symplectite (Table 4) and its mode (44 vol%).

To investigate the metamorphic conditions recorded by the garnet mantle, we removed the garnet core components from the bulk composition (inside “1” in Fig. 7). The core components were estimated using the garnet mode (~44 vol%), core fraction (~0.72), and the average core compositions. Uncertainties arose from performing area to volume extrapolation and from variation between garnet grains, but the resulting effective bulk composition was not very sensitive to these uncertainties. The effective bulk composition was similar to the starting bulk composition except for lower FeO and MnO (Table 8). The bulk Mn content was negligible, and so was the modeled garnet X_SPS. The low Mn content made the garnet solid solution effectively ternary. Since ternary systems can be constrained with just two degrees of freedom, we plotted the garnet Mg# isopleth instead of the X_ALM and X_PPR isopleths [X_ALM = Fe²⁺/(Fe²⁺ + Mg + Mn + Ca), ALM for almandine; X_PRP = Mg/(Fe²⁺ + Mg + Mn + Ca), PRP for pyrope], as well as the X_GRS isopleth in the phase diagrams.

Cao et al. (2021) systematically described the partial melting textures in the retrogressed eclogite. The extent of melt fractionation posed uncertainties in the bulk composition. However, the outcrop showed a very low degree of melt segregation with minor veining. The felsic veinlets showed no external connectivity; thus, the bulk composition should have remained effectively untouched with such a low degree of melt segregation. The massive samples were collected far from these leucosome veins.

The water content was derived from the modified LOI. The bulk H₂O content made water barely saturated at the solidus at ~1.8 GPa (H₂O 4.67 mol%; Fig. 8, inset), where the isopleth of the garnet Ca peak (X_GRS = 0.22–0.24; “1” in Fig. 7A) intersect the solidus (Fig. 8, inset). As subsolidus prograde metamorphism is accompanied by dehydration, the assemblage should be saturated with H₂O at the solidus before free water is absorbed by the partial melting. As long as water is present along the solidus, the solidus temperature, garnet isopleths, and phase relations are largely insensitive to the bulk H₂O content (Fig. 8B). Although the presence of calcite suggests that the fluid might have been more complicated than pure water, the H₂O activity is primarily controlled by the equilibrium between hydrous melt and hydrous minerals, so other impurities like carbon species do not affect the supersolidus phase relations (Chu and Ague, 2013) and thus the peak-pressure condition.

We assumed a closed system, so the bulk-rock Fe³⁺/(Fe²⁺ + Fe³⁺) ratio [Fe³⁺/Fe_T] and the amount of H₂O would not change during the metamorphic evolution. The application of a pseudosection also assumed equilibrium in the P-T range of interest. Both are strong assumptions that are unlikely to have been true. Nonetheless, we note that the focus of the pseudosection thermobarometry was the metamorphic peak pressure, where chemical equilibrium is most likely archived (Powell and Holland, 2008). The variations of bulk H₂O and Fe³⁺/Fe²⁺ ratios within their reasonable ranges do not significantly impact the general estimate of the peak pressure (see below). On the other hand, the assemblage tended not to evolve in global equilibrium far past peak conditions (Powell et al., 2019), so the retrograde exhumation and cooling information from the P-T pseudosection should be viewed as qualitative.

The calculated phase equilibrium diagram for sample XCG6D is presented in Figure 8A. Most assemblages in the calculated P-T range are melt-bearing assemblages. The phase relations are heavily dependent on pressure, and the boundaries of stability fields intersect the temperature-sensitive solidus and SiO₂-saturation boundary. The inferred peak assemblage L (melt) + o + g + bi + ru + hb (see the Fig. 8 caption for mineral abbreviations) dominates a large P-T range (1.7–2.2 GPa, 700–900 °C). With decreasing pressure, the breakdown of omphacite is accompanied by the formation of plagioclase, the conversion from rutile to ilmenite, and the appearance of orthopyroxene. Garnet is a stable phase in the assemblages at >1 GPa, although its modal abundance is strongly dependent on pressure. Amphibole is stable in the entire calculated P-T field. All amphibole is denoted as hornblende (hb) for simplicity, but at high-P and high-T conditions, the modeled amphibole is tschermakitic and Ti-rich and occurs in negligible amounts. The observed retrogressed subsolidus assemblage of di + g + bi + ilm + hb + pl ± ru ± q without orthopyroxene suggests that the pressure does not drop to an opx-stable condition.

Relevant compositional isopleth contours are plotted in the P-T fields of interest. The Na/(Na + Ca) isopleths of clinopyroxene (omphacite or diopside) are subparallel to the temperature axis (Fig. 8A). The ratio of reintegrated symplectite (0.49 ± 0.02; Fig. 4A) corresponds to ~2.1 ± 0.2 GPa at 850 °C. From the solidus to the metamorphic peak condition, the garnet growth is accompanied by an increase in Mg# and a decrease in X_GRS (Fig. 8A). The garnet X_GRS isopleths are largely temperature-sensitive features. At ~1.8 GPa, the modeled temperature decreases from near the solidus with modeled X_GRS of 0.23 to 850 °C with modeled X_GRS of 0.175. In a plagioclase-free assemblage, clinopyroxene and garnet are the two major Ca carriers. While the Na/(Na + Ca) ratio of clinopyroxene is more dependent on pressure (Fig. 8A), the Ca content of garnet (X_GRS) is diluted as garnet grows with increasing temperature. The garnet Mg# increases toward higher temperatures and pressures. The isopleth of the highest Mg# of the garnet rim (~0.35; Figs. 5H and 5I) intersects that of X_GRS (0.175; Fig. 5H) at 2.1 GPa, 830 °C (“2” in Figs. 7A and 8). This P-T condition agrees with the clinopyroxene Na/(Na + Ca) ratio of the reintegrated symplectite (0.49; Figs. 4A and 8). The highest X_GRS of the garnet reaches 0.22–0.24 in the broad high-Ca zone (Fig. 7A1); these isopleths broadly intersect the solidus at ~700 °C, 1.8 GPa (“1” in Figs. 7A and 8). The mode (vol%) of garnet increases significantly (32%→39% in a garnet core–removed assemblage) from “1” to “2,” which explains the euhedral peritectic overgrowth mantle. The corresponding Mg# isopleth (~0.3) crosses the solidus ~40 °C higher, which may be the result of faster diffusion of Mg and Fe compared to Ca (Chu and Ague, 2015).

The decompression path sees a monotonic decrease in the garnet mode (right panel of Fig. 8A), which indicates the partial resorption of garnet. We note that garnet could still be stable in the assemblage because (1) the phase diagram uses a garnet core–removed bulk composition, (2) garnet resorption could be hindered by local equilibrium, as suggested by minor plagioclase coronae, and (3) the decompression path does not cross the garnet-out boundary. The replacement of omphacite with diopside and a decrease in the cpx Na/(Na + Ca) ratio also take place during decompression, corresponding to the formation of di + pl symplectite. To assess how the variation in bulk H₂O and Fe³⁺/Fe_T ratio affected the pressure estimate, we also calculated P-H₂O and P-Fe³⁺ pseudosections at 850 °C (Figs. 8C and 8D). The phase relations remained largely insensitive to the variations, except that the assemblage in a drier system would be subsolidus (Fig. 8C). The garnet X_GRS isopleths were the most sensitive to the bulk H₂O and Fe³⁺/Fe_T, and the bulk composition used (black arrows in Figs. 8C and 8D) corresponded to the observed garnet X_GRS in the peak assemblage. In the P-H₂O pseudosection, the garnet Mg# and clinopyroxene j(o) isopleths are parallel. If the system had been hydrated during retrogression, or the peak assemblage had been drier, the garnet Mg# and integrated clinopyroxene compositions would lead to a lower pressure 1.8 GPa within the stability field of the peak assemblage. Similarly, if the system had been oxidized during retrogression, the peak pressure would be overestimated by 0.3 GPa at most. In summary, the general estimate of the peak pressure is not significantly affected by the uncertainties in the bulk H₂O and Fe³⁺/Fe_T ratio.

The P-T conditions (Kohn, 2020) corresponding to the high-Zr group of rutile (1539 ± 46 ppm; Fig. 6A) are plotted in Figure 8 (~840–870 °C, >1.5 GPa). The independent thermometer yielded slightly higher temperature or lower pressure than the garnet isopleth intersection but was largely consistent with the temperature of the peak pressure constrained from garnet and clinopyroxene compositional isopleths (Fig. 8A). If Zr were undersaturated, the temperature calculated with the Zr-in-rutile thermometer would be underestimated.

Additionally, we used a Ti-in-biotite geothermometer empirically calibrated by Wu and Chen (2015). Although incorporating Ti in biotite cannot be attributed to specific cation exchange reactions, Ti shows a clear affinity to biotite. As biotite decomposes at high temperatures, Ti is passively enriched in biotite. Thus, Ti content in biotite is strongly dependent on temperature (Henry and Guidotti, 2002; Henry et al., 2005). Sample XCG6 is suitable for the application of this geothermometer since it is TiO₂-saturated, as evidenced by the presence of rutile and ilmenite. Because the matrix biotite flakes overprint the peak-pressure assemblage and symplectite, they grew after or near the end of decompression. At 1.0–1.5 GPa, the Ti-in-biotite thermometer yields 800–900 °C. This condition agrees well with the temperature peak estimated by Cao et al. (2021) (C21 in Fig. 8A). Compared with the Zr-in-rutile geothermometer and phase equilibria analysis, the result of the Ti-in-biotite geothermometer is less robust. First, the Ti-in-biotite geothermometer by Wu and Chen (2015) was calibrated for metapelite. Applying the thermometer to a metabasite leads to unquantified uncertainties. Second, the resulting P-T condition marginally overlaps the calibrated P-T range (450–840 °C, 0.1–1.9 GPa), which arguably weakens the statistical significance of the result of the Ti-in-biotite thermometer.

We attempted to constrain the time scales of decompression and cooling using multicomponent diffusion simulation on the garnet profiles near the rims. We used the diffusion model provided in Chu and Ague (2015), which factors in temperature, pressure, unit cell dimension, and oxygen fugacity. We modeled along the decompression path from 2.0 GPa to 1.2 GPa at 850 °C, followed by cooling to 650 °C at 1.0 GPa, at a constant carbon-carbon oxide (CCO) oxygen fugacity (f_O2) buffer. This P-T path was inferred from phase equilibria analysis and Cao et al. (2021) (Fig. 8A), but neither segment is quantitatively constrained by thermobarometry. Given that the activation energies of these four diffusing species are 200–300 kJ/mol, an uncertainty in temperature of ±100 °C leads to a variation of the diffusion coefficient of ±1 log₁₀ unit (Chu and Ague, 2015), and accordingly ±1 order of magnitude in the time scales.

Each P-T segment was discretized to 20 isothermal steps, and diffusive relaxation was modeled at each successive step. The multicomponent diffusion coefficient matrix was constructed using the formalism of Lasaga (1979) under the assumption of ideal solid solution; such simplification has been shown to be valid for natural garnet compositions (Borinski et al., 2012). We utilized a Crank-Nicolson finite-difference scheme (Matlab code from Chu and Ague, 2015) to solve the diffusion equations numerically, and we included the spatial averaging effect by electron probe microanalysis beams (Borinski et al., 2012). We assumed the garnet grains had spherically symmetric compositions to the first order, and the modeled garnet profiles were sectioned through the center of symmetry.

For internal profiles of garnet grains XCG6B (Fig. 9) and XCG6D (Figs. S2 and S3), we were primarily interested in the diffusion modification on the X_GRS profile. Fe, Mg, and Mn were more homogenized, so the assumptions for initial profiles entailed large uncertainties, and the observed profiles could be fit using a series of initial profiles. This also illustrates the difficulty of constraining the prograde history recorded by the garnet core. Therefore, although multicomponent diffusion modeling was conducted on four components, the fitting of X_GRS was prioritized. We assumed the initial compositional zoning to be step functions. Profile 1 (Fig. 5D) crossed a thinner, high-Ca mantle. The highest X_GRS (“1” in Figs. 7 and 8) corresponded to the solidus condition, and the outer rim (“2” in Figs. 7 and 8) recorded growth en route to the peak pressure. Thus, we used a higher X_GRS (0.225) as the initial condition for the high-Ca mantle (“1” in Fig. 9A), with the assumption that the X_GRS of the narrower mantle (up to 0.175; Fig. 5H) had been modified by diffusion from both sides. Another exercise, without the extrapolation of X_GRS, was conducted for the garnet in XCG6D (Fig. S2F).

For the resorption profiles at the rims (Fig. 9A; Fig. S3), we used the diffusion profiles from the previous step (green-dotted line in Fig. 9A) or simplified flat profiles (Fig. S3) as the initial conditions. We used the Dirichlet boundary condition, assuming the garnet rim stayed in equilibrium locally with the matrix as diffusion proceeded. We note that diffusion rate increased exponentially with temperature, and so diffusion modification of the initial cooling (close to ~850 °C) was responsible for most interior smoothing, while the rim was continuously modified by moving boundary conditions. The extent of resorption was estimated by assuming Mn conservation (Spear, 2014), where all Mn components of the incrementally resorbed garnet rim were piled up at the garnet resorption front due to strong affinity of Mn to garnet. The moving boundary conditions for Fe, Mg, and Ca were linearly discretized for simplicity, from the initial boundary conditions to the apparent rim compositions of resorbed garnet (white arrows in Fig. 9A). For a conservative order-of-magnitude estimate, we assumed that rim resorption took place during cooling only, at a uniform cooling rate. A nonlinear or segmental cooling history could fit the profiles better (e.g., Chowdhury and Chakraborty, 2019; Zou et al., 2020).

Best-fit profiles and associated modeling durations were selected by visual inspection. We prioritized Ca profiles in the fitting because Ca diffuses the slowest in garnet. Compared with other diffusion models, like those of Carlson (2006) and Borinski et al. (2012), the program we used leads to slower diffusion of Ca in garnet. The modeling results are shown in Figure 9 and in Supplemental Material Figures S2F, S3F, and S4 (see ^{footnote 1}). The X_GRS and X_ALM profiles across the boundary between garnet core and thick mantle (profile 2; Figs. 5 and 8B) could be fitted perfectly using step-function initial conditions. If diffusion modification of profile 2 of XCG6B had taken place solely during decompression, the duration would be 80 k.y. (Fig. 9B). Profile 1 of garnet in sample XCG6B (Fig. 9A) yielded 62 k.y. of decompression followed by 88 k.y. of cooling. For Garnet_1 of XCG6D as a comparison (Fig. S2), we assumed all diffusion modification took place during decompression and no X_GRS extrapolation. The Ca profile (X_GRS) produced a 48 k.y. duration; if a higher-X_GRS mantle is assumed, the duration would be longer. The resorption profiles of Garnet_3 of XCG6D (Fig. S3) record a 240 k.y. cooling history from 850 °C to 650 °C. In sum, decompression from the peak pressure to a granulite-facies condition lasted <100 k.y., followed by cooling for ~100–200 k.y.

The inherent uncertainties from the diffusion model are ±0.5 orders of magnitude (Chu and Ague, 2015). The modeled temperature is close to or overlaps the temperatures at which diffusion experiments were conducted (summarized in Chu and Ague, 2015), so the potential error caused by down-temperature extrapolation from experimental data is minor. Additional uncertainties in the diffusion modeling results arise from uncertainties in the input P-T conditions and the validity of assumptions made about the initial and boundary conditions. The input P-T path has a significant impact on the diffusion time scales. In particular, a hotter temperature condition would result in a shorter duration to reproduce the observed profiles and vice versa. The P-T condition input was constrained from compositional isopleths in phase equilibria analysis, independent thermometers, and literature on other HP rocks associated with the western Grenville Allochthon Boundary thrust (Fig. 8A and summarized in the “Regional Geology” section). Conservative uncertainties associated with phase equilibria analysis are around the order of ±50 °C and ±0.1 GPa for strongly temperature- and pressure-dependent boundaries (Powell and Holland, 2008; Palin et al., 2016). Even if the temperature had been overestimated by 100 °C, the modeled diffusion time scale of decompression would be longer by about one order of magnitude but still brief (<1 m.y.).

The diffusion rate is proportional to f_O2^1/6; a higher oxygen fugacity by six orders of magnitude leads to one order-of-magnitude faster diffusion (Chakraborty and Ganguly, 1991). The f_O2 of eclogite-to-granulite-facies metabasite commonly varies between CCO + 1 and CCO + 2.5 (e.g., Tao et al., 2017). Thus, the diffusion rate is underestimated at CCO, and, accordingly, the diffusion time scale is overestimated.

We simulated diffusion profiles in a spherical model, and the relative distances to the core were estimated. A small shift in the radial coordinate causes negligible differences in the modeling results (Fig. S4 [see ^{footnote 1}]). Additionally, the garnet profiles used in the modeling were assumed to have been sectioned through the center of the garnet grain. Any off-center analysis and inclined sectioning would result in the direction of diffusion not being parallel to the profile, yielding longer time-scale estimates (Ague and Baxter, 2007; Chu and Ague, 2015).

We assumed sharp step functions as initial growth zoning developed in supersolidus conditions and investigated diffusion modification during the subsequent decompression. If the initial condition were gradational, the diffusion modeling results would yield shorter durations (Ague and Baxter, 2007). If the high-Ca mantle (“1” in Figs. 7A and 9A) had been higher, for example, X_GRS = 0.25, then diffusion to lower the central content of an ~0.01-mm-thick sheet to X_GRS = 0.225 would require an additional ~60 k.y. Thus, this different initial condition would not affect the order-of-magnitude estimate, and such overestimation does not affect the implication of transient metamorphism.

Synthesizing all the P-T information we gathered through phase equilibria modeling, independent geothermometers, and literature on rocks from similar geologic settings, we formulated the P-T evolution history of the retrograde eclogite from the Mattawa area.

The prograde path is poorly recorded, and it is beyond the scope of this study for a few reasons. First, the euhedral cores that recorded the prograde growth have been heavily modified by diffusion along the early prograde path and subsequently at or near the peak temperature. Mg# zonation would have been a valuable indicator of the increasing temperature, but the Mg# profile has been nearly completely smoothed except at the mantle and rim (Figs. 5H and 5I). Thus, the high-P conditions (“1”→“2” in Fig. 8A) are largely constrained by the Ca zonation, which diffuses the slowest in garnet. In addition, the phase equilibria analysis was conducted with an effective bulk composition that removes the euhedral garnet cores from the modeled system. Finally, the presence of melt and its potential removal pose additional uncertainty, as the preserved bulk composition is only a snapshot of the final stage.

The P-T path entered supersolidus conditions at ~690 °C and 1.8 GPa (“1” in Fig. 8A), where Ca- and Ti-rich garnet formed an overgrowth on the euhedral cores (Fig. 7A). These conditions (~680 °C, 1.8 GPa) are broadly consistent with the “pressure peak” revealed in Cao et al., 2021 (C21 in Fig. 8A) and other high-P metabasites close to Allochthon Boundary thrust (Marsh and Culshaw, 2014; MC14 in Fig. 8A). The subsequent heating to ~850 °C at 2.0 GPa had not been documented in previous studies of other early Ottawan high-pressure rocks associated with western Grenville Allochthon Boundary thrust. In this study, this peak pressure at a much higher temperature was supported by the garnet Ca and Mg# isopleths and the reintegrated omphacite Na/(Na + Ca) ratio (Fig. 8A). The Mg# peak of the garnet (“2” in Fig. 7), closer to the rim than the Ca peak, suggests prograde growth and distinguishes this compositional pattern from retrograde resorption. The garnet rim that records the peak-pressure conditions is susceptible to resorption during decompression. If retrograde resorption had been responsible for modifying the rim, then the Mg#, X_GRS, and X_SPS upturns should align on the measured profile (e.g., Fig. S3; fig. 7 in Marsh and Kelly, 2017; fig. 7 in Cao et al., 2021). In garnet grains with a higher degree of diffusion and resorption modification, more smeared Mg# and X_SPS profiles could potentially lead to alternative interpretations.

The eclogite-facies stage was followed by decompression to granulite facies. The formation of symplectite indicates a failure to reach equilibrium in a relatively anhydrous condition (Gaidies, 2021), and the phase relations should be viewed as qualitative. If garnet had been resorbed to a greater extent, it would not record the pressure peak at 2.1 GPa and ~850 °C. The specific garnet that we used to estimate the peak condition (Fig. 5) was only weakly resorbed, likely due to the sluggish reaction kinetics and a shrinkage of equilibration volume (Stüwe, 1997; Guevara and Caddick, 2016). The P-T condition at the end of the decompression path is recorded as the “temperature peak” in Cao et al. (2021) (Fig. 8A). We note that the predicted garnet composition at the “temperature peak” (1.2–1.3 GPa, 850–900 °C) differs from the pressure peak in its lower Mg# (~0.35 vs. ~0.27). Mg# is more susceptible to modification than Ca content, so the equilibrated garnet rim (e.g., Mg# = 0.27 at Grt_3 rim; Fig. S3) yields lower pressure and/or temperature. At another locality (Marsh and Kelly, 2017), the garnet rim compositions have been completely modified, so the equilibrium among the garnet rim and matrix phases cpx + hb + pl + q yielded a retrograde P-T condition of 1.2 GPa, ~860 °C (Fig. 8A). The high TiO₂ content in the retrograde biotite flakes broadly supports this P-T condition (Fig. 8A).

The peak pressure was followed by rapid (<100 k.y.) near-isothermal decompression (Fig. 8A). The exact decompression and cooling paths are not as robust as the peak-pressure condition, and they were simplified for the diffusion simulation. Diffusion is a cumulative effect, so a different or more complex retrograde history would not change the overall time scales of cooling significantly. The Zr contents in rutile suggest that the system was heated to ~850 °C at the peak pressure or during decompression. If decompression took place at lower temperatures, then the duration would have been longer by one order of magnitude only if the temperature was overestimated by ~100 °C; even in this case, the decompression time scale is still short (<1 m.y.). Decompression and/or cooling from the peak pressure eventually led to the crystallization of the partial melt and the retrograde assemblage observed in the field. The metamorphic peak condition was not revealed in previous studies, likely due to the modification of garnet by diffusion and resorption at rims.

Both outcrop and microscopic textures showed that a small amount of melt was present in the retrogressed eclogite, corroborated by the supersolidus assemblages from the modeled phase diagram (Fig. 8A). Microtextures indicating partial melting were described in great detail in Cao et al. (2021), including polymineralic inclusions near garnet rims and sparse felsic veinlets with no external connectivity. The lack of external connectivity in the felsic veinlets suggests that they did not originate from fluid flux. Utilizing the documented melt-related assemblages and modeled phase diagrams, Cao et al. (2021) proposed that phengite-dehydration melting with the involvement of amphibole and omphacite led to initial melting of the metabasite, and amphibole-dehydration melting involving garnet, omphacite, and biotite generated additional melt in the mafic system. Additionally, Feng et al. (2021) proposed that the omphacite breakdown reaction is also an important contributor to producing partial melt in eclogite through dehydroxylation (Wang et al., 2017, 2020). Omphacite breakdown is strongly dependent on pressure (Fig. 8A), so the rapid decompression likely contributed to the production of partial melt also through omphacite breakdown.

Cao et al. (2021) modeled the maximum modal abundance of melt to be ~10% at 850 °C. The critical porosity for melt segregation in the lower crust is estimated to be anywhere from 10% to 50%, and the lower crust can hold up to 30% interstitial melt without draining or eruption (e.g., Wickham, 1987; Yu and Lee, 2016). The value of 10% partial melt is close to the lower limit for melt segregation. The subsequent in situ recrystallization left behind little or no textural or mineralogical evidence for partial melting except inclusions in refractory phases. This small fraction of melt, however, could have significantly weakened the host rock and altered the rheological structure of the lower crust (Rosenberg and Handy, 2005). The presence of partial melt facilitates exhumation because it can increase buoyancy, lower the bulk viscosity, and enhance deformation; had the partial melt segregated, it would have lubricated the shear zones and facilitated exhumation (Labrousse et al., 2011). Although the minor amount of interstitial melt was not significant enough to form a composite melt transport system, the presence of melt still could have contributed to the fast decompression by assisting exhumation.

Assuming the thermodynamic pressure recorded by the mineral assemblage is lithostatic, the decompression rate converts to an exhumation rate of decimeters per year. The uncertainties associated with this estimation are mainly inherited from the pressure in the phase equilibria analysis and the time scale in the diffusion simulation. Using ±0.3 GPa uncertainty for the pressure drop and 0.5 orders of magnitude uncertainty for time scale (±50 °C), the lower bound of the range (0.08 m/yr) is still rapid. Some of the fastest exhumation rates previously reported range from 0.01 to 0.1 m/yr (e.g., Rubatto and Hermann, 2001; Hacker et al., 2003; Parrish et al., 2006). Most of these rates were estimated based on isotopic dating techniques. In comparison, exhumation rates estimated from diffusion speedometry are typically higher because rates calculated with isotopic methods are limited by the time resolution of isotope geochronology (Viete and Lister, 2017).

The reconstructed P-T path suggests that cooling occurred after the rock was exhumed to midcrustal levels (35–40 km). If there were significant residence time at high temperature following the decompression path, the garnet compositional profiles would have been further smeared out. Diffusion rates at temperatures <650 °C are slow enough that the Mattawa eclogite could have remained moderately hot (e.g., ~600 °C) over an extended period while minor diffusion modification took place, especially if the rock was not rehydrated (Zhang et al., 2019). Fast cooling and sluggish reaction kinetics at low temperatures also explain the relatively pristine anhydrous assemblage and minor amphibolite-facies overprinting. The short residence time at granulite-facies condition contrasts with the protracted high-temperature conditions inferred for the long-lived Ottawan orogeny (McLelland et al., 2010; Hynes and Rivers, 2010; Volkert and Rivers, 2019; Indares, 2020). For example, monazite records >70 m.y. of supersolidus conditions in the Grenville hinterland (Turlin et al., 2018). U-Pb depth profiling of rutiles in eclogite–HP granulite from the Shawanaga domain indicates an extended period (60 m.y.) of high-temperature metamorphism from 850 °C to 750 °C, and subsequently even slower cooling (Smye et al., 2018).

On the other hand, Grenvillian crusts also feature unique high-temperature signatures (Moecher et al., 2014), anorthosite massifs (Corrigan and Hanmer, 1997), and massive granulite-facies terranes (Indares, 2020, and references therein). Protracted heating with the lower crust and heat upwelling from the underlying mantle would significantly destabilize the orogenic plateau by gravitational collapse in the weak, hot, overthickened crust (Rivers, 2012), thermal erosion of the lower crust (Corrigan and Hanmer, 1997), and delamination (Lieu and Stern, 2019). Indeed, repeated thinning events punctuating the long-lived Ottawan phase are recorded by bulk-rock and zircon proxies for crustal thickness (Brudner et al., 2022). Crustal-scale syn- and postcollisional extension is well documented in the western (e.g., Schwerdtner et al., 2014, 2016) and central (e.g., Corrigan and Breemen, 1997; Soucy La Roche et al., 2015) Grenville Province, and in the Adirondacks (e.g., Wong et al., 2012; Baird, 2020). In the early Ottawan phase particularly, the activity of Midcontinent rifting (1110–1085 Ma; Swanson-Hysell et al., 2019) overlapped the metamorphic ages of eclogites along the Allochthon Boundary thrust (ca. 1090 Ma; Marsh and Culshaw, 2014; Marsh and Smye, 2017). The transient and local extension was not restricted to the Midcontinent Rift, as contemporaneous intraplate magmatism was widespread across Laurentia (e.g., Guitreau et al., 2016; McLelland et al., 2010; Bright et al., 2014). A small rift basin in Arctic Laurentia yielded a 1087.1 ± 5.9 Re-Os depositional age, and series of such isolated rift basins of similar ages demonstrate that localized extension at the margins of Laurentia was pervasive (Greenman et al., 2021, and references therein).

The fast exhumation and cooling revealed by the Mattawa eclogite might correspond locally to extension or collapse events shortly after or coeval with the culmination of burial. In response to extension in the upper- to midcrustal sections, the spread of partially molten lower-crustal material to the foreland thickens the orogen edge and results in eclogite-facies metamorphism at high temperatures; this is closely followed by upward migration of hot, low-viscosity migmatites that host metabasites (Whitney et al., 2015; Rey et al., 2017). As soon as rifting takes place in the upper crust, the gap created accommodates rapidly ascending, hot lower-crustal material. Deep crustal rocks beneath a localized region of extension first experience near-isothermal decompression, and then the still-hot rocks are emplaced into cold mid- to upper crust, which helps to rapidly cool the rocks (Rey et al., 2009). In western Grenville, such midcrustal metamorphic core complexes are found juxtaposed against the upper crust (Rivers, 2012). Given that metabasite dikes or pods in the Mattawa domain are ~500 m wide (Fig. 1C), the time scale of conductive cooling would have been ~10 k.y. The activity of the Allochthon Boundary thrust and reworked extensional detachment faulting might have facilitated subsequent exhumation at a slower rate (Rivers et al., 2002).

This study adds to the few records of Precambrian short-duration metamorphism (e.g., Guevara et al., 2017; Zou et al., 2020). The previous lack of such records left the question open-ended: Was short-duration metamorphism genuinely absent in the Precambrian, or did the lack of evidence result from other biases (Viete and Lister, 2017)? The recent discoveries of Precambrian fast metamorphism suggest that additional occurrences could be expected, and diffusion speedometry is a suitable tool with which to identify them because its resolution is independent of age. The emerging new records of fast metamorphism shed light on the transition of tectonic regimes to modern plate tectonics (Chowdhury et al., 2021) and add another dimension to petrotectonic processes in the well-studied Precambrian Grenville orogen.

The retrogressed eclogite near Mattawa, Ontario, Canada, records a P-T path of supersolidus isobaric heating to eclogite-facies conditions followed by decompression and cooling. Garnet diffusion speedometry constrains the durations of decompression to be <100 k.y. and cooling to within 100–240 k.y. The short time scales imply transient metamorphism within a long-lived hot orogen. The fast decompression rate estimated from the short duration indicates rapid exhumation to <40 km, likely in response to syn- or postcollisional crustal thinning in the early Ottawan phase. The discovery of Precambrian short-duration metamorphism is a potential indicator for operative modern-style plate tectonics in the late Mesoproterozoic.

This project was supported by the Undergraduate Student Research Award (USRA, award 551809-2020 to X. Fan) and Discovery Grant (RGPIN-2018-03925 to X. Chu) provided by the Natural Sciences and Engineering Research Council of Canada (NSERC) and the National Natural Science Foundation of China (grant 42002211 to Y. Zou). We thank Yanan Liu at the University of Toronto for her help with electron microprobe analyses, and Adam Brudner and Junxing Chen for their discussion, comments, and proofreading. We sincerely thank Renaud Soucy La Roche, Fred Gaidies, and Aphrodite Indares for their constructive and careful reviews, and Nancy Riggs and Dawn Kellett for editorial handling and comments.