Deciphering crustal growth in the southernmost Arabian Shield through zircon U-Pb geochronology, whole rock chemistry and Nd isotopes

ABSTRACT New U-Pb zircon geochronology using high-spatial resolution secondary ion mass spectrometry fills a data gap and provides crystallization ages for granitoids from the Asir composite terrane in the southernmost Arabian Shield of Saudi Arabia. Ages of c. 810–685 Ma, c. 663–636 Ma, and 625–610 Ma reflect oceanic island arc genesis, subduction-related arc accretion (syn-collisional), and post-collisional stabilization, respectively. All samples have juvenile εNd(t) compositions with no evidence of older material being involved in their genesis, indicating that this part of the Arabian Shield grew through juvenile magmatic addition and that assimilation by syn- and post-tectonic magmatism involved an isotopically juvenile component(s). The crustal thickness derived from the (La/Yb)N proxy indicates significant thickening from 10–20 km to c. 70 km at c. 650 Ma, consistent with timing of orogenic uplift and increasing crustal thickness post-dating peak Nabitah orogeny. The age of an intrusion cross-cutting the Atura formation, when combined with other data, provides a well-constrained depositional age of c. 646–625 Ma for the Atura formation and indicates that erosion of the orogenic edifice in this part of the Arabian Shield began at latest by 625 Ma. Our new data indicate that denudation occurred 80–100 m.y. before the development of the prominent sub-Cambrian peneplain, consistent with previous assertions that major pulses of denudation occurred prior to the waning stages of Nabitah orogenesis.


Introduction
The amalgamation of Gondwana occurred through numerous collisional events, including that which generated the well-known East African Orogen (Stern 1994;Jacobs and Thomas 2004;Collins and Pisarevsky 2005). The Arabian-Nubian Shield (ANS) occupies a large part of the East African Orogen ( Figure 1) and represents one of the greatest volumes of Neoproterozoic juvenile crust preserved on Earth (Patchett and Chase 2002). The ANS is therefore an ideal place to study Neoproterozoic crustforming processes and the tectonic mechanisms associated with the formation of Gondwana.
The protracted evolution of the ANS involved the accretion of oceanic island arcs, remnants of ophiolites, and oceanic plateaus from within the Mozambique Ocean, as well as microcontinental fragments (Stern 1994;Stein and Goldstein 1996;Johnson et al. 2011). Island-arc formation and accretion were followed by shield-wide emplacement of post-tectonic calc-alkaline and alkaline plutons during final consolidation of the shield (c. 630-550 Ma) (Bentor 1985;Be'eri-Shlevin et al. 2009a;Johnson et al. 2011;Fritz et al. 2013;Robinson et al. 2014;Yeshanew et al. 2015). While the core of the ANS is a collage of juvenile island-arc terranes, its margins were influenced by either pre-Neoproterozoic crust reworked during Neoproterozoic amalgamation or by the involvement of older continental material (Agar et al. 1992;Windley et al. 1996;Teklay et al. 1998;Whitehouse et al. 1998;Be'eri-Shlevin et al. 2009b;Yeshanew et al. 2015Yeshanew et al. , 2017. Across the ANS terranes and their boundaries reveal Neoproterozoic crust-forming processes and elucidate the tectonic mechanisms associated with terrane assembly. Terrane correlation, however, requires detailed knowledge of rock types, ages, metamorphic grade, etc. Unfortunately, key locations within the ANS lack the critical data that allow such correlations to be made. For example, near the border of Saudi Arabia and Yemen where extensive Tonian to Cryogenian basement of the southern Arabian Shield is well exposed (Figure 2a), rugged topography combined with an unstable political situation along the border means that this region remains little studied. Elsewhere, exposure is poor. While a terrane framework for Saudi Arabia exists, there is a paucity of data on which to extrapolate this framework south (Johnson and Kattan 2012 and references therein). A synthesis characterizing the terranes of northwestern Yemen has recently been published (Al-Khirbash et al. 2021), but given the general lack of data in the region the correlation to the terranes within Saudi Arabia remains speculative.
We present new whole-rock chemistry, zircon U-Pb ages, and Nd isotopic data from granitoids in the southernmost part of the Asir terrane near the Saudi Arabia-Yemen border (Figure 2b). These new data are used to characterize the age of juvenile intra-oceanic island-arc granitoids, subduction-related syntectonic granitoids associated with arc accretion, as well as the postcollisional granitoids associated with consolidation of the southernmost Arabian Shield. Integrating these data with previously published work, we are able to refine the tectonomagmatic evolution of the southern ANS in relation to the timing of crustal thickening and erosion of the orogenic edifice associated with Nabitah orogenesis.   Stern et al. 2010;Johnson et al. 2011). Terrane boundaries and terrane names shown in red. Inset shows location of study area within the Asir composite terrane, southern Saudi Arabia (after Johnson and Kattan 2012). The Ablah belt is a younger (post-Neoproterozoic) feature.

Geological background
The juvenile Neoproterozoic crust of the ANS comprises over 610,000 km 2 across northeast Africa and the western part of the Arabian Peninsula ( Figure 2a) and is best exposed along the uplifted flanks of the Red Sea (Fleck et al. 1980;Stern 1994;Johnson and Woldehaimanot 2003;Hargrove et al. 2006). It occupies the northern segment of the East African Orogen (Stern 1994) and is believed to have formed by juvenile magmatic addition between c. 850-570 Ma (Stoeser and Camp 1985;Stoeser and Frost 2006;Stern et al. 2010). The core of the Arabian Shield consists of Tonian to Cryogenian island-arc terranes formed in the Mozambique Ocean and later accreted/obducted during the assembly of Gondwana (Stoeser and Camp 1985;Johnson and Woldehaimanot 2003;Stoeser and Frost 2006;Johnson et al. 2011). The vast tracts of juvenile island arcs of the Arabian Shield comprise numerous, distinct terranes separated by ophiolitebearing suture zones (Stoeser and Camp 1985;Stoeser and Frost 2006;Johnson et al. 2011) and document distinct tectonomagmatic cycles. Thick and extensive oceanic tholeiite generation began at c. 900 Ma and lasted c. 20 m.y. -this marks the earliest magmatism in the Arabian Shield and has been interpreted as the shallow expression of a 'plume head' (Bentor 1985;Stein and Goldstein 1996;Stein 2003). Andesitic volcanism, calc-alkaline plutonism with island-arc affinity, and syntectonic plutonism occurred between c. 870-650 Ma (Fleck et al. , 1980Bentor 1985;Stern 1994;Johnson and Woldehaimanot 2003;Stein 2003). Post-collisional calc-alkaline and alkaline batholiths were emplaced at c. 640-550 Ma and are volumetrically more important in the northern and eastern Arabian Shield. Post-tectonic magmatism was integral to the final consolidation of the shield (Bentor 1985;Black and Liegeois 1993;Stern 1994;Be'eri-Shlevin et al. 2010;Robinson et al. 2014;Yeshanew et al. 2015).
The Arabian Shield in southern Saudi Arabia is represented by the Asir composite terrane (Stoeser and Camp 1985), which comprises north-trending structural belts of island-arc material separated by shear zones (Figure 2b). The Asir composite terrane is typically divided into the 'western' and 'eastern' terranes along the Nabitah fault zone (NFZ; also known as the Nabitah suture). The NFZ formed during Nabitah orogenesis at c. 680-640 Ma and was historically considered to be the suture associated with the assembly of Gondwana between western juvenile terranes and the more evolved eastern arc and continental terranes (Johnson et al. 2011). However, it has been clearly demonstrated that there is little difference between terranes adjacent to the NFZ and therefore the suture between west and east Gondwana must be east of the NFZ (Flowerdew et al. 2013).
Following the nomenclature of Johnson and Kattan (2012), the western terranes are represented by the Al Lith assemblage, the Bidah belt, the Shwas and Tayyah belts (equivalent to the An Nimas terrane of Stoeser and Frost 2006), and the Khadrah belt (equivalent to the Al Qarah terrane of Stoeser and Frost 2006) (Figure 2b). The eastern terrane is represented by the Malahah belt (partially correlative with the composite Tathlith-Malalah volcanic arc terrane of Stoeser and Frost 2006) (Figure 2b). The rocks within each belt have distinct compositions, depositional facies, metamorphic grades, and deformational styles (Fleck et al. 1980;Greenwood et al. 1982). At present it is not known if each belt represents distinct arcs that collided to form the composite terrane or if they constitute a smaller number of arcs repeated through folding and faulting and superposed by N-S shear.
The structural fabric of the Asir composite terrane clearly strikes into northwest Yemen, but the correlation of units is tenuous due to the lack of data, access and/or exposure. For example, the location of the NFZ in Saudi Arabia south of ca. latitude 20°N is obscured by intrusion of the Wadi Tarib batholith (Johnson and Stewart 1995) and in Yemen, much of the projected trace of the NFZ is hidden beneath younger rocks (Windley et al. 1996). The Shwas and Tayyah belts in Saudi Arabia continue into Yemen and should be directly linked to the Haydan terrane of Yemen, but there are no data to confirm this correlation (Al-Khirbash et al. 2021). The At Talh composite terrane is the purported extension of the Khadrah and Malahah belts in Yemen, but most of the At Talh is hidden beneath younger rocks. Detailed geochronological and isotopic investigations are clearly needed to verify the correlation of these terranes.
The focus of this study is on rocks from the southernmost Tayyah and Khadrah belts, which are separated by the Junaynah fault zone (Johnson and Kattan 2012) ( Figure 2b). In both belts the rocks are generally metamorphosed to greenschist facies. When proximal to large batholiths or gneiss domes, granulite facies may be reached confirming a genetic relationship between emplacement and metamorphism (Fleck et al. 1980;Greenwood et al. 1982).
The Tayyah belt is an early Cryogenian assemblage of metavolcanic and metasedimentary arc rocks intruded by numerous plutonic bodies and gneiss domes. The An Nimas batholith is one of the largest calc-alkaline intrusions in the southern Arabian Shield and consists of gabbro, diorite, tonalite, trondhjemite and granodiorite. It intrudes the Tayyah belt north of the study area; previously reported Rb-Sr ages for the An Nimas batholith range from 837 ± 50 Ma to 818 ± 95 Ma (Fleck et al. 1980). Cooper et al. (1979) present a more robust U-Pb age for the main phase of the An Nimas batholith of 816 ± 4 Ma (95% confidence upper-intercept age, zircon population method), which also provides a minimum age for the Tayyah belt. Fairer (1985) identified two generations of granitic gneiss, which also intrude the older layered volcanic arc assemblages and gneissic rocks of the Tayyah belt, but which post-date the An Nimas batholith. His 'older' generation includes the Wadi Baqarah dome and the Wadi Tarib gneiss dome (both north of the study area). Cooper et al. (1979) provided U-Pb zircon (population method) ages for the Wadi Baqarah dome (c. 797-763 Ma) and the Wadi Tarib gneiss dome (714 ± 28 Ma; 95% confidence). Flowerdew et al. (2013) reported a similar U-Pb zircon (secondary ion mass spectrometry, SIMS) age of 716 ± 8 Ma for a tonalite gneiss from the Tarib arc. Fairer (1985) grouped the granitic Harisi dome and the Shada pluton ( Figure 3) into this 'older' generation because they have a well-developed fabric and on the basis of a poorly constrained 2-point Rb-Sr isochron for the granitic Harisi dome of 773 Ma (Fleck et al. 1980). The 'younger' generation is represented by the highly sheared intrusion of the Khamis Mushayt complex with an age of 663 ± 4 Ma (2σ) (our weighted average of four ages reported by Cooper et al. 1979) (Figure 3).
In the southern part of the Tayyah belt, greywacke, quartzite, marble and lesser amounts of chert are overlain by thick-pillow basalts. However, the metasedimentarymetavolcanic contact is not obvious and they may be contemporaneous (Greenwood et al. 1982). Pillow basalts are widespread throughout the Tayyah belt and the most voluminous of these is best preserved in the south close to the Yemen border (Greenwood et al. 1982). Basaltic and andesitic pyroclastic rocks and sub-areal welded tuffs are more abundant in the northern Tayyah belt (Greenwood et al. 1982). The earliest plutonism in the Tayyah belt post-dating the An Nimas batholith includes the tonalitic to granitic gneisses of Wadi Baqarah (north of Abha) at c. 797-763 Ma (Cooper et al. 1979).
The Khadrah belt is slightly younger and consists of strongly deformed, greenschist and amphibolite facies volcanic and sedimentary rocks dated to 746 ± 16 (Rb-Sr isochron; Fleck et al. 1980). The rocks are discontinuously exposed because they are extensively intruded by syntectonic granitoids associated with Nabitah orogenesis, such as the Tarib batholith ( Figure 3). The Tarib batholith consists of diorite, tonalite, and granodiorite and has two reliable U-Pb upper-intercept ages (zircon population method, 95% conf.) that are within uncertainty of each other at 732 ± 3 Ma and 729 ± 3 Ma (Stoeser et al. 1984). The Tarib plutonic suite is correlated with the Hulayfah group supracrustal rocks, which are together known as the Tarib arc assemblage. The Tarib arc existed in the north, central and southern shield from 785-680 Ma (Stoeser et al. 1984 and references therein). Cooper et al. (1979) dated a foliated tonalite intruding the Tarib batholith to 663 ± 8 Ma (95% confidence upper-intercept age, U-Pb zircon population method), an age similar to the Khamis Mushayt complex (c. 663 Ma) presented above (Figure 3), suggesting that the syntectonic magmatismlasted a bit longer.
The arc assemblages are mostly comprised of metamorphosed andesites (flows and pyroclastics) and dacites (Greenwood et al. 1982;Johnson and Kattan 2012). As in the Tayyah belt, basalts predominate in the southern part of the Khadrah belt; in the northern Khadrah belt, the upper sections of these sequences contain rhyolite (Greenwood et al. 1982). Rb-Sr ages of 785 ± 96 Ma to 746 ± 16 Ma are reported for these volcanic rocks (Fleck et al. 1980). The oldest reported age for plutonic rocks intruding the Khadrah belt is c. 730 Ma (Fleck et al. 1980).

Units and sampling
Nine samples were collected from the southernmost parts of the Tayyah and Khadrah belts for geochemical and geochronological investigations ( Figure 3). These include five samples from the Tayyah belt and four from the Khadrah belt, providing a total of eight plutonic rocks and one igneous clast from a polymict conglomerate of the Atura formation in the Khadrah belt (Supplementary Figure 1). Sample numbers, coordinates, rock types, units, and simplified petrography are summarized in Supplementary  Table 1. The samples represent a variety of units, span a compositional range, and are presumed to have formed at different stages during the tectonic evolution of the region.
Tayyah belt samples. These five samples include a homogeneous, weakly foliated sample (FGY12-04; Supplementary Figure 1) collected from the northern part of the Harisi dome (Jabal Harisi). This sample represents the older generation of deformed, twofeldspar granitoids (unit gs) that intrude the early basaltic and andesitic assemblages of the Tayyah belt (Fairer 1985). A highly sheared quartzofeldspathic sample (FGY12-16; Supplementary  Figure 1) within the Sabya Formation represents the gneissic country rocks (unit sa) into which the Shada pluton intrudes. Two samples of deformed granitoids at Jabal Banī Mālak (FGY12-05) and Jabal Fayfā (FGY12-07) were also collected; these are mapped as small syenite bodies (unit sy) by Fairer (1985) and intrude the Sabya Formation. Field relations suggest that sample FGY12-05 is from a bimodal pluton with mafic and felsic zones, with the latter containing mafic enclaves. This sample is from the mafic part, with distinctive biotite clots surrounded by plagioclase, and cut by small granitic veins (Supplementary Figure 1). FGY12-07 at Jabal Fayfā has numerous mafic enclaves and 'ghostly' (resorbed) mafic 'stringers'; the enclaves include altered pyroxenes mantled by amphibole (Supplementary Figure 1). This sample represents homogeneous (enclave-free) host material. These two samples may represent the younger intrusions of Fairer (1985). An undeformed, two-feldspar granite (FGY12-08) was also collected at Jabal al Qarahthis sample (unit gdn) cross-cuts the local rock units and structures, and may be a post-tectonic intrusion (Supplementary Figure 1).
Khadrah belt samples. These four samples were collected near the southwestern part of the Tarib batholith. They include an undeformed, hornblendebearing two-feldspar granite (FGY12-09; unit gm) that cross-cuts the Atura formation (unit ma), as well as a rounded trondhjemitic clast (FGY12-12) collected from a poorly sorted conglomerate within the Atura formation. The Atura formation is interpreted as the deposits of a small intermontane basin and these two samples should constrain the age of the Atura formation. A deformed granodiorite (FGY12-11; unit gdn) intrudes the tonalite-granodiorite phase (tdg) of the Tarib batholith and was also collected; this sample should confirm whether the two generations of intrusions identified by Fairer (1985) occur in both belts or are restricted to the Tayyah belt. The final sample (FGY12-10; unit gm) is a late, fine-grained and undeformed granodiorite that cross-cuts the local units and structures of the Khadrah belt and is, therefore, presumed to be a post-tectonic intrusion.

Whole-rock major and trace element geochemistry
X-ray fluorescence (XRF) analyses for major elements were performed using a Rigaku Primus II at the PetroTectonics analytical facility, Stockholm University. Unaltered portions of the samples were powdered using a stainless steel swing mill. Sample preparation procedures were similar to those described in Yeshanew et al. (2015). The XRF analytical protocols are similar to other laboratories (e. g. Johnson et al. 1999). Major element analytical precision is generally 1-2%. The homogeneous fused glass disks prepared for major element analyses were also used to analyse trace elements via laser ablation inductively coupled plasma mass spectrometry (LA-ICP-MS) at the PetroTectonics analytical facility, Stockholm University, where a New Wave 193 excimer laser is coupled to a Thermo X-series 2 quadruple ICP-MS. A spot size of 120 µm was ablated using a laser pulse frequency of 10 Hz with a laser energy density of 7 Jcm -2 . Each sample was analysed three times for adequate precision. The samples were bracketed by analyses of the international reference glass NIST 612 for external calibration and elemental quantification was made using Si. Accuracy was monitored by analyses of the international reference material BCR-2.

Zircon U-Pb geochronology
Zircon was separated and concentrated from c. 2 kg of sample using standard hydrodynamic and heavy liquid (methylene iodide) techniques. After rinsing with deionized H 2 O, zircon crystals were mounted in epoxy together with the international reference zircon standard 91,500 (Wiedenbeck et al. 1995). The epoxy mount was polished to expose the crosssectional area of the zircon. Cathodoluminescence (CL) imaging of zircon was performed using an Hitachi S4300 Scanning Electron Microscope (SEM) fitted with an ETP Semra Robinson CLdetector at the Swedish Museum of Natural History.
U-Th-Pb analyses of zircons were performed using a large geometry CAMECA IMS-1280 Secondary Ion Mass Spectrometer (SIMS), Swedish Museum of Natural History following the methodology described by Whitehouse et al. (1999) and Jeon and Whitehouse (2015). The O 2 − primary ion with an intensity of ~10 nA was accelerated at −13 kV. Pre-sputtering, centring of the secondary ion beam, mass calibration and secondary ion energy optimization were automated for all analyses. The elliptical spot size from the analytical beam was c. 18 × 15 µm. In the secondary ion beam optics, a 45 eV energy window was used. A mass resolution of ~5400 was used to discriminate the Pb isotope peaks from those of interfering species. Common Pb corrections were made using the measured nonradiogenic 204 Pb and when significant (>3x detector background) corrections were made using the present-day model terrestrial Pb composition (Stacey and Kramers 1975). Analyses of the reference zircon 91,500 (1065 Ma; Wiedenbeck et al. 1995) were interspersed with analyses of the unknown zircon to monitor instrumental drift. Decay constants are those of Steiger and Jäger (1977). Data reduction was done using an in-house spreadsheet. Age calculations and data plotting were performed using Isoplot 3.75 (Ludwig 2012). Individual ages in the data table are reported with 1σ (percentage) uncertainties and the results of age calculations are given at 2σ error limits, with concordia ages presented using the combined MSWD of concordance and equivalence following the recommendation of Ludwig (1998).

Whole-rock Nd isotopes
Whole-rock powders were digested in a mixture of concentrated HF and HNO 3 within sealed pressurized containers and heated at 205°C for 72 hours. Ten drops of concentrated HNO 3 were added to each sample capsule and evaporated. This procedure was repeated 3-4 times to break down the fluorides, then 6 N HCl was added to each capsule and the samples were re-heated at 205°C for another 24 hours. Standard column chemistry was employed to separate Sm and Nd. Determinations of Sm and Nd isotope ratios were made on a Thermo Scientific TRITON thermal isonization mass spectrometer using the total spiking method with a mixed 147 Sm/ 150 Nd spike at the Swedish Museum of Natural History. Detailed sample preparation and analytical conditions closely resemble those described in Yeshanew et al. (2017). Samarium concentrations were determined in multi-collector static mode, neodymium was run in static mode, and both elements were run using rhenium double filaments. A ratio 149 Sm/ 148 Sm = 1.22937 was used to normalize Samarium ratios. A normalization value of 146 Nd/ 144 Nd = 0.7219 was used for the calculated Nd isotope ratios. The external precision for 143 Nd/ 144 Nd, evaluated using the La Jolla reference material, was 11 ppm. An accuracy correction was not necessary as the mean 143 Nd/ 144 Nd ratio (0.511861 ± 55, n = 12) was within error of the accepted value of 0.511859 (Pier et al. 1989).

Element geochemistry
Major and trace element data for these nine samples are given in Supplementary Table 2. The samples represent different units and are not cogenetic.

Major elements
Major elements for these samples range in SiO 2 from 54-72 wt%, in Al 2 O 3 from 14.4 to 16.8 wt%, and combined K 2 O+Na 2 O = c. 8 wt% (Supplementary Table 2). Due to the coarse grain size of the samples, the samples are named on the basis of their normative An, Or, and Ab mineralogy (after Barker, 1979) Table 2); samples FGY12-05 and FGY12-07 are the exceptions and these are strongly metaluminous (Supplementary  Table 2).

Trace elements
Trace element (including rare earth elements) data of the ferroan and magnesian compositional groups are plotted separately and normalized to the N-MORB values of Sun and McDonough (1989) in Figure 5. Relative to N-MORB, the large ion lithophile elements (e.g. Rb, Ba, K) and high field strength elements (e.g. Nb, Hf, Zr) are generally enriched (e.g. Ba N = 58-365 and Nb N = 0.7-18; Figure 5a). Most samples have negative Sr and P anomalies, whilethe ferroan samples also have small negative Ti anomalies (Figure 5a). Ta-Nb troughs also occur in most samples, but are more pronounced in the magnesian samples (Figure 5a).
With respect to N-MORB-normalized rare earth elements (REE), both the ferroan and magnesian samples are enriched in light (L)REE relative to heavy (H)REE (up to 35x N-MORB); a single exception is FGY12-16 which has LREE similar to, or less than, N-MORB. The REEs are strongly fractionated (Figure 5b). The HREEs are generally depleted relative to N-MORB (Ho N to Lu N < 1), with flat to slightly negative slopes ([Ho/Lu] N = 1.0 to 1.7) (Figure 5b), but there are two exceptions (FGY12-04 and FGY12-16) -these samples represent the 'older' gneisses of Fairer (1985) and both samples have increasing HREE from Ho to Lu (Figure 5b). Sample FGY12-04 is unique in having a pronounced positive Eu anomaly [Eu* = 4.7] (Eu* = Eu N /√(Sm N x Gd N ), whereas Eu anomalies associated with the other samples are nominal (Eu* = 0.64 to 1.59) (Figure 5b).

Zircon U-Pb geochronology
The analytical data are presented in Supplementary Table 3, CL images of representative grains with analytical locations are shown in Supplementary Figure 2, and concordia diagrams are compiled in Figure 6.

Sample FGY12-04: gneissic granite (gs)
Zircons from this sample are yellowish-brown, elongated and euhedral with sizes ~100-220 µm. They often have a dark CL response, but when core-rim structures are evident, the cores are usually brighter than the rims. Most grains show simple oscillatory zoning. A total of 30 analyses were performed and overall show a slight reverse discordance, likely reflecting overcorrection of common Pb. Of the seven discordant analyses, two are reversely discordant (n4924-4, n4924-7) and these are excluded from the final age determination. The five younger discordant analyses are interpreted to represent Pb-loss. Excluding these seven discordant analyses, the remaining 23 analyses yielded a concordia age of 772 ± 4 Ma (MSWD = 1.5; Figure 6a), which is interpreted to be the igneous crystallization age of the sample. The Pb-loss recognized in this and the remaining samples is attributed to a regional post-crystallization resetting due to Nabitah orogenesis (c. 680-640 Ma) and/or Najd faulting (c. 620-525 Ma) (Johnson et al. 2011;Robinson et al. 2014).

Sample FGY12-05: granite (sy)
This sample yielded only a few zircons that were colourless euhedral prisms, ~ 100-135 µm long, with aspect ratios of 2:1 to 3:1. Grains show typical sector zoning with variable CL response. All four zircons from this sample were analysed with one of them being reversely discordant due to the analytical spot partially overlapping with epoxy. The remaining three analyses yielded a concordia age of 648 ± 9 Ma (MSWD = 1.5; Figure 6b) and this is interpreted to be the crystallization age of the sample.

Sample FGY12-07: granite (sy)
Most zircons from this sample are translucent, euhedral, and stubby. These zircons are c. 100 to 215 µm long, with aspect ratios of 1:1 to 4:1. The majority of grains display sector zoning, but some have oscillatory zoned rims. Thirteen grains were analysed and a single discordant analysis was excluded. The 12 remaining analyses yielded a concordia age of 641 ± 3 Ma (MSWD = 1.08; Figure 6c) which is interpreted as the crystallization age of the sample.

Sample FGY12-08: granite (gdn)
The zircons are light pink to light brown, euhedral and elongate, with aspect ratios of 3:1 to 8:1. In CL the crystals are generally medium-bright with weakly zoned central portions, mantled by thin and oscillatory zoned outer portions. Oscillatory growth zoning shows smooth transitions without notable truncations and suggests continuous growth during the same crystallization sequence rather than rim overgrowth. Of 17 analyses, a single reversely discordant analysis (reflecting overcorrection of common lead) and five analyses defining recent Pb-loss were excluded from the final age calculation. The remaining 11 analyses yielded a concordia age of 636 ± 6 Ma (MSWD = 1.7; Figure 6d).
modern Pb-loss array. The remaining eight analyses combine to yield a concordia age of 625 ± 5 Ma (MSWD = 1.6; Figure 6e).

Sample FGY12-10: granodiorite (gm)
Zircons from this granodiorite sample are elongate, euhedral prisms. They are 90 to 150 µm long and 35 to 60 µm wide, with aspect ratios of 2:1 to 4:1. These zircons show variable CL response with generally irregular zoning and some grains appear bleached. Of 13 analyses, four are reversely discordant and these are excluded from the final age determination. The remaining nine analyses yield a concordia age of 613 ± 3 Ma (MSWD = 0.98; Figure 6f).

Sample FGY12-11: granodiorite (gdn)
Zircons extracted from this sample are euhedral and range from stubby to elongate; they are 110 to 190 µm long and 30 to 95 µm wide, giving aspect ratios of 2:1 to 5:1. In CL, most zircons show simple oscillatory zoning with a thin dark mantle. A total of 17 analyses were performed and eight of these define a horizontal array interpreted to reflect modern Pb-loss. The remaining nine analyses yielded a concordia age of 663 ± 6 Ma (MSWD = 1.5; Figure 6g).

Sample FGY12-12: trondhjemite clast (in ma)
Zircons from this sample are euhedral, equant to elongate, 65 to 160 µm long and 20 to 48 µm wide, giving aspect ratios of 1.5:1 to 7:1. In CL, they show variable internal structures and luminescence. Both zoned and unzoned grains are present; internal domains are often weakly zoned and mantled by oscillatory zoned outer portions/rims. A total of 27 analyses were performed. Rejecting four analyses that are reversely discordant and 13 that represent a combination of ancient and modern Pb-loss, the remaining 10 analyses yielded a concordia age of 685 ± 3 Ma (MSWD = 1.8; Figure 6h).

Sample FGY12-16: trondhjemite (in sa)
Zircons extracted from this sample are euhedral, 78-206 µm long and 22-80 µm long, giving aspect ratios of 2:1 to 5:1. In CL, these grains are bright to mediumbright and most grains show simple oscillatory zoning with some core-rim relationships evident. Twenty analyses were performed and two of these were excluded from the final age calculation: one is interpreted to be xenocrystic and the other is discordant. The remaining 18 analyses yielded a combined concordia age of 811 ± 5 Ma with MSWD = 1.4 (Figure 6i) which is interpreted to be the crystallization age of the sample.  Pearce (1983) with maximum mobility at Ba-Th; both plots are normalized with the N-MORB values of Sun and McDonough (1989).

Whole-rock Nd isotope data
Whole-rock powders from the nine samples were used for Nd isotopic analysis (Supplementary Table 4). Their Sm and Nd concentrations, determined by isotope dilution, range from c. 1 to 16 ppm and c. 6-91 ppm, respectively. Their 147 Sm/ 144 Nd measured isotopic values range from 0.0881 to 0.1416 and for 143 Nd/ 144 Nd from 0.512437 to 0.512702. Initial epsilon Nd (εNd (i) ) values were calculated for these samples using Figure 6. U-Pb inverse concordia diagrams for southern Asir granitoids. Individual analysis errors are indicated by ellipses at 2σ (black outlines = analyses included in the final age determination; grey outline = excluded due to inferred Pb-loss, discordance, or large errors). All ages represent concordia ages (filled grey ellipses) and are quoted at 2σ. MSWD (mean square weighted deviation) = con-cordance+ equivalence for concordia ages. Horizontal arrows suggest modern Pb-loss. their respective U-Pb zircon crystallization ages. All samples have juvenile εNd (i) values ranging from +3.8 to +6.9.

Discussion
Our new data from the Asir composite terrane represent about 200 m.y. of tectonomagmatic evolution in the southernmost Arabian Shield. These data provide new and/or more precise age constraints and the integration of these data with elemental and isotopic geochemistry elaborates the petrogenic evolution of these juvenile arc granitoids. The data further allow us to constrain the tectonic and crustal evolution of the Asir terrane within the context of Nabitah orogenesis and the ultimate assembly of Gondwana.

Tayyah belt
Our re-sampling of gs at Jabal Harisi (FGY12-04) defines an age of 772 ± 4 Ma for island arc diorite-trondhjemite plutonism, which is similar to, but better constrained than, a 2-point Rb-Sr isochron upper-intercept age c. 763 ± 53 Ma (Fleck et al. 1980). Similar ages reported from the Khadrah belt are discussed below.
Our new data confirm a Cryogenian age for a number of small, deformed syentitic bodies (sy) mapped in the region, including the Banī Mālik pluton (FGY12-05) with an age of 648 ± 9 Ma and the Fayfā pluton (FGY12-07) with a similar age of 641 ± 3 Ma; these new data define the age for this phase of magmatism and are in agreement with the age of 646 ± 2 Ma (U-Pb zircon SIMs; Robinson et al., 2018) for a late granitic dike that also intrudes the Tayyah belt. The foliated granodioritic gneiss (gdn) at Jabal al Qahrah (FGY12-08) south of the Khamis Mushayt complex has a similar age of 636 ± 6 Ma, so although gdn is supposedly the mafic phase of the Khamis Mushayt complex, at Jabal al Qaharah, gdn post-dates the Khamis Mushayt complex by c. 30 Ma (dated by Cooper et al. 1979 to 663 ± 4 Ma -our weighted average of their four ages). We conclude that not all foliated granodioritic rocks are related to the Khamis Mushayt complex. In the southern part of the Asir composite terrane, units mapped as sy and gdn have similar ages (within error) and therefore foliated granodiorites seem to be coeval with Banī Mālik and Fayfā magmatism. These ages are consistent with emplacement ages of c. 660-640 Ma for syn-collisional plutonism within the Khadrah belt north of the study area and from the adjacent Tathlith-Malahah terranes (Flowerdew et al. 2013).
In addition, our new data from FG12-16 constrain the trondhjemitic quartzo-feldspathic gneiss (presumably intruding unit sa) to 811 ± 5 Ma. Notably, this is within error of the An Nimas batholith (816 ± 4 Ma; Cooper et al. 1979) and confirms that age-equivalent arc-plutonism extends to the border with Yemen.

Khadrah belt
The age of 685 ± 3 Ma for a trondhjemitic conglomerate clast (FGY12-12) from the Atura formation constrains its maximum age of deposition (MAD; c. 685 Ma), while the cross-cutting granitic intrusion (gm, FGY12-09) defines its minimum age (625 ± 5 Ma). Cooper et al. (1979) reported an age of 663 Ma for the quartz diorite upon which the Atura formation is unconformably deposited and the youngest detrital zircon age peak in the Atura formation defines an MAD of 646 Ma (Robinson et al. 2018). These data combined provide a well-constrained depositional age range for the Atura formation of 646-625 Ma.
The other gm sample (FGY12-10) that also cross-cuts the local units and structures of the Khadrah belt is younger (613 ± 3 Ma) than FGY12-09 (625 ± 5 Ma). The deformed granodiorite (gdn, FGY12-11; unit) has an age of 663 ± 3 Ma, consistent with the fact that it must be younger than the c. 730 Ma Tarib batholith (tdg) which it intrudes.
Our oldest island arc sample (trondhjemite, FGY12-16) yields a U-Pb crystallization age of 811 ± 5 Ma and refines the age estimate from previous studies which suggested that the earliest plutonism represented by diorite-trondhjemite batholiths commenced prior to c. 800 Ma (Greenwood et al. , 1982Fleck et al. 1980). Previous mineral dating (K-Ar and whole-rock Rb-Sr) from the southwestern and western Arabian shield identified c. 900-800 Ma diorite-trondhjemite batholiths, with younger ages of c. 750-680 Ma in south-central and eastern parts of the shield (Fleck et al. 1980;Fairer 1985). North of the study area, a tonalite gneiss from the Khadrah belt with an age of 716 ± 8 Ma (95% conf.; SIMS U-Pb zircon; Flowerdew et al. 2013), records island arc magmatism at this time. Our new data indicate that the initial phase of intra-oceanic arc magmatism began by at least c. 811 Ma and lasted until c. 685 Ma (the youngest sample in the ferroan, alkali-calcic group). This implies that magmatism associated with island arc (diorite-trondhjemite-tonalite) and subduction-related accretion (granite-granodiorite) genesis was contemporaneous in the region from c. 770 Ma, consistent with the overlap in timing and periods of magmatism suggested by Stoeser and Frost (2006) (Table 1). In addition, these data together confirm that the two generations of intrusions recognized by Fairer (1985), as well as the post-collisional intrusions, occur in both the Tayyah and the Khadrah belts, suggesting that the Junaynah fault zone is a within-terrane structure, rather than a suture between terranes (see discussion in Johnson and Kattan 2012).

Tectonomagmatic evolution of the ANS
In combining the geochemical, geochronological, and isotopic data, we establish an improved petrogenetic evolution. Relative to Nabitah orogenesis, three coherent tectonomagmatic phases are defined  Pallister et al. (1988).
in the southern Asir composite terrane on the basis of age and major element geochemistry (Table 2): 1) intra-oceanic island arc genesis at 811-685 Ma characterized by slightly ferroan, alkali-calcic, and peraluminous geochemistry; 2) subduction-related arc accretion (syntectonic phase) at 663-636 Ma characterized by magnesian geochemistry; and 3) a postcollisional phase at 625-613 Ma characterized by ferroan, alkali, and metaluminous geochemistry. None of our samples reflect the later generation of anorogenic granitic plutonism as defined by Robinson et al. (2014). These three stages fit within the established tectonic framework of the ANS (Bentor 1985;Stoeser and Frost 2006;Johnson et al. 2011;Robinson et al. 2014). Although the early island arc phase and the syn-tectonic phase are consistent with subductionrelated arc magmatism and the younger alkalic granites reflect post-collisional magmatism (Frost et al. 2001 and references therein), most samples have subduction-related trace element geochemistry (Figures 6, 7a) reflected by the decoupling of N-MORB-normalized LILE, and HFSE, with enriched LIL and Ta-Nb troughs (e.g. -Hawkesworth et al. 1993;Pearce and Peate 1995). Over time, the 'subduction zone' signature becomes more pronounced, with average (LREE/HREE) N increasing from (La/Lu) N = c. 12 for island arc phase to (La/Lu) N = c. 48 for accretion phase to (La/Lu) N = c. 60 for post-collision phase (Figure 7b; Table 2). This tripartite geochemical evolution is also seen in the subtle but discernible secular variation of εNd (i) to less juvenile (more radiogenic) compositions through time (Figure 7c; Table 2) -the intra-oceanic island arc magmatic phase is the most isotopically juvenile (εNd (i) >6) with Nd model ages (T DM ) that approximate their U-Pb zircon crystallization ages, the subduction-related arc accretion (syn-collisional) phase is less juvenile (εNd (i) c. 5-6) with Nd model ages that are slightly older than the U-Pb zircon crystallization ages, and the post-collisional phase has less juvenile/more radiogenic signatures (εNd (i) c. 4-5) with Nd model ages older than their associated U-Pb zircon crystallization ages ( Table 2).
The Tayyah Belt contains some of the oldest (> 815 Ma) signatures of the Arabian Shield (Johnson 2006). Robinson et al. (2018) obtained strongly evolved εHf (t) signatures in c. 990-810 Ma detrital zircons extracted from a mature, well-sorted metasandstone within the Tayyah Belt; because these data were difficult to reconcile with juvenile isotopic signatures of the local magmatic record within the Asir composite terrane, a distal source with continental affinity was inferred for the metasandstone. In contrast, our new data combined with existing geochemistry and isotopic data from the region show little indication of an evolved 'continental' source. Pb isotopes, which are sensitive indicators of crustal contamination (e.g. Stacey and Stoeser 1983), show no indication of continental input (Stoeser and Frost 2006). The juvenile εNd(t) of our samples, combined with the Pb isotopes of Stoeser and Frost (2006), does not support the significant involvement of continental sediments during the intra-oceanic island arc stage for the Asir composite terrane. The juvenile, albeit less radiogenic, isotopic compositions of the post-collisional plutons and their arc-like signatures suggest that their genesis involved recycling of earlier formed arc crust with minor (if any) continental sedimentary input. Therefore, we interpret the observed geochemical variation with time (increasing REE and radiogenic Nd, and the absence of crustal Pb) to be most consistent with the reworking (assimilation) of juvenile arc rocks. 6.2 6.9 6.6 5.9 5.6 4.9 5.8 3.8 5.1 Notes: 1, Fe* = FeO*/(FeO*+MgO); 2, Modified alkali-lime index (MALI) = Na2O + K2O -CaO; 3, Alumina saturation index (ASI) = Molecular ratios Al/(Ca -1.67P + Na +K); 4, NMORB normalization values of Sun and McDonough (1989), c.f.   Pearce et al. 1984) shows most samples represent arc and syncollisional granitoids (symbols as in Figure 4). B) REE geochemistry shows increasing REE from older to younger samples (N-MORB normalization values of Sun and McDonough 1989). Island arc samples in grey, syn-collisional samples in yellow, and post-collisional samples in green. C) Age (Ma) versus εNd (t) illustrating the juvenile nature of the Asir granitoids, but with a slightly increasing radiogenic component through time (depleted mantle, DM, after DePaolo 1981). Island arc samples = triangles, syn-collisional samples = squares, and post-collisional samples = circles (reference fields for the western island arc granitoids, WAG, and eastern arc granitoids, EAG, after Stoeser and Frost 2006). D) Age (Ma) versus crustal thickness (km) showing pronounced increase at c. 650 Ma from 20-30 km to c. 70 km (crustal thickness calculated using (La/Yb) N proxy of Profeta et al. 2015). Western and eastern arc terrane (WAT and EAT, respectively) databases for the northern ANS (black squares) from Robinson et al. (2015), Eyal et al. (2010), Eliwa et al. (2014) and Be'eri-Shlevin et al. (2011); the southern ANS database (grey squares) is from Robinson et al. (2015); circles represent data from this study. Associated errors (1 SD) typcially ± 6-10 km; grey bar represents the lower limit of the method (see Supplementary  Table 5).

Crustal thickening and post-collisional denudation
Our tripartite evolution may also reflect crustal thickening over time. Whole-rock trace element ratios have been successfully applied as proxies for the pressure/ depth of melting and used to determine crustal thickness in subduction-related magmatic arcs, e.g. , Chiaradia (2015), Chapman et al. (2015). These proxies exploit the differences in elemental partition coefficients in residual phases and intermediate melts as a function of pressure (a proxy for depth; see Profeta et al. 2015 and references therein). The ratio of light-to heavy-REEs is used to distinguish between deep and shallow fractionation processes as, for example, Y is incompatible at low pressures but readily partitions into amphibole or garnet with increasing pressures. Low La/ Yb ratios are expected in melts of subducted oceanic crust, reflecting the LREE-depleted nature of oceanic crust. Elevated La/Yb ratios are predicted in rocks involving fractionation from garnet and/or amphibole-rich residues in thick arcs. Thus, La/Yb ratios in arc rocks can be used as proxies to estimate the thickness of the crust.
Applying the method of Profeta et al. (2015) to dated samples for which appropriate geochemistry exists (restricted to Mg<4 wt% and intermediate SiO 2 = 55-68 wt%), we calculate crustal thickness associated with the ANS Western Arc Terranes from c. 850-600 Ma using the normalized element ratio (La/Yb) N from 50 arc granitoids ( Figure 7d). Excluding all ages <600 Ma, which is approximately the time of post-orogenic and anorogenic activity (Table 1) and despite the fact that a larger dataset would provide a more robust statistical regional average (beyond the time and length-scale of single intrusions), crustal thickness for the Western Arc Terranes shows a pronounced increase beginning at c. 650 Ma from 20-30 km to c. 70 km (Figure 7d). While additional pre-650 Ma data will refine this result, it nonetheless suggests an increase in crustal thickness at this time.
This increasing crustal thickness is consistent with the change in magmatic cycle from the early intra-ocean arc stage to the arc-accretion/syn-collisional stage (Figure 8). The overall effect of crustal thickening would likely result in a reduction of the melt column(s), in turn producing lower degrees of partial melting (e.g. Whattam and Stern 2016), and is consistent with the evolution of the chemical signatures documented here (Figure 7). Robinson et al. (2014) suggested that anorogenic magmatism at 590-550 Ma was due to lithospheric delamination and certainly if the WAT had a crustal thickness of c. 70 km, it is entirely plausible that the delamination of overthickened lithosphere could have occurred post-650 Ma.
The Atura formation also offers insights into the denudation processes and geodynamics that shaped the palaeogeography of the southwestern-most Arabian Shield. The Atura formation unconformity is dated to c. 645-625 Ma. Nabitah orogeny was active at c. 680-640 Ma and thickening of the Arabian Shield crust began at c. 650 Ma (Figure 8) as the result of collision between the accreted ANS terranes with the Saharan Metacraton to the west and the docking of exotic continental fragments to the east (Avigad and Gvirtzman, 2009). Intense denudation affected a large segment of the ANS shortly after the post-collisional plutonic phase, evidence for which is preserved across the northern Arabian Shield (Garfunkel 1999;Johnson 2003;Avigad and Gvirtzman, 2009;Johnson and Kattan 2012). Throughout the northern Arabian Shield, there are  Profeta et al. 2015). Crustal thickening, uplift, and denudation were likely episodic until regional Precambrian peneplanation occurred. Data for granitoids from the shield in gold (Stoeser and Frost 2006) and in brown (Robinson et al. 2015). Data for the southern Asir terrane granitoids in orange (this study). conglomerate and volcanics with unconformities pierced by post-tectonic plutons of c. 590 Ma (Stern and Hedge 1985;Jarrar et al. 1993;Garfunkel 1999;Yaseen et al. 2013). While geological evidence indicates that throughout the northern Arabian Shield denudation started prior to emplacement of post-collisional plutons (Figure 8), the major phase of denudation began immediately after the consolidation of post-collisional calcalkaline granitoids, c. 70 m.y. prior to formation of the prominent sub-Cambrian peneplain (Garfunkel 1999;Johnson 2003;Avigad and Gvirtzman, 2009). In the southern Arabian Shield, the Atura formation unconformity at c. 645-625 Ma developed during the waning stages of orogenesis, c. 80-100 m.y. prior to the formation of the sub-Cambrian peneplain, and provides additional evidence that multiple phases of denudation occurred before the termination of Nabitah orogenesis.

Conclusions
The granitoids in the Asir composite terrane fill a regional data gap and record about 200 m.y. of tectonomagmatic evolution. Island arc subduction-related plutonism occurred at c. 810 Ma to 685 Ma, with arcaccretion and syn-collisional plutonism following from c. 680-640 Ma. Post-collisional alkalic plutonism from c. 635-610 Ma marks the end of orogenesis in the region, as well as the termination of island arc terrane accretion phase at c. 635 Ma. This implies that the closure of the Mozambique Ocean associated with the assembly of Gondwana occurred by c. 640 Ma in the region.
Geochemical variation in magmatism (increasing REE, less radiogenic Nd, juvenile Pb) with time likely reflects minor contamination/assimilation associated with thickening juvenile arc-related crust. The accreted crust reached its maximum thickness of c. 70 km at c. 650 Ma. By the time of post-collisional plutonism sometime after c. 640 Ma, arc accretion ceased and multiple pulses of uplift/denudation occurred during the waning stages of orogenesis.
The Atura formation unconformity constrains one denudation pulse to 645-625 Ma. The unconformity records erosion of the uplifted orogenic edifice of the southern Arabian Shield, consistent with a buoyant, thickened crust, and indicates that denudation in the Arabian-Nubian Shield began 80-100 m.y. before the sub-Cambrian peneplanation of the Gondwana supercontinent.