The Pre-Grenvillian assembly of the southeastern Laurentian margin through the U–Pb–Hf detrital zircon record of Mesoproterozoic supracrustal sequences (Central Grenville Province, Quebec, Canada)

Abstract

The detrital zircon perspective on the pre-collisional crustal evolution of the Grenville Province remains poorly explored. In this study, we conducted in situ laser ablation U–Pb–Hf isotopic microanalysis on detrital zircon grains from three pre-orogenic (>1 Ga) supracrustal sequences that crop out in the Central Grenville Province (Lac Saint-Jean region, QC, CA). Detrital zircon grains from vestiges of these sequences record three dominant age peaks at c. 1.46 Ga, 1.62 Ga, 1.85 Ga, and a subordinate peak at 2.7 Ga. The 1.46 Ga and 1.62 Ga age peaks are recorded in detrital zircon grains from a quartzite associated with a metavolcanic sequence (i.e. Montauban Group) with a maximum depositional age of c. 1.44 Ga. In contrast, the c.1.85 Ga age peak is observed from recycled zircon grains in metasediments with maximum depositional ages between 1.2 and 1.3 Ga. The suprachondritic Hf isotope composition in detrital zircon grains of the 1.46 Ga and 1.62 Ga age populations records juvenile crustal growth during peri-Laurentian accretionary orogenesis related to the Pinwarian (1.4–1.5 Ga) and Mazatzalian–Labradorian (1.6–1.7 Ga) events. The detrital zircon grains associated with Penokean–Makkovikian (1.8–1.9 Ga) source rocks record reworking of c. 2.7 Ga continental crust derived from a near-chondritic mantle reservoir. Overall, crust-forming and basement reworking events associated with accretionary orogenesis in southeastern Laurentia are retained in the detrital zircon load of Precambrian basins even after the terminal Grenvillian collision and assembly of Rodinia.

1. Introduction

The Grenville orogen (c. 1090–980 Ma; Rivers, 2008, 2021) is a global-scale collisional orogenic system that terminated the long-lived accretionary growth of Laurentia, leading to the formation of the supercontinent Rodinia (Hoffman, 1988; Whitmeyer and Karlstrom, 2007). In the southeastern Canadian shield (Grenville Province; Fig. 1a), c. 1.15 Ga (Shawinigan) and 1.05 Ga (Ottawan) high-grade metamorphism and penetrative deformation concealed the traits of accretionary tectonism, leaving behind gneiss terranes of ambiguous origin. Regional-scale Sm–Nd isotope mapping along the Grenville Province has attempted to untangle the pre-collisional crustal formation history of these large tracts of gneisses (Dickin & McNutt, 1989; Dickin, 2000; Dickin et al., 2009). However, the perspective of detrital zircon on the pre-collisional crustal growth of the Grenville Province remains underexplored despite zircon grains of Laurentian affinity being documented in pre-, syn- and post-collisional metasedimentary sequences (Sager-Kinsman and Parrish, 1993; Van Breemen and Corriveau, 2005; David et al. 2010; Kamo et al. 2011; Lasalle et al. 2013; Spencer et al. 2015; Groulier et al. 2018; Peck et al. 2019). Together, these studies have documented, through detrital zircon U–Pb geochronology, the presence of 1.4–1.9 Ga sources but, further insights about their crustal production and evolution are missing.

Fig. 1. Simplified geological map of the Lac-Saint-Jean region. The inset (a) depicts the crustal terranes of the Grenville Province sensu Dickin, 2000). The white dashed rectangle depicts the study area. The red circles in (b) denote the four sampling sites.

In this regard, the interrogation of the detrital zircon record, from the Hadean to Phanerozoic, using U–Pb–Hf ± O isotope systematics has placed critical new constraints on the secular evolution and growth of the continental crust (Stevenson & Patchett, 1990; Kemp et al. 2006; Belousova et al. 2010; Hawkesworth et al. 2010; Parman, 2015). Crucially, the initial Hf isotopic composition of single zircon grains monitors crust–mantle differentiation events while remaining resistant to high-grade metamorphism and/or hydrothermal fluid flow that can compromise the whole-rock Sm–Nd isotope system (Moorbath et al. 1997; Vervoort & Blichert-Toft, 1999; Hoffmann et al., 2011). The advantage of the U–Pb and Lu–Hf isotopic couplet has been demonstrated in many crustal evolution studies (e.g. Condie et al., 2005; Kemp et al. 2015; Fisher & Vervoort, 2018) but remains underexploited in the metasedimentary sequences of the Grenville Province. The scarcity of zircon U–Pb–Hf isotopic data from these sequences has inhibited our ability to ground-truth crustal evolution models, which are heavily based on whole-rock Sm–Nd isotopic data from felsic orthogneisses, and in some cases interpreted without U–Pb age constraints from the same samples (e.g. Dickin, 2000).

In this contribution we combine, for the first time, detrital zircon U–Pb (n = 623) and Lu–Hf (n = 366) isotopic data from rare outcrops of quartzites in three Mesoproterozoic supracrustal sequences, from the Lac Saint-Jean region in the central part of the Grenville Province (Fig. 1b). The main aims of our study are to: (a) provide new constraints on the pre-collisional geodynamic evolution and crustal growth of the Grenville Province through examination of the U–Pb–Hf detrital zircon record, (b) compare the zircon U–Pb–Hf isotope signatures with the existing whole-rock Nd isotopic record of crustal growth in the central Grenville Province and (c) better constrain the depositional ages and sources of these poorly documented supracrustal sequences.

2. Pre-collisional accretionary terranes of the Grenville Province

The Grenville Province is divided into two major lithotectonic belts, the Parautochthonous and the Allochthonous Belts (Rivers, 1997). The Parautochthonous consists of basement rocks with Laurentian isotopic and age characteristics, whereas the Allochthonous Belt comprises translated nappes of amphibolite to granulite-facies gneisses (Rivers, 2008). Biotite gneisses in the western part of the Parautochthonous belt exhibit depleted-mantle Nd model ages from 2.5 to 2.9 Ga (Guo & Dickin, 1996) although Archaean zircon U–Pb crystallization ages have rarely been reported from these areas (Krogh, 1994; Jordan et al. 2006).

The Allochthonous Belt has been characterized, through Nd isotope mapping of granitoid orthogneisses, as an amalgamation of Palaeo- to Mesoproterozoic accretionary terranes (Dickin, 2000). These terranes, sensu Dickin, (2000), are the Makkovikia–Barilia, Labradoria and Quebecia terranes (Fig. 1a) with crustal formation ages (i.e. depleted-mantle model ages) of c. 1.5 Ga, 1.7–1.8 Ga and 1.9 Ga respectively. In more detail, the Makkovikia has been characterized as a terrane of reworked crust located at the easternmost part of the Grenville Province (i.e. Labrador), and the Barilia as a 1.8–2.0 Ga juvenile Palaeoproterozoic arc in Ontario (Dickin et al., 2009). These terranes are regarded as products of an arc–continent collisional tectonothermal event known as the Penokean–Makkovikian orogeny, which is well documented in the Great Lakes region and took place between c. 1.8 and 1.9 Ga (Van Schmus & Bickford, 1981; Barovich et al. 1989). U–Pb age-peaks in detrital zircon spectra at c. 1.8–1.9 Ga, from different metasedimentary sequences, have been reported along the Grenville Province and are commonly ascribed to Penokean–Makkovikian sources (Kamo et al. 2011; Moukhsil et al. 2013; Culshaw et al. 2016; Peck et al. 2019).

Labradoria was formed during the Mazatzalian–Labradorian orogeny between 1.6 and 1.7 Ga (Wasteneys et al. 1997; Papapavlou et al. 2017) and crops out mostly in Labrador, in the eastern part of the Grenville Province (Connelly & Heaman, 1993; Gower, 2019). Labradorian-age metaplutonic and mafic/ultramafic suites have also been documented in the central and western parts of the province (Dunning & Indares, 2010, Culshaw et al. 2016). In eastern Labrador, c. 1.7 Ga tonalitic and granitic rocks exhibit near-chondritic isotope signatures (ϵNd = −0.13 to +0.42; Scharer, 1991). Later studies suggested that the Labradorian crust-forming event was mostly juvenile but also involved magmatic reworking (i.e. remelting) of older Makkovikian crust in the eastern part of the province (Dickin, 2000). As with detrital zircon grains of Penokean age, grains of Labradorian age are ubiquitous in the pre- and syn-Grenvillian metasedimentary record (Corrigan & Van Breemen, 1997; Lasalle et al. 2013; Culshaw et al. 2016; Peck et al. 2019).

The youngest accretionary event in the pre-collisional history of the Grenville Province is the c. 1.46–1.52 Ga Pinwarian orogeny (Tucker & Gower, 1994; Gower & Krogh, 2002) that is expressed through titanite U–Pb dating in reactivated mylonitic shear zones in adjacent orogens (e.g. Southern Province; Papapavlou et al. 2018). In the eastern part of the Grenville Province, Pinwarian magmatism reworked Labradoria whereas further west it led to the construction of Quebecia (Dickin, 2000). Recent studies focusing on the Pinwarian magmatism suggest that it took place above two subduction zones with opposing polarities: a north-dipping subduction zone below Laurentia and a south-dipping intra-oceanic subduction zone, which led to the construction of a nascent magmatic arc (Martin & Dickin, 2005; Vautour & Dickin, 2019).

The Mid- to Late-Proterozoic Elzevirian (c. 1.2–1.3 Ga) and Shawinigan (c. 1.19–1.14 Ga) orogenic pulses constitute the final accretionary events (Rivers, 2008) before the collisional Grenville orogeny that lasted until c. 960 Ma (Jannin et al. 2018, Turlin et al. 2019) and was ended by late collisional alkaline magmatism (Groulier et al. 2020b).

3. Regional geology of the Saguenay – Lac Saint-Jean region

The Lac Saint-Jean region is located in the central part of the Grenville Province within the Allochthonous Belt (Fig. 1). The oldest lithostratigraphic unit in the area is the Pinwarian Cap de la Mer amphibolite, with a U–Pb zircon age of 1506 ± 13 Ma, which is included within a group of basement units called the Saguenay Gneiss Complex (Hébert & Van Breemen, 2004). This area belongs to the Portneuf – St Maurice lithotectonic domain that consists of two main units, the La Bostonnais complex, which is a 1.4–1.37 Ga suite of metamorphosed calc-alkaline plutonic rocks, and the Montauban group, a c. 1.45 Ga metamorphosed volcano-sedimentary sequence commonly observed as enclaves within the La Bostonnais complex (Fig. 1b; L Nadeau et al., unpub. data, 1992; Corrigan & Van Breemen, 1997). The area is dominated by anorthosite–mangerite–charnockite–granite (AMCG) magmatic suites and hosts the >20 000 km² Lac Saint-Jean Anorthositic Suite (Woussen et al. 1981). U–Pb zircon dating (isotope dilution thermal ionization mass spectrometry (ID-TIMS)) has shown that the Lac Saint-Jean anorthosite suite was emplaced between c. 1160 and 1140 Ma and was followed by two less voluminous AMCG magmatic pulses at c. 1082–1050 Ma and 1020–1010 Ma (Higgins & Van Breemen, 1996). Regional Grenvillian metamorphism in the area attained granulite facies with P–T conditions of 6.2 kbar and 780 °C followed by a retrograde amphibolite facies event reaching 4.7 kbar and 690 °C , respectively (Hébert & Van Breemen, 2004).

3.a. Supracrustal lithological complexes

In this study we explore the detrital zircon record from three of the four supracrustal complexes that crop out in the wider area of the Lac Saint-Jean region. These are the Montauban, Wabash, Barrois and Saint-Onge complexes (Fig. 1b). The Montauban complex consists of quartzofeldspathic gneisses interlayered with amphibolite, pillowed tholeiitic metabasalts with boudinage structure, volcanic and volcaniclastic rocks (i.e. lappili tuff), metalliferous calc-silicate rocks, garnetite, and minor quartzite which was sampled in this study (sample 17GC1086; Stamatelopoulou-Seymour & Maclean, 1977; Rondot, 1978). Supracrustal lithologies of the Montauban group are of economic importance since they host metamorphosed exhalative massive sulphide (Zn–Pb–Cu–Au–Ag) and stratiform gold occurrences (Bernier et al. 1987). The Wabash group is a metamorphosed volcano-sedimentary sequence that consists of four tectonostratigraphic units that mainly comprise migmatitic paragneisses, calcitic and dolomitic marbles, minor impure quartzite layers (sample 18YD2075), amphibolites and quartzofeldspathic rocks. The Barrois group comprises chiefly Opx-bearing orthogneisses, Bt-bearing paragneisses, charnockitic to enderbitic gneisses, amphibolite and garnetite. The Barrois group is divided into the Barrois I unit that contains metasedimentary lithologies (i.e. paragneisses) and the Barrois II unit that consists of orthogneisses which likely constitute the basement onto which the Barrois I unit was deposited. The biotite-bearing paragneisses consists of interbedded, decimetre- to metre-scale, quartzite horizons, which were sampled here (samples 19AM02 and 19FS7085). Aeromagnetic maps suggest the juxtaposition of the Barrois and Wabash supracrustal packages most likely along a stratigraphic contact. The Saint-Onge group is composed of quartzofeldspathic paragneisses, calc-silicate rocks, marbles and minor quartzite. Lithological members of this group show extensive alteration in outcrop scale and were not sampled in this study. Four quartzite samples from the Lac Saint-Jean region were collected: these samples are the 17GC1086 (Montauban Group), 18YD2075 (Wabash Complex), 19FS7085 and 19AM02 (Barrois Complex). Sample coordinates and petrographic descriptions are provided in supplementary material file 1 (in the Supplementary Material available online at https://doi.org/10.1017/S0016756821001023).

4. Analytical techniques and methodology

Sample processing and microanalytical work were conducted at the facilities of GEOTOP-UQAM (Montreal, QC, CA) using conventional zircon separation procedures involving crushing with a jaw crusher and disc mill followed by Wilfley table, heavy liquid, and magnetic separation using a Franz isodynamic separator. The zircon grains were annealed at 1000 °C for 48 hours and were then mounted in epoxy for imaging and microanalysis. Cathodoluminescence (CL) imaging of the zircon grains was conducted using a Hitachi S-3400N scanning electron microscope (SEM). The U–Pb and Lu–Hf isotopic analyses were performed during separate sessions on overlapping or adjacent analytical spots within the same CL growth domain assuming no differences in the isotopic composition or age of the different volumes. The U–Pb and Lu–Hf analytical sessions were performed using a Photon Machines Analyte G2 excimer laser (193 nm) coupled with either a single-collector high-resolution inductively coupled plasma mass spectrometer (Nu Attom, HR-ICP-MS, for U–Pb microanalysis) or a multicollector ICP-MS (Nu Plasma II, for Hf isotope microanalysis). The ablation took place in a two-volume HELEX II sample cell with He as the carrier gas which was mixed with the Ar make-up gas downstream before entering the ICP torch. Instrument parameters of the analytical sessions and secondary reference material data are provided in supplementary material file 1 (in the supplementary material available online at https://doi.org/10.1017/S0016756821001023). The U–Pb and Lu–Hf time-resolved isotopic data were reduced and processed using iolite 3.63 and IsoplotR (Paton et al. 2011; Vermeesch, 2018).

4.a. U–Pb isotopic analysis

The U–Pb isotopic data were collected using a spot size of 25 μm at a repetition rate of 5 Hz and fluence of 3 J cm⁻². The 91500 zircon (²⁰⁶Pb/²³⁸U age of 1062.4 ± 0.4 Ma; Wiedenbeck et al., 1995) was used as the primary reference material for the correction of laser-induced elemental fractionation and instrument drift. Secondary zircon reference materials were interspersed every 10–15 unknowns and yielded the following weighted mean ²⁰⁷Pb/²⁰⁶Pb ages: BB9 (567 ± 16 Ma, 2s, MSWD = 2, n = 8), Mudtank (713 ± 17 Ma, 2s, MSWD = 1.5, n = 35), GJ-1 (604 ± 11 Ma, 2s, MSWD = 2.7, n = 38), UQZ8 (1119 ± 9 Ma, 2s, MSWD = 5, n = 15) and Pentecote (1367 ± 39 Ma, 2s, MSWD = 3, n = 10). These dates corroborate within uncertainty the published reference values (Martignole et al. 1993; Machado & Gauthier, 1996; Jackson et al. 2004; Santos et al. 2017; Gain et al. 2019).

The U–Pb isotopic data (see also U–Pb supplementary material available online at https://doi.org/10.1017/S0016756821001023) that are presented from the cores of zircon grains have ≤10% discordance and Th/U > 0.1 following the rationale of previous studies (Papapavlou et al. 2021). The filtered data are presented as ²⁰⁷Pb/²⁰⁶Pb ages in Kernel Density Estimation (KDE) diagrams with an adaptive bandwidth using the Java application DensityPlotter (Vermeesch, 2012). Maximum depositional ages (MDAs) were calculated using the maximum likelihood age (MLA) estimation algorithm of Vermeesch (2021). Grains that do not show igneous zoning in CL were not considered during the calculation of MDAs but are presented in the U–Pb supplementary material (available online at https://doi.org/10.1017/S0016756821001023).

4.b. Lu–Hf isotopic analysis

The Lu–Hf isotopic data were collected using a 50 μm spot with a laser frequency of 15 Hz and fluence of c. 9 Jcm⁻². A typical Lu–Hf analysis consisted of 25 s of background measurement and 20 s of data acquisition. The ¹⁷¹Yb, ¹⁷³Yb, ¹⁷⁴Yb, ¹⁷⁵Lu, ¹⁷⁶(Yb + Lu + Hf), ¹⁷⁷Hf, ¹⁷⁸Hf, ¹⁷⁹Hf, ¹⁸⁰Hf, ¹⁸¹Ta and ¹⁸²W beams were collected using 11 Faraday cups in fixed positions attached to 10¹¹ Ω feedback resistors, with the mass separation being controlled by an optical zoom lens between the magnetic sector and the detectors. The 91500 zircon was used as the primary reference material to correct instrument drift during the analytical session. Mass bias and isobaric interference corrections of ¹⁷⁶Yb and ¹⁷⁶Lu on ¹⁷⁶Hf were made following the rationale and guidelines of previous studies (Woodhead et al., 2004; Fisher et al. 2014). Specifically, the β _Yb (Yb mass bias) was determined by measuring the ¹⁷³Yb/¹⁷¹Yb ratio relative to ¹⁷³Yb/¹⁷¹Yb = 1.132685 (Chu et al., 2002). Following the calculation of β _Yb, the contribution of ¹⁷⁶Yb on ¹⁷⁶Hf was subtracted using the intensity of the interference-free ¹⁷³Yb and a ¹⁷⁶Yb/¹⁷³Yb = 0.79618 (Chu et al., 2002). The contribution of the ¹⁷⁶Lu on the ¹⁷⁶Hf signal was estimated assuming that β _Lu = β _Yb using the intensity of the interference-free ¹⁷⁵Lu with a ¹⁷⁶Lu/¹⁷⁵Lu = 0.02656 (Blichert-Toft et al. 1997). The mass bias correction on the ¹⁷⁶Hf/¹⁷⁷Hf ratios was estimated using ¹⁷⁹Hf/¹⁷⁷Hf = 0.7325 (Patchett & Tatsumoto, 1981). Mass bias corrections for Yb, Hf and Lu were made using an exponential law. For the Yb mass bias corrections of grains with low Yb/Hf ratios (e.g. Mudtank and Munzirc02a, which have ¹⁷³Yb beam intensities below 5mV), the β _Hf was used instead of the poorly characterized β _Yb. For quality control on the mass bias corrections, the ¹⁷⁸Hf/¹⁷⁷Hf ratio on the analysed reference materials yielded a value of 1.46731 ± 52 (1s), in agreement with the published value of 1.46735 (Thirlwall & Anczkwiewicz, 2004).

Secondary zircon reference materials with variable ¹⁷⁶Yb/¹⁷⁷Hf ratios were analysed during the Lu–Hf session to ensure that our ¹⁷⁶Yb correction yielded accurate results (see also Lu–Hf supplementary material file, available online at https://doi.org/10.1017/S0016756821001023). The following weighted mean ¹⁷⁶Hf/¹⁷⁷Hf values were obtained: BB9 zircon (0.281666 ± 26, MSWD = 1.5, n = 9, 2s, ϵHf value of −27.2 ± 1.2), Maniitsoq (0.280865 ± 15; MSWD = 0.9, n = 10, 2s, ϵHf value of −0.28 ± 0.9), Munzirc12c (0.282129 ± 12; MSWD = 1.3, n = 17, 2s), Munzirc02a (0.282160 ± 8, MSWD = 1.5, n = 18, 2s), Munzirc32c (0.282150 ± 19; MSWD = 2, n = 14, 2s), Munzirc42b (0.282164 ± 14; MSWD = 0.3, n = 11, 2s), Plešovice (0.282475 ± 6; MSWD = 0.6, n = 30, 2s, ϵHf value of −3.5 ± 0.6), SA01 (0.282296 ± 24; MSWD = 1.5, n = 11, 2s, ϵHf value of −5.7 ± 1.5), GJ-1 (0.282006 ± 16; MSWD = 0.6, n = 17, 2s, ϵHf value of −14.3 ± 1.1) and Mudtank (0.282528 ± 11; MSWD = 1.5, n = 16, 2s, ϵHf value of +7.1 ± 1) (Morel et al. 2008; Sláma et al., 2008; Fisher et al. 2011; Santos et al. 2017; Gain et al. 2019; Marsh et al. 2019; Huang et al., 2020). The reproducibility (2s) of the secondary reference materials varies from 0.004 to 0.029 %. The analytical errors are quadratic additions of the in-run error (2SE) and the reproducibility (2s) of the Plešovice standard and are typically ∼± 2 ϵHf units. The initial ¹⁷⁶Hf/¹⁷⁷Hf isotope ratios and ϵHf values were calculated using the ¹⁷⁶Lu decay constant of 1.867 × 10⁻¹¹ (Söderlund et al. 2004) and the Lu–Hf Chondritic Uniform Reservoir (CHUR) parameters of Bouvier et al. (2008). The initial ¹⁷⁶Hf/¹⁷⁷Hf isotope ratios and ϵHf values for each detrital zircon were calculated using the ²⁰⁷Pb/²⁰⁶Pb date of each analytical spot.

5. Results

5.a. Sample 17GC1086 (Montauban Group)

The examined zircon grains are rounded to sub-rounded with aspect ratios 1:1 to 4:1, exhibiting oscillatory, sector-zoned or homogeneous zones in CL with brighter cores commonly overgrown by thin (<15 μm) darker rims (Fig. 2; spot 1086R7.3). Planar banding and domains with faint variations in CL response are also observed. Four dominant age peaks at 1.46 Ga, 1.62 Ga, 1.87 Ga and 2.7 Ga are recorded in this sample (Fig. 3; sample 17GC1086). The dominant age components at c. 1.46 Ga and 1.62 Ga indicate the prevalence of sources with Pinwarian (34 % of the analyses) and Labradorian (29 % of the analyses) affinities, with subordinate amounts of grains with Penokean ages (22 % of the analyses). The MLA algorithm yields a date of 1442 ± 8 Ma (Fig. 4a) that is interpreted as the MDA of the sample. This date agrees with the youngest cluster of nine U–Pb analyses from grains with oscillatory zoning in CL and Th/U values from 0.3 up to 1.1 that yield a weighted mean ²⁰⁷Pb/²⁰⁶Pb date of 1436 ± 11 Ma (n = 9, MSWD = 0.9, 2s). A cluster of younger ²⁰⁷Pb/²⁰⁶Pb dates from 1399 Ma to 1334 Ma is also observed, but the zoning patterns in CL of the grains are not diagnostic of igneous growth and were not considered in the estimation of the MDA (see also Supplementary Material, available online at https://doi.org/10.1017/S0016756821001023).

Fig. 2. Representative cathodoluminescence images from detrital zircon grains of the four quartzite samples. Red and blue circles show the location of the 25 μm and 50 μm beam location for U–Pb and Lu–Hf microanalysis, respectively.

Fig. 3. Kernel density estimation (KDE) ²⁰⁷Pb/²⁰⁶Pb age diagrams and weighted mean ²⁰⁷Pb/²⁰⁶Pb age diagrams from the youngest zircon grains in the four quartzite samples. The numbers adjacent to the peaks denote the age maxima in Ga. The weighted mean plots next to the KDE diagrams represent MDA estimates based on the youngest acceptable dates for each supracrustal complex.

Fig. 4. Radial plots and MDA estimates of the (a) Montauban Group quartzite, (b) Wabash Complex quartzite and (c, d) Barrois Complex quartzites using isoplotR (Vermeesch, 2018) and the Maximum Likelihood Age (MLA) algorithm (Vermeesch, 2021). The radial plots display heteroscedastic data (i.e. data points with unequal uncertainties) where each data point represents a ²⁰⁷Pb/²⁰⁶Pb date. These plots are scatter plots of standardized estimates (y axis) against precisions (x axis).The angular position of each data point denotes the date, whereas the horizontal distance from the point of origin varies with the precision of the analysis (i.e. imprecise measurements on the left side and more precise on the right side of the plot). For further details on age calculation see Vermeesch (2021) and Galbraith (1988).

The ϵHf values of the analysed zircon grains vary from −18.4 to +10.7 (Fig. 5). The majority of zircon grains with ²⁰⁷Pb/²⁰⁶Pb ages between c.1.4 and 1.56 Ga (92 % of the analyses) exhibit suprachrondritic ϵHf values (ϵHf > 0) consistent with the juvenile signature of Pinwarian sources. Similarly, 84 % of the analyses from the 1.6–1.7 Ga age population also yield suprachondritic ϵHf values (ϵHf = +0.2 to +5.8) with 16 % recording ϵHf values between −1.1 and −8. The third major age component, with ages between c. 1.8 and 1.9 Ga, shows an increasing contribution from sources that define a vertical mixing array between sub- and suprachondritic components (ϵHf between −1.2 and +7.0). The minor 2.7 Ga age component consists chiefly of grains with ϵHf values that vary from −10.3 to +3.8.

Fig. 5. Epsilon Hf versus age with Lu–Hf isotopic data from the detrital zircon grains of the four quartzites reported in this study. Black and white coloured are Lu–Hf data points from the dataset of Spencer et al. (2015). The black-filled data points represent data from the pre-orogenic Siamarnekh Formation (SF). The white-filled data points depict, for comparative purposes, Lu–Hf isotopic data from pre-, syn- and post-orogenic formations from Labrador and Scotland. The lines with different slopes (DM = depleted mantle, MC = mafic crust, CC = continental crust) correspond to different Lu–Hf isotope evolution curves for the depleted mantle (¹⁷⁶Lu/¹⁷⁷Hf = 0.03976; Vervoort et al. 2018), continental crust (¹⁷⁶Lu/¹⁷⁷Hf = 0.0113; Rudnick & Gao, 2003) and mafic crust (¹⁷⁶Lu/¹⁷⁷Hf = 0.022; Nebel et al. 2007) reservoirs. The CHUR reference line (ϵHf = 0) involves an uncertainty of ±0.4 ϵHf units (Bouvier et al. 2008). The grey band is drawn parallel to the ¹⁷⁶Lu/¹⁷⁷Hf = 0.022 line and covers c. 70 % of the 1.8–1.9 Ga zircon grains.

5.b. Sample 18YD2075 (Wabash complex)

Two grain-size populations consisting of euhedral zircon grains with aspect ratios varying from 1:1 up to 3:1 are observed in this sample. The examined grains exhibit bright cores in CL, with oscillatory or sector zoning, overgrown by darker euhedral rims of possible metamorphic origin (Fig. 2). Zircon grains with homogeneous grey to dark grey response in CL are also common. The region between the brighter rounded to sub-rounded cores and the darker rims often has a brighter response in CL, indicating a separate growth event (Fig. 2; 2075R1.3). The four age components recorded in the detrital spectrum are c. 1.61 Ga, 1.80 Ga, 1.90 Ga and 2.68 Ga (Fig. 3; 18YD2075). The MLA algorithm yields a date of 1309 ± 38 Ma (Fig. 4b) that is interpreted as the MDA of the sample, corroborating previous MDA constraints of this complex (Moukhsil et al. 2015). In comparison, the two youngest dates with Th/U values of 0.14 and 0.19 agree with the MDA yielding a weighted average ²⁰⁷Pb/²⁰⁶Pb date of 1316 ± 150 (2s).

The ϵHf values from the analysed zircon grains vary from −26.9 to +10.5 (Fig. 5). Specifically, the Lu–Hf analyses from the minor 1.4–1.5 Ga and c.1.6–1.7 Ga age populations have ϵHf values varying from +1.2 to +10.5 and −12.4 to +8.1, respectively. The majority of the analysed grains belong to the 1.7–2.0 Ga age range (56% of the analyses) and exhibit a spread in ϵHf values from −26.9 to +7.2. Zircon grains with ages >2 Ga (30% of the analyses) have ϵHf values that vary from −10.6 to +7.5.

5.c. Sample 19AM02B (Barrois complex)

The zircon grains of the 19AM02 quartzite are rounded to sub-rounded with aspect ratios 1:1 to 2:1. The examined grains exhibit bright cores in CL, commonly with sector zoning, overgrown by rims (≤20 μm thick) that are dark grey in CL (Fig. 2; spot 02R10.15). Two dominant age peaks are recorded in the detrital spectrum at 1.85 Ga and 2.65 Ga (Fig. 3; 19AM02B). The MLA algorithm yields a date of 1202 ± 40 Ma (Fig. 4c) that is interpreted as the MDA of the sample. In comparison, the three youngest dates from zircon grains with Th/U values between 0.4 and 1, that overlap at 2s, yield a weighted average ²⁰⁷Pb/²⁰⁶Pb date of 1303 ± 11 Ma.

The ϵHf values on the examined grains vary from −12.3 to +10. In more detail, two analyses of the youngest age population at c. 1.2 Ga and 1.3 Ga yield ϵHf values of −11.3 and −12, respectively. The Lu–Hf analyses on grains of the c. 1.8–2.0 Ga and 2.3–3.0 Ga age populations exhibit ϵHf values that vary from −12.3 to +10.0 and −8.6 to +7.5 (Fig. 5).

5.d. Sample 19FS7085 (Barrois complex)

The examined zircon grains are rounded to sub-rounded with aspect ratios 1:1 to 3:1. The imaged grains exhibit bright cores in CL overgrown by oscillatory zoned mantle domains (Fig. 2; 7085R3.11); domains with dark grey response in CL that exhibit oscillatory or no zoning are also observed. A dominant age peak is recorded in the sample at 1.83 Ga, with two subordinate peaks at 1.44 Ga and 2.63 Ga (Fig. 3; 19FS7085). The MLA algorithm yields a date of 1281 ± 48 Ma (Fig. 4d) that is interpreted as the MDA of the sample. The three youngest U–Pb dates with Th/U values from 0.47 to 0.84 yield a weighted average ²⁰⁷Pb/²⁰⁶Pb date of 1294 ± 64 (2s) which is consistent with the calculated MDA.

The ϵHf values on the examined grains vary from −17.8 to +10.5 (Fig. 5). Specifically, the Lu–Hf analyses in the 1.7–1.9 Ga age population show the widest spread in ϵHf values from −17.8 to +10.4 with a prevalence of negative ϵHf values (84% of the analyses). Four U–Pb analyses that fall within the 1.6–1.7 Ga age population show a smaller spread in ϵHf values from −3.66 to +5.56. Seven Lu–Hf analyses from grains of the c. 1.4–1.5 Ga population vary from −17.5 to +8.8, with a prevalence of positive ϵHf values (72 % of the analyses), whereas those with ages >2.3 Ga yielded ϵHf values varying from −15.3 to +3.7.

6. Discussion

6.a. Timeline of crustal evolution

6.a.1 Juvenile crustal growth at c. 1.4–1.5 Ga and 1.6–1.7 Ga

The suprachondritic Hf isotopic values in the c. 1.4–1.5 Ga age population support juvenile crustal additions during the accretionary growth of southeastern Laurentia or reworking of a juvenile source with short crustal residence time. In the case of the 1.6–1.7 Ga age population the presence of subchondritic ϵHf values is consistent with an older crustal component assimilated during mantle-derived magmatism. Juvenile whole-rock Nd isotopic compositions of Pinwarian and Labradorian orthogneisses exhibit ϵNd values that vary from +2.2 to +5.8 and +2.3 to +3.6, respectively (Dickin, 2000; Zelek & Dickin, 2013). Geochemical and zircon Lu–Hf isotopic data from other Pinwarian metaplutonic suites in the Grenville Province also show juvenile signatures with ϵHf values varying from 6.3 ± 0.5 to 10.7 ± 0.7 (Augland et al., 2015), supporting the Lu–Hf detrital zircon dataset of this study.

It has been proposed that juvenile plutonic rocks have low preservation potential during supercontinent assembly (Hawkesworth et al. 2010; Spencer et al. 2015). However, in our study the age peaks at c. 1.45 Ga and 1.65 Ga in the detrital spectra record accretionary orogenesis that survived the effects of the protracted collisional Grenvillian orogenesis. These peri-Laurentian accretionary events are members of the Great Proterozoic Accretionary Orogen of Laurentia that contains the largest fraction of juvenile continental crust produced globally during the Palaeo- and Mesoproterozoic (Condie, 2013). Thus, here we show that, in Central Grenville Province, the record of juvenile crust production prior to collision is preserved in vestiges of pre-orogenic basins, supporting the few U–Pb and Lu–Hf zircon studies from pristine meta-igneous lithologies (Augland et al., 2015).

6.a.2 Crustal reworking at 1.8–1.9 Ga

The vertical mixing array in ϵHf-time space is consistent with the assimilation of mantle-derived magmas with pre-existing crust during the Penokean–Makkovikian orogeny (Fig. 5). Alternatively, the suprachondritic ϵHf values at 1.8–1.9 Ga age could represent a distinct juvenile isotopic source crystallized between 1.8 and 1.9 Ga, with the subchondritic ϵHf data reflecting pure crustal reworking. Neither scenario can be excluded since juvenile intra-oceanic and continental arcs that reworked pre-existing crust were constructed between 1.8 and 1.9 Ga (e.g. Labrador; Culshaw et al. 2000) and could have shed detritus into the examined pre-collisional basins.

In our dataset, a ¹⁷⁶Lu/¹⁷⁷Hf evolution trajectory with slope of 0.011 between a near-chondritic 2.7 Ga reservoir and the end-member ϵHf value of the sub-chondritic source at 1.8–1.9 Ga (arrow line; Fig. 5) is consistent with the reworking of Archaean (i.e. Superior Province) continental crust. The ϵHf-time pattern could alternatively be interpreted as the reworking of variably contaminated mafic sources extracted from a depleted-mantle reservoir between c. 2.8 Ga and 3.6 Ga (¹⁷⁶Lu/¹⁷⁷Hf = 0.022; grey band in Fig. 5). We lean towards the former scenario, considering the development of an Andean-type convergent margin during the Penokean orogeny at the southern margin of Laurentia (Schulz & Cannon., 2007) and the presence of 2.7 Ga chondritic ϵHf detrital zircon values. The only previous Lu–Hf isotopic data from detrital zircon grains of a pre-orogenic formation come from the Siamarnekh formation with MDA at c.1.25 Ga (Spencer et al. 2015) that crops out in northern Labrador. The zircon grains of the Siamarnekh formation exhibit two dominant age peaks at 1.9 Ga and 2.15 Ga and show a vertical spread in ϵHf-time space. Interestingly, the 2.15 Ga grains fall along the crustal reworking trajectory (¹⁷⁶Lu/¹⁷⁷Hf = 0.011; Fig. 5), but further constraints about their origin and provenance have not been reported. Also, the most proximal Penokean-affinity sources, from the Western Grenville Province, that could have shed detritus towards the sequences examined here are gneisses of the Barillia terrane with Nd model ages in the 1.8–2.0 Ga range. These gneisses contain no Nd isotopic evidence of reworked Archaean crustal components and have been characterized as juvenile (Dickin, 2000). However, lead isotopic data from leached feldspar crystals from various lithologies in the western Grenville Province suggest mixing of Archaean crust with 1.8 or 1.5 Ga juvenile magmatic rocks, with the latter preferred in order to be compatible with the 1.8 Ga Nd model ages in the area (DeWolf & Mezger, 1994). More distal terranes, in the Great Lakes region and Eastern Labrador, exhibit stronger evidence for Penokean–Makkovikian reworking of Archaean basement crust. In the Great Lakes region, synorogenic Penokean plutons that intrude the Archaean Marshfield terrane show ϵNd values that vary from −1 to −7 resulting from the contamination of magmas with Archaean crust, supporting previous isotopic studies (Van Schmus & Anderson, 1977; Van Wyck & Johnson, 1997). Further eastwards, in the Makkovik Province of eastern Labrador, the Island Harbour Bay Plutonic Suite constitutes a calc-alkaline magmatic arc that reworked Archaean crust between c. 1.89 and 1.87 Ga (Ketchum et al. 2002). Similarly, Makkovikian granitoids in this area show negative ϵNd values that range from −7.2 to −2.5 as a result of variable assimilation of granitic magmas with Archaean crust of the Nain Province (Barr et al. 2001). Taking into consideration that the 1.8–1.9 Ga grains are not first cycle and there is no U–Pb isotopic evidence for 1.8–1.9 Ga magmatism proximally to the study area, we attribute the parental sources of the 1.8–1.9 Ga grains to reworked Archaean basement terranes of the Great Lakes area or Makkovik Province.

6.b. Provenance and tectonic constraints

6.b.1. Montauban quartzite

The suprachondritic Hf isotopic compositions in the zircon grains of the c.1.46 Ga age population indicate they were derived from juvenile sources. Potential sources are the c. 1.45 Ga Montauban arc (L Nadeau et al., unpub. data, 1992), an intraoceanic arc that accreted to Laurentia at c. 1.39 Ga (Sappin et al. 2009), and juvenile Pinwarian granitoids (Augland et al., 2015). The detrital spectrum in Montauban quartzite and the c. 5–10 Ma age gap between the dominant age peak and the MDA would be consistent with deposition in a convergent margin (possibly a back-arc) setting (e.g. Cawood et al. 2012). The mingling of clastic, pyroclastic and pillow basalts in the Montauban Group likely records the early stages of back-arc basin formation at c. 1.44 Ga (Fig. 6a). The Escoumins supracrustal belt is a laterally equivalent supracrustal sequence to the Montauban Group, with a similar MDA of 1476 ± 10 Ma (Groulier et al. 2018). The two sequences have age spectra with dominance from Pinwarian sources, but in the Montauban quartzite there is also a significant contribution from Labradorian sources (Fig. 6b). The minor age peak in both sequences between c. 1.8 and 1.9 Ga indicates the limited contribution from Penokean–Makkovikian sources. Comparatively, the detrital signal of Montauban and Escoumins sequences is markedly different from quartzites of the Pinwarian, Wakeham intracontinental rift basin, in eastern Grenville Province (Van Breemen & Corriveau, 2005). The unimodal detrital spectra with an age peak at 1.8–1.9 Ga in Wakeham basin quartzites (Fig. 6c) indicate provenance most probably from the Palaeoproterozoic Torngat orogen or the Makkovikian terrane. The disparate detrital age distributions demonstrate the different tectonic setting and catchment areas of the Pinwarian basins along the orogen.

Fig. 6. (a) Tectonic setting (not in scale) of the examined quartzite from the supracrustal sequence of the Montauban Group; (b) KDE plots of U–Pb dates from the Montauban and Escoumins supracrustal group (Groulier et al. 2018); (c) KDE plots of U–Pb dates from the Wakeham basin (U–Pb data from Wodicka et al. 2003 and Van Breemen & Corriveau., 2005); (d) tectonic setting (not in scale) of the Wabash–Barrois supracrustal groups; (e) KDE plots of U–Pb dates from the Wabash–Barrois and Saint-Boniface supracrustal groups (U–Pb data from this study and Peck et al. 2019); (f) KDE plots of U–Pb dates from supracrustal sequences of the Adirondack lowlands (U–Pb data from Chiarenzelli et al. 2015; and pers. comm. with J. Chiarenzelli).

6.b.2. Wabash and Barrois quartzites

The age gap between the MDA and the major age peaks at 1.8–1.9 Ga in Barrois and Wabash complexes is more than 500 Ma. Specifically, c. 97 % of the analyses record an age gap of more than 150 Ma between the dominant age peak and the MDA (Fig. 6d), suggesting that both sequences were deposited in intracratonic extensional sedimentary basins (e.g. Cawood et al. 2012). In terms of proto-sources, Penokean–Makkovikian (c. 1.8–1.9 Ga) lithologies to our knowledge are not exposed in the study area. Therefore, the 1.8–1.9 Ga detrital zircon grains, hosted in the Wabash and Barrois quartzites, either constitute distal detritus deposited in intermediate sedimentary reservoirs and recycled during Elzevirian rifting (1.2–1.3 Ga) or were transported via orogen-parallel drainage (present coordinates), possibly in large fluvial systems. U–Pb detrital zircon spectra with dominant age peaks between c. 1.8 and 1.9 Ga and MDA at 1.2–1.3 Ga, similarly with the Wabash and Barrois complexes, are observed within different crustal domains of the Grenville Province. Specifically, southwestwards of the study area, in the Morin terrane, the Saint-Boniface metasedimentary sequence (Corrigan & Van Breemen, 1997) with MDA at c. 1.25 Ga (Peck et al. 2019) shows dominant age-peaks at c. 1.25 and 1.82 Ga (Fig. 6e). These age peaks also indicate the feeding of these basins dominantly from Penokean and Elzevirian sources. Further west, in the Adirondack region, three quartzites from metasedimentary sequences with MDAs at c. 1.25–1.3 Ga contain a higher proportion of Elzevirian and a lower proportion of Penokean detrital zircon grains than those found in Quebec (Fig. 6f; Chiarenzelli et al. 2017). These 1.2–1.3 Ga basins support inferences from bimodal volcanic and AMCG plutonic suites which indicate that extensional tectonism was occurring during that period in the Grenville Province and on a global scale at the core of Nuna supercontinent (Evans & Mitchell, 2011; Indares & Moukhsil, 2013).

7. Conclusions

From the examination of three poorly exposed Mesoproterozoic supracrustal complexes from Precambrian basins of the Central Grenville Province we show that:

• Detrital zircon U–Pb data from a quartzite hosted within a metavolcanic sequence (Montauban complex) indicates the development of a Pinwarian back-arc basin with maximum depositional age of 1442 ± 8 Ma (2s).

• The detrital zircon grains of the Montauban quartzite record two dominant age peaks at c. 1.46 Ga and 1.62 Ga that exhibit suprachondritic ϵHf values recording juvenile crust formation during the Pinwarian (c. 1.4–1.5 Ga) and Labradorian (c. 1.6–1.7 Ga) orogenies.

• The U–Pb detrital zircon spectra from two supracrustal complexes in the Lac Saint-Jean region (Wabash and Barrois complexes) with MDA of c. 1.2–1.3 Ga suggest deposition of 1.8–1.9 Ga recycled or far-travelled detrital zircon grains in an intracontinental extensional basin.

• The sub- and suprachondritic ϵHf values found in the 1.8–1.9 Ga detrital zircon grains manifest the reworking of a near-chondritic 2.7 Ga continental basement during the Andean-type, Penokean–Makkovikian orogeny.

• The accretionary stages of supercontinent assembly are retained in the detrital zircon archive of pre-collisional basins in Central Grenville Province even after the terminal Grenvillian collision at c. 1 Ga.