High-precision zircon geochronology and geochemistry of evolved magmatic centres in the Paraná-Etendeka LIP: Temporal placement and tectono-magmatic origin of the Damaraland complexes, Namibia

In the African section of the Paran ´ a-Etendeka Large Igneous Province (LIP), several early Cretaceous intrusive complexes with highly variable bulk rock compositions ranging from subalkaline to alkaline ‑ carbonatitic are strikingly aligned on a SW-NE trend and are clearly related to the large igneous province that straddles Africa and South America based on their 40 Ar/ 39 Ar and Rb – Sr ages. We use high-precision zircon U – Pb geochronology, zircon Hf isotope compositions, and bulk-rock and zircon geochemistry of these intrusions to decipher the spatial pattern of magmatism and the mechanisms of magma formation and evolution. High-precision CA-ID-TIMS zircon ages show that the plutons were emplaced without clear spatial systematics over at least 1.4 Myr from 133.2 to 131.8 Ma and likely postdate the peak of volcanic activity in the southern Etendeka (Namibia) and the Paran ´ a (Brazil), estimated to be between 135 and 133 Ma. The youngest Damaraland zircon ages overlap with previously reported baddeleyite ages for mafic dykes and sills in the southern Paran ´ a, confirming that mag-matism in the whole LIP continued beyond 132 Ma. In the whole Damaraland province, intrusive rocks follow two distinct evolution paths towards both the granite/rhyolite and the syenite/phonolite minima. The near-synchronous emplacement of silicic to highly alkaline bodies with no clear spatial or temporal pattern is best explained by lithospheric extension driven by mantle upwelling along reactivated lineaments accommodating active N-S rifting of the South Atlantic. Negative zircon ε Hf(t) values in silicic complexes ( (cid:0) 5 to (cid:0) 12) likely stem from crustal contamination during ascent of primitive mafic melts via “ path contamination ” and/or from a heterogenous mantle source already contaminated by recycled crustal material ( “ source contamination ” ). Similarly variable, but on average less negative zircon ε Hf(t) in alkaline and carbonatitic intrusions (0 to (cid:0) 6) reflect source contamination and preferential low degree melting of fertile portions in a heterogenous mantle source. These mantle heterogeneities likely originate from recycling of ancient lithosphere during subduction in the Neoproterozoic, or from fertile material brought up by the Tristan-Gough mantle plume via intrusion and/or exposed by lithospheric erosion, fertilizing the lithospheric mantle under the Damaraland prior to and during the main flood basalt event.


Introduction
Large Igneous provinces (LIPs) are produced during some of the most extensive magmatic events occurring on Earth, emplacing igneous volumes of over 750,000 km 3 within very short timescales of 1-5 Myr, and are often directly related to global extinction events and climatic turmoil episodes, the break-up of continents, and the formation of important economic resources (e.g., Bryan and Ernst, 2008;Burgess et al., 2017;Pirajno and Hoatson, 2012).Despite their immense impact for the Earth's surface, the formation processes and evolution of many LIPs are still under debate, especially under the light of the general question, whether the thermomechanical effects of mantle plume impingement, or lithospheric processes such as crustal extension and thinning, are the main driving forces behind the emplacement of LIPs and the initiation of continental rifting (Matton and Jébrak, 2009;Niu, 2020;Sheth, 1999).A major challenge is that LIP magmatism may develop over very short timescales which are only resolvable by highest-precision geochronology (Kasbohm et al., 2021).
The Early Cretaceous Paraná-Etendeka Large Igneous Province LIP, one of the largest continental LIPs known, is related to the break-up of Gondwana and opening of the South Atlantic Ocean, and it is suspected to have significantly contributed to Early Cretaceous climatic changes such as the Weissert anoxic event (Cavalheiro et al., 2021).While the major magmatic volume of this bimodal LIP was emplaced within the Paraná basin in South America as mafic (tholeiites to basaltic andesites) to silicic (dacite to rhyolites) lava piles, the volumetrically smaller Etendeka province in southwestern Africa also contains a series of remarkable, well-preserved and -exposed plutonic centres in addition to the mafic to silicic volcanic rocks (Harris, 1995;Trumbull et al., 2000).These intrusions range in composition from granites to carbonatites, have similar 40 Ar/ 39 Ar ages (see supplemental table S1 for compilation of literature ages) to the nearby Etendeka volcanic rocks and are collectively known as the Damaraland igneous province or Damaraland complexes (e.g., Harris, 1995;Milner et al., 1995;Trumbull et al., 2000).The tectono-magmatic origin of the Paraná-Etendeka continental flood basalts and the associated Damaraland intrusions remains controversial.In particular, the Damaraland complexes have been interpreted as a continental hotspot track of the Tristan-Gough plume based on their striking NW-SE trending alignment, their close spatial-temporal proximity to the age-progressive oceanic seamount ridges in the Southern Atlantic (Fig. 1A), and the often near bulk-silicate earth/asthenospheric mantle Sr-Nd-Pb isotope compositions (e.g., Trumbull et al., 2003 and references therein).However, a systematic spatial age progression between individual complexes has not yet been demonstrated (Harris, 1995).
Reliable radiometric age information is not only crucial for the tectono-magmatic interpretation of Damaraland magmatism, but also for assessing the potential environmental impact of the Paraná-Etendeka LIP magmatism on Early Cretaceous life by providing correlation points and estimations of magmatic durations and fluxes (Bryan and Ernst, 2008).For the Etendeka region, previous attempts to unravel its temporal evolution and to correlate magmatic events with the Paraná region have been hampered by the lack of a reliable geochronological framework covering the whole LIP, as earlier 40 Ar/ 39 Ar ages of hornblende, plagioclase, and mica were often insufficiently precise and accurate to establish the exact temporal relationships within the LIP (Gomes and Vasconcelos, 2021).Here, we present the first high-precision chemical abrasion-isotope dilution-thermal ionization mass spectrometry (CA-ID-TIMS) U-Pb zircon ages from the African side of the Paraná-Etendeka LIP and use them to explore the timing and the emplacement mechanisms of the Damaraland igneous suite.Bulk rock and zircon geochemistry, and zircon Hf isotope compositions aid in tracking the magmatic origins and evolution paths of the melts.

Geological setting
The Damaraland igneous province is located in northern to central Namibia (Fig. 1A-B) and comprises at least 13 Early Cretaceous intrusive centres which are exposed to the south of the main Etendeka lava field (e.g., Harris, 1995;Trumbull et al., 2004).Based on the lithologies they contain, the complexes can be divided into three main types: 1) subalkaline complexes consisting mainly of granites, 2) bimodal complexes containing both subalkaline (e.g., granite) and alkaline or mafic components (e.g., nepheline syenite, gabbro), and 3) alkaline complexes consisting of silica-undersaturated rock types such as nepheline syenite, often including carbonatites (Harris, 1995).Silica-oversaturated granitoids are the most common rock types in the Damaraland complexes (ca.84% of outcrop area; Fig. 2).Mafic rocks (~9%) and syenites (~7%) are less abundant, and carbonatites are only locally present and volumetrically negligible (Fig. 2).The Damaraland complexes all intruded along the Damara orogenic belt, which formed in the Pan-African (580-500 Ma) orogeny by collision of the Precambrian Congo and Kalahari cratons (Trumbull et al., 2004).Country rocks comprise folded metasedimentary sequences intruded by syn-and post-collisional granite plutons (e.g., Trumbull et al., 2000).All Damaraland complexes follow a 350 km long, 120 km wide zone parallel to the SW-NE strike direction of major lineaments and the ancient continental margins of the Congo and Kalahari cratons at the time of collision (e.g., Raab et al., 2002).From N to S, four major lineaments are present (Fig. 1B): the Khorixas-Gaseneirob thrust, the Autseib Fault-Otjohorongo Thrust, the Omaruru Lineament-Waterberg Thrust, and the Okahandja Lineament (Clemson et al., 1997;Kramm et al., 2017).Cross-cutting relationships and thermochronological data show that these major fault zones were repeatedly reactivated throughout the Phanerozoic, including during the Cretaceous (Clemson et al., 1997;Raab et al., 2002).

Bulk rock geochemistry
Igneous rocks including granite, syenite and rhyolite/comendite were sampled from the major rock units of 11 different Damaraland complexes: Cape Cross, Messum, Brandberg, Spitzkoppe, Erongo, Okenyenya, Otjihorongo, Ondurakorume, Etaneno, Paresis, Okorusu (Fig. 1B; sample list in table S2 in the supplementary material).Bulk rock preparation and analyses were carried out at the Department of Earth Sciences, ETH Zurich.Representative sample material (ca.150-300 g) was freed from alteration crusts and veins with a blade saw and crushed in a hydraulic press down to <4 mm grain size.An aliquot of ~50 g per sample was milled into fine powder using a tungsten carbide vibratory disc mill.Loss on ignition (LOI) was determined on dried rock powder (110 • C, overnight) by weighing in a ceramic crucible before and after firing (2 h, 900 • C).Fused glass beads were prepared by mixing ~1.6 g of sample powder with lithium tetraborate/lithium metaborate (66%/34% mix) fluxing agent in a ratio of 1:5.The samplefluxer mix was homogenised for 4 min in an agate mortar and cast into glass disks using a PANalytical Claisse Eagon 2 automated fusion instrument.Measurements of major and minor elements in the prepared glass beads were done by X-ray fluorescence spectrometry (XRF) using a PANanalytical AXIOS WD-XRF spectrometer.The glasses were shattered afterwards, and trace elements were measured at ETH Zurich on the shards by LA-ICP-MS with a Coherent ComPex Pro 102F 193 nm ArF laser system coupled to a PerkinElmer NexION 2000 quadrupole mass spectrometer.Measurements were carried out with a spot size of 115 μm, a laser repetition rate of 10 Hz and on-sample fluence of 8-10 Jcm − 2 .The carrier gas consisted of high-purity He and Ar, which were both added with a flow rate of 1.1 lmin − 1 .The measurement routine consisted of 30 s of gas blank acquisition and ~40 s of sample ablation, following standard primary reference material (RM) -flux blanksecondary RMsample bracketing.As the primary RM, NIST SRM610 glass was measured at a spot size of 40 μm.Flux blank and the secondary RM BCR-2 glass were measured at 115 μm spot size to check the accuracy and reproducibility of the measurements.Reduction of the mass spectrometry data was done using the MATLAB-based program SILLS (Guillong et al., 2008).For most samples, the CaO content measured by XRF was used as the internal standard, and, where CaO concentrations are below 1 wt.-%,Al 2 O 3 was used instead.

Zircon U-Pb geochronology and trace elements by LA-ICP-MS
Zircon and baddeleyite were separated from the rock samples by standard techniques including fragmentation with the Selfrag device, heavy mineral separation using a Wilfley shaking table, and hand picking under the binocular microscope at the Department of Earth Sciences, ETH Zurich.All zircon and baddeleyite grains were thermally annealed at 900 • C for 48 h prior to analysis to recrystallize radiationdamaged domains and were embedded in epoxy mounts using Struers   (Rohde et al., 2013) and the present-day locations of the Tristan da Cunha and Gough islands.(B) Locations of sampling sites, major Precambrian lineament zones, and important Early Cretaceous volcanic outcrops in the Damaraland province.Each complex is marked with an individual symbol with symbols in red denoting exclusively subalkaline complexes, symbols in green representing bimodal complexes, and symbols in blue referring to exclusively alkaline complexes.(C) Locations of the alkaline complexes of the Kalkfeld group with the emplacement ages shown for the older generation of intrusions.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)EpoFix resin.Mounts were ground with SiC polishing paper for maximum surface exposure, polished using diamond suspension (0.25 μm grain size) and coated with a 20 nm conductive carbon coating.
Cathodoluminescence (CL) images were taken at ETH Zurich using a JEOL JSM-6390 LA scanning electron microscope (SEM) equipped with a Deben Centaurus panchromatic CL detector to identify strongly metamict zircon domains, inclusions or other inhomogeneities.
Uranium-Pb ages and trace element compositions of thermally annealed zircon and baddeleyite grains were measured by laser ablation inductively coupled plasma mass spectrometry (LA-ICP-MS).Analyses were carried out using an ASI (Resonetics) RESolution S-155 193 nm ArF excimer laser ablation system equipped with a Laurin Technic dualvolume ablation cell at a laser repetition rate of 4 Hz, on-sample fluence of ~2 Jcm − 2 , and a spot size of ~19 μm.Each spot was ablated for 30 s during measurement.In the cell, 100% He was used as carrier gas for the sample aerosol, and Ar was added as make-up gas with a flow rate of 0.90-0.99lmin − 1 .For the U-Th-Pb measurements GJ-1 zircon was used as primary RM, and 91500, AUSZ7-1, and Plesovice zircon served as secondary RM.Standardization of trace element measurements was done on NIST 610 glass using Si for internal standardization, and four natural zircon RM (91500, GJ-1, Plesovice, Temora-II) were measured to vet the accuracy of the analyses.Data reduction was done in the Iolite 4 program (Paton et al., 2011), using the VizualAge add-on (Petrus and Kamber, 2012) and down-hole fractionation corrections following Paton et al. (2010).For the baddeleyite measurements, primary standardization of U-Pb isotope measurements was done with Phalaborwa baddeleyite, and the secondary RM used were FC-1 and Mogok baddeleyite.All ages are reported with uncertainties stated as 2σ.
Prior to the calculation of weighted mean ages, the geochronological dataset was filtered using predefined cut-off criteria: discordance <20%, calculated from 207 Pb/ 235 U vs 206 Pb/ 238 U ratios, equivalency of 207 Pb/ 235 U and 206 Pb/ 238 U ages within 2σ uncertainties, 35 ppm Al as proxy for silicate/melt inclusions, 2000 ppm P as proxy for apatite inclusions, and high LREE indicative of a hydrothermal origin of the zircon.For Erongo zircon, an Al content of 100 ppm was used as exclusion criterion, as the sample show higher Al contents in general, and a cut-off at 35 ppm would have excluded most measured grains.Concordia ages and weighted averages were then calculated using Iso-plotR (Vermeesch, 2018).

High-precision zircon geochronology by CA-ID-TIMS
Between 5 and 8 concordant zircon grains previously analysed by LA-ICP-MS were selected from seven complexes (Cape Cross, Messum, Brandberg, Okenyenya, Etaneno, Paresis, Okorusu) for high-precision dating by chemical abrasion isotope dilution thermal ionization mass spectrometry (CA-ID-TIMS) based on their lack of inheritance, their homogeneity under BSE and CL imaging, and grain size.
Zircon from Spitzkoppe, Erongo and Otjihorongo show textural features typical for U-rich zircons crystallized from highly differentiated hydrothermal fluids, such as high degree of radiation damage and dissolution-reprecipitation-related textures with abundant microinclusions and were therefore not dated by CA-ID-TIMS (see CL images in the supplementary material figures Fig. S1-S13).The preliminary U-Pb obtained by LA-ICP-MS for Ondurakorume revealed that the complex is not related to the ~132 Ma Early Cretaceous magmatism, but to an older, likely Permian phase of magmatism, and its temporal difference to the other Damaraland complexes is readily resolved at the analytical uncertainties of the LA-ICP-MS measurements.Hence, the Ondurakorume complex was also not dated by CA-ID-TIMS.
Sample preparation and TIMS analyses were done at ETH Zurich.Zircon grains were extracted from the epoxy mounts and chemically abraded (Mattinson, 2005) for 12 h at 190 • C in 200 μl Savillex PFA microcapsules inside a steel pressure vessel using concentrated HF.Afterwards, the zircon grains underwent two ultrasonic cleaning steps in 3 ml Savillex PFA hex beakers with 6 N HCl and 7.5 N HNO 3 .Cleaned zircon grains were loaded into 200 μl microcapsules with 6-15 mg of the EARTHTIME 202 Pb-205 Pb-233 U-235 U tracer (Condon et al., 2015;McLean et al., 2015) and dissolved in 40 μl of concentrated HF in a steel bomb at 210 • C (60 h).The dissolved samples were evaporated on a hot plate (140 • C) and the fluoride species were converted into chloride by redissolving the sample in 6 N HCl at 190 • C (12 h) in the steel bomb.After conversion, the samples were evaporated again and loaded in 3 N HCl into 50 μl columns filled with Bio-Rad™ AG1-X8 resin.Matrix elements were washed out with 3 N HCl before eluting Pb and U into the sample beaker with 6 N HCl and water, respectively.The purified Pb and U solution was mixed with a micro-drop of H 3 PO 4 and dried down on the hotplate (140 • C).
Measurements of U and Pb isotopes were conducted with a Thermo Scientific Triton Plus TIMS.Samples were loaded onto outgassed, zonerefined Re filaments using 1 drop of Si-gel emitter to stimulate ionization of U and Pb.Isotope ratios of Pb were measured in static mode using Faraday cups equipped with 10 13 Ω resistors (von Quadt et al., 2016).The following Pb species were detected: 202 Pb (spike), 205 Pb (spike), 206 Pb, 207 Pb, 208 Pb.Isotope 204 Pb was detected on the axial ion counter due to the expected low abundance of this isotope.Mass fractionation was corrected online using the known 202 Pb/ 205 Pb of the spike.
Uranium was measured as oxide species (UO 2 ) using a static mode   et al., 2017).Oxide interferences were corrected online using 18 O/ 16 O derived from the analysis of mass 272 whenever permitted by beam size, or alternatively using a fixed 18 O/ 16 O = 0.002064 ± 20.Correction for the instrumental mass fractionation was done offline using the double spike composition and assuming sample 238 U/ 235 U of Hiess et al. (2012).Data reduction and age calculations were done with TRIPOLI and ET redux (Bowring et al., 2011;McLean et al., 2011) using decay constants of Jaffey et al. (1971), and all ages were corrected for initial Th/U disequilibrium using a fixed Th/U magma ratio of 3.5 as best estimate.All ages are reported with 2σ uncertainties.Further information about the data reduction and age calculations are available in the supplementary table S6 in the appendix.

Zircon Hf isotopes by LA-MC-ICP-MS
The remaining zircon grains with ages between 140 and 124 Ma (LA-ICP-MS) were ablated for Hf isotope measurements by Multi-Collector Inductively Coupled Plasma Mass Spectrometry (MC-ICP-MS) at ETH Zurich with an ASI (Resonetics) RESolution S-155 193 nm ArF excimer laser ablation system connected to a Nu plasma 2 multicollector sectorfield mass-spectrometer.A laser spot size of ~50 μm, a repetition rate of 5 Hz and an on-sample fluence of ~4 Jcm − 2 were chosen during data acquisition.High purity He and Ar were used as carrier gas.Data were reduced offline using the Iolite software and following an in-house data reduction scheme.Instrumental mass bias for Hf and Yb isotopes were corrected using natural abundance ratios for 179 Hf/ 177 Hf and 173 Yb/ 171 Yb of Chu et al. (2002), and radiogenic 176 Hf/ 177 Hf isotope ratios were corrected for isobaric interferences by 176 Yb and 176 Lu with the mass-bias corrected ratios of 176 Yb/ 173 Yb = 0.796218 and 176 Lu/ 175 Lu = 0.02656.The effects of the isobaric interferences of 176 Yb on 176 Hf was empirically monitored and corrected for 176 Yb/ 177 Hf ratios up to 0.1 using the measured values of reference materials GHR-1 and Temora 2. Spot measurements affected by high Yb ( 176 Yb/ 177 Hf higher than 0.1) were discarded.The following zircon reference materials were measured for controlling the accuracy and external reproducibility of the method: Mud Tank, GHR-1, RAK17, 91500, GJ-1, Plesovice and Temora 2. Initial Hf isotopic ratios εHf(t) were calculated using crystallization ages for the same spots previously measured by LA-ICP-MS.

Zircon U-Pb geochronology
Zircon grains are 100-300 μm in size and display oscillatory or sector zoning (Fig. S1-S13 in the supplementary material).For each sample, between 11 and 67 concordant zircon/baddeleyite analyses were obtained by LA-ICP-MS.All complexes gave Early Cretaceous (ca.134-130 Ma) ages (Fig. 3; Fig. S14-S15 in the supplementary material), except for an Ondurakorume syenite, which shows two zircon age populations dated at 275.3 ± 7.9 Ma and 510.3 ± 5.1 Ma (Fig. 3D, F).While the older population are most likely inherited zircon grains scavenged from the Pan-African basement, we interpret the younger population around 275 Ma as the magmatic emplacement age of the complex.Ondurakorume is thus unrelated to the Early Cretaceous Paraná-Etendeka LIP and the other Damaraland complexes but formed during an earlier episode of magmatism, similar to the nearby Kalkfeld alkaline complex, which was dated to 242 ± 6.5 Ma using Rb-Sr isochron dating (Kramm et al., 2017).The two youngest ages for Cape Cross were detected as statistical outliers using the Chauvenet criterion and excluded from the weighted average age, as well as the youngest grain from the Amis satellite intrusion at Brandberg.Clearly inherited zircon grains were absent in most samples, except for the Erongo granite which shows three inherited grains (430 Ma, 440 Ma, 1.07 Ga) out of 32, and the syenites from Etaneno (3/70) and Okorusu (1/25).The youngest grain from Okenyenya was rejected as an outlier by IsoplotR.For Otjihorongo, the youngest grain was excluded from the mean age, as well as two older grains at ~150-160 Ma, which we interpret as mixed ages with inherited zircon parts.
Overlapping LA-ICP-MS dates for baddeleyite and zircon from the Okorusu syenite (131.5 ± 1.7 Ma and 131.3 ± 1.4 Ma, respectively) suggest that Okorusu zircon grains were not trapped in the alkaline melt as xenocrysts, as baddeleyite gives primary crystallization ages even in strongly silica undersaturated magmas (e.g., Heaman and LeCheminant, 1993).Secondary zircon rims over baddeleyite indicative of later assimilation processes or magma mixing with more silicic melts, which could have provided zircon xenocrysts, were not observed under the SEM (see Fig. S13 in the supplementary material).
Pristine zircon crystals devoid of inclusions and inherited cores were used for CA-ID-TIMS dating.In total, 39 high-precision single-crystal U-Pb dates were obtained.Single-crystal 206 Pb/ 238 U dates range from 133.282 ± 0.060 to 131.776 ± 0.062 Ma but form a single age population in most individual samples.Due to the small degree of observed in-sample age dispersion and because zircon crystals were generally extracted from the most voluminous units in the complexes, weighted average ages were calculated for all complexes and interpreted as best estimates for the emplacement of the complexes.For Cape Cross, an abnormally young grain at 128 Ma was attributed to unmitigated Pb loss and rejected from the weighted average.For the Okenyenya complex, a single unusually old grain which was interpreted as an antecryst was excluded from the average value.For Okorusu, the three oldest grains were excluded from the weighted average.
Weighted mean intrusion ages range from 133.2 to 131.8 Ma (Fig. 4) and do not exhibit any spatial trend from the NE (Okorusu) to the SW (Cape Cross).Instead, the ages testify to a short period of magmatism across the entire province over ~1.4 Myr.The Etaneno syenite, Paresis comendite, Brandberg granite and Cape Cross granite (weighted mean ages of 133.219 ± 0.017, 132.963 ± 0.017, 132.875 ± 0.075 and 132.858 ± 0.052 Ma, respectively) show the oldest emplacement ages for individual intrusive phases, while the Messum quartz syenite, Amis granite (peralkaline intrusion within Brandberg) and Okenyenya syenite show the youngest (132.017± 0.022, 131.879 ± 0.020 and 131.811 ± 0.023 Ma, respectively).

Bulk rock and zircon geochemistry
In both the W and E part of the province, the Damaraland complexes follow two different melt evolution paths (Fig. 5): the Etaneno (Central) and Okorusu (E) complexes evolve only towards the phonolite minimum.The Messum (W) and Okenyenya (Central) complexes mostly evolve towards the phonolite minimum as well, but locally granitic rocks are formed.Rocks of the Cape Cross (E) and Paresis (W) complexes mostly progress towards the granite minimum, but locally also towards the phonolite minimum.All other complexes evolved only towards the granite minimum.Therefore, no clear systematic variation in magma types from E to W is recognizable.
Damaraland zircons show variations in Hf concentration and isotope composition which most likely depend on differences in host lithology and/or magma source (Fig. 6).In the subalkaline and bimodal complexes, zircon from the Cape Cross granite, Messum quartz syenite, Brandberg granite, Spitzkoppe granite, Erongo granite, Okenyenya syenite and Otjihorongo granite show overlapping Zr/Hf of 30-70.(Vermeesch, 2018).Error ellipses are plotted as 2σ, and reported MSWD and chi-square probabilities were calculated in IsoplotR for combined concordance and equivalence.

Timing and duration of Damaraland magmatism
Our new high-precision zircon ages allow us to place the Damaraland complexes in the frame of the Paraná-Etendeka LIP evolution (Fig. 7).Emplacement ages of the oldest intrusions (Cape Cross granite, Brandberg granite, Etaneno syenite and Paresis comendite) between 133.2 and 132.8 Ma show that the Damaraland magmatism dominantly occurred immediately after the 135-133 Ma volcanic activity peak in the southern Etendeka region and in the southern to central Paraná region, Brazil ( 40 Ar/ 39 Ar data, Gomes and Vasconcelos, 2021).
The oldest emplacement ages for the Damaraland intrusions likely postdate the 134.7 ± 0.3 to 133.9 ± 0.7 Ma zircon/baddeleyite U-Pb CA-ID-TIMS ages recorded for mafic dykes of the Florianópolis dyke swarm, which are interpreted as the magmatic plumbing system of the high-Ti mafic lavas in the southern Paraná (Urubici lava) and northern Etendeka (Khumib lava) provinces (Florisbal et al., 2014).Similar CA-ID-TIMS zircon ages to the Florianópolis dyke swarm were reported for alkaline and mafic, tholeiitic intrusions in the Ponta Grossa Arch, Central Paraná, which mostly crystallized from 133.9 to 133.4 Ma, with a single mafic intrusion yielding an older best-fit zircon age of 134.5 ± 0.1 Ma (Almeida et al., 2018).
The stratigraphically oldest, early erupted low-Ti silicic sequence (Palmas-type dacite to rhyolite) in the southern Paraná dated at 133.6 Ma by CA-ID-TIMS on zircon is slightly older, or synchronous with the oldest Damaraland intrusions (Rocha et al., 2020).The oldest Damaraland intrusions are also indistinguishable in age from the 132.9Ma eruption age reported for the late-erupted high-Ti silicic lava (Chapecótype trachydacite) in the central Paraná basin (Rocha et al., 2020).Our zircon crystallization ages for the younger Damaraland intrusions (Messum quartz syenite, Amis granite, Okenyenya and Okorusu syenites) at 132.2 to 131.8 Ma clearly postdate the published zircon ages for all intrusive and extrusive magma types of the Paraná province.They also postdate or overlap the young baddeleyite TIMS ages reported for tholeiitic high-Ti-Sr dykes and dolerite intrusions (Taió sill) belonging to the Florianópolis dyke swarm in southern Paraná, Brazil (132.53 ± 0.40 and 132.07 ± 0.27 Ma; Fig. 4), suggesting that these young baddeleyite ages may be geologically meaningful, rather than being caused by secondary Pb loss as previously suspected (Rocha et al., 2023).
The cessation of Paraná-Etendeka LIP magmatism was previously constrained based on eruption ages of the youngest volcanic units in the Paraná lava pile, a 132.906 ± 0.072 Ma silicic high-Ti rock (Guarapuava trachydacite) in the south and a 133.2 ± 0.3 Ma basalt (Paranapanema high-Ti basalt) in the central/northern section (Gomes and Vasconcelos, 2021 and references therein).Our new data indicate that all complexes in the ~350 km Damaraland igneous province intruded within a shortlived (~1.4 Myr) pulse during the late evolution stage of the Paraná-Etendeka magmatism and after the emplacement of the main lava pile in Paraná and, most likely, the Etendeka.Although the 1.4 Myr duration for the Damaraland intrusive phase is much shorter than the previous estimates of 3-5 Myr based on 40 Ar/ 39 Ar dating (e.g., Renne et al., 1996;Zhou et al., 2020), the new ages are significantly younger than previously reported U-Pb zircon ages from the Paraná basin and confirm that major magmatic activities in the Paraná-Etendeka LIP lasted at least until ~132 Ma.The relatively short lifetime of the Damaraland province magmatism is comparable to the <2 Myr timescales observed for evolved intrusive complexes in other continental LIPs, including the evolved intrusions associated with the Deccan Traps, India (Basu et al., 2020;Schoene et al., 2015) and the alkaline layered intrusions of the Siberian Traps, Russia (Augland et al., 2019).

Magma source and the role of crustal contamination
The interplay of isotopic and geochemical information allows us to speculate on the type(s) of mantle that experienced partial melting in the Early Cretaceous to generate the parental melts of the Damaraland complexes.Zircon grains from the Damaraland province show variable Hf isotope compositions ( 176 Hf/ 177 Hf at 130 Ma from 0.28237 to 0.28286), which plot around or below the Hf isotope evolution line for the chondritic uniform reservoir (CHUR) at 130 Ma, and the negative zircon εHf(t) values observed from several of the units (Fig. 8) indicate that some of the complexes interacted with old crustal materials during their petrogenesis.Previous experimental results and geological observations suggest that A-type granitoids are likely not produced by melting of crustal material, but instead originate from magmatic differentiation of mantle-sourced primary liquids (Bonin, 2007; see also Harker plots in Fig. S16-S17 in the supplementary material), which can attain crust-like isotope compositions only through crustal contamination either in the mantle source by metasomatism ("source contamination"), or by assimilation during magma ascent ("path contamination"; see discussion in Cornet et al., 2022).In the latter case, any relevant assimilation of quartz-bearing crustal rocks or crustally-derived melts would shift the melts to silica-richer compositions and therefore force the magma evolution towards Si-rich side of the feldspar thermal divide and to the granite minimum, even if the original primary melt was strongly silicaundersaturated (Fig. 5).Conversely, it is unlikely that melts that evolved towards the phonolite minimum suffered any significant path  (Rocha et al., 2020;Rocha et al., 2023).Dates shown in grey are excluded from the calculation of the weighted mean age.

contamination.
The two Damaraland complexes which exclusively evolve towards the phonolite minimum (Etaneno, Okorusu) display clearly negative εHf (t) values of − 3 to − 7 (Fig. 6).Highly variable zircon εHf(t) of +6 (Messum) to − 5 (Okenyenya) are shown by syenites in predominantly alkaline complexes (Messum), which are silica-oversaturated but nevertheless show absence of significant crustal assimilation during melt ascent and evolution as demonstrated by their whole rock oxygen isotope data, and by the narrow clustering of Sr-Nd-Pb isotopic compositions near bulk-earth values (Mingram et al., 2000;Trumbull et al., 2000;Trumbull et al., 2003;Trumbull et al., 2004 and references therein;Fig. 8).This suggests that the primary melts forming the alkaline complexes arose from an isotopically heterogeneous mantle source, the negative εHf(t) values testifying to heterogeneous crustal-related sources within the lithospheric mantle.The evolution of Hf isotope compositions for global carbonatites and carbonatite-hosted zircon/ baddeleyite crystals over the past 3 Gyr shows that variations towards more crust-like values as observed in the Damaraland region became more common over the last 1 Gyr, which likely reflects the increasing contamination of mantle source regions by recycled crustal materials following the establishment of colder and deeper subduction mechanisms during the Late Proterozoic (Yaxley et al., 2022 and references therein).Similar conclusions have also been drawn from mafic to silicic rocks in convergent plate margins, e.g., the Famatinian arc in Argentina, where increasing zircon O-Hf isotope variations over geologic time show an increasing incorporation of subducted crustal materials into the mantle source regions over the past 2.5 Gyr (Cornet et al., 2022).In the Paraná basin, the role of mantle source contamination in the petrogenesis of volcanic and plutonic rocks has been suspected based on highly variable Sr-Nd isotope compositions of alkaline and ultramafic rocks of the Jacupiranga complex (Ponta Grossa Arch Alkaline Province, SE Paraná), including lithologies unaffected by crustal contamination such  Bowden et al., 1990;Haapala et al., 2007;Harris et al., 1999;Miller, 2008;Mingram et al., 2000.as clinopyroxenite and carbonatite, were interpreted as being inherited from a heterogenous lithospheric mantle source (Chmyz et al., 2017).The possibility of mantle source contamination was also discussed for mafic to intermediate lavas from the Paraná, in which variable and negative zircon εHf(t) of − 1 to − 8 were recorded, although, in that case, the role of crustal contamination in producing this crust-like isotopic signature could not be excluded (Hartmann et al., 2019).
In contrast to the alkaline complexes, the negative εHf(t) values reported from the intrusions which evolved towards the granite minimum can result from either source or path contamination or a combination of 87 Sr/ 86 Sr (132 Ma) ratios of 0.713-0.760(Trumbull et al., 2000 and references therein).The crustally contaminated Damaraland plutonic rocks either fall on a mixing line between pristine mantle values and these Panafrican granites/metasediments, or, in extreme cases like the highly contaminated Erongo granite, overlap in Sr-Nd isotope compositions with the basement lithologies (Fig. 8).
The comendite from the bimodal Paresis complex, which evolved both towards the granite and the phonolite minimum, shows highly negative zircon εHf(t) around − 20.Due to the isobaric interferences of 176 Yb on 176 Hf, which cannot be reliably corrected for the Paresis zircon grains due to their extremely high Yb content (up to 2320 ppm), these values are considered unreliable and therefore were omitted from the plot (Fig. 6).The crust-like Sr-Nd isotope compositions of Paresis rocks published in the literature (Fig. 8), however, indicate that these anomalously Hf values do not entirely stem from the 176 Yb interference, but it is very likely that significant path contamination during magma evolution also contributed to these negative values.The extremely non-radiogenic Nd ratios of εNd = − 21 and radiogenic 87 Sr/ 86 Sr i = 0.7117-0.7138are consistent with a genesis of the silica-oversaturated Paresis rhyolite and comendite by extensive path assimilation of Precambrian crust, most likely gneisses originating from the nearby Congo craton, from originally strongly undersaturated primary melts (Mingram et al., 2000).

Possible origin of mantle heterogeneities in the Damaraland and the role of the Tristan-Gough plume
The involvement of mantle heterogeneities in the petrogenesis of the Damaraland complexes raises the question which geologic events and mechanisms are responsible for inducing fertile crustal materials into the mantle regions below the Damaraland crust.The Damara mobile belt experienced a long period of relative tectonic stability after its assembly in the Late Neoproterozoic to Cambrian Pan-African orogeny (e.g., Kramm et al., 2017;Trumbull et al., 2000).Therefore, one possibility is that mantle fertilization far predated the Early Cretaceous magmatic activity and occurred already during the assembly of the belt in the Pan-African orogeny by the subduction of Precambrian crust, which is further indicated by the early formation of the alkaline complexes Kalkfeld and Ondurakorume in the Damaraland during the Permian to Triassic.
Another possibility is that the Tristan-Gough plume was responsible for the crust-derived signature in the melt source region, and plumerelated enrichment contributed to modify the lithospheric mantle source of the melts (Pilet et al., 2008), imprinting the crustal Hf isotope signatures into it and producing at least part of the observed mantle heterogeneities under the Damaraland crust.The presence of high 3 He/ 4 He (>24 R A ) olivine crystals in mafic dykes of the Etendeka province is an indicator for the involvement of deep mantle-derived plume material in the generation of melts during the Etendeka magmatism (Stroncik et al., 2017).Presence of the Tristan-Gough plume head under southwestern Africa during the Early Cretaceous is supported by stratigraphy and thermochronology studies, which show surface uplift in the Damaraland prior to and after the Paraná-Etendeka flood basalt event, indicating that the area was located above the plume centre region (Krob et al., 2020).Reconstructions of magma flow directions based on magnetic susceptibility surveys around the Ponta Grossa and the Florianópolis dyke swarms (Salomon et al., 2017, and references therein) and by seismic and magnetotelluric data support the localization of the incipient Tristan-Gough plume head under African plate near the landfall of Walvis ridge (Heit et al., 2015;Jegen et al., 2016).While the exact time of plume impingement is not wellconstrained, geophysical evidence such as the limited extent of underplating around the Walvis ridge suggest an early arrival of the plume head before the opening of the South Atlantic and the formation of the first oceanic crust (Fromm et al., 2015).The earliest phase of extension related to the South Atlantic rift was active around 155 Ma and is most likely related to the opening of the Weddell Sea between southern Africa and Antarctica (Fig. 9E; Franke, 2013).Seafloor-spreading in the southernmost section of the South Atlantic rift, located near the Falkland-Agulhas Fracture zone (Fig. 9F), was established before the start of the Paraná-Etendeka flood basalt magmatism as shown by the Valanginian age (~137 Ma) of breakup unconformities in offshore sediment basins near the Cape region, South Africa and the Falkland Islands, and gradually progressed towards the north (Franke, 2013 and references therein).The South Atlantic rifting reached the area affected by the Tristan-Gough plume by ~114 Ma at the latest as indicated by oldest 40 Ar/ 39 Ar ages from the Walvis ridge seamount track, which manifested when the newly formed South Atlantic oceanic crust moved over the location of the plume head (Rohde et al., 2013).
Mantle fertilization by a deep plume can occur through different geological processes.A rising mantle plume can bring up recycled, subducted lithosphere material, which modifies the lithospheric mantle through the intrusion of melts or metasomatic veins (Pilet et al., 2008).Carbonatite melts, for example, are viable vectors of mantle refertilization (Chmyz et al., 2022), and super-chondritic bulk-rock Zr/Hf ratios positively correlated with the degree of undersaturation are reported from intraplate basalts such as the carbonatite-associated Cape Verde Islands basalts, and are often interpreted as stemming from mantle heterogeneities caused by carbonatite metasomatism in the mantle source regions of the melts (Dupuy et al., 1992), although other authors explain Zr/Hf variations with varying degree of partial melting, or fractional crystallization of phases like clinopyroxene, titanite and amphibole (David et al., 2000).Clearly super-chondritic Zr/Hf ratios which increase with the degree of bulk rock undersaturation are also observed in the Damaraland zircon crystals, and bulk rock Zr/Hf are consistently super-chondritic for the Damaraland (Fig. 8; Fig. S18 in the supplementary material), possibly indicating the involvement of carbonatite metasomatism in the genesis of the host rocks.

Tectono-magmatic model of the Damaraland igneous province
Hawaii-style, narrow hotspot magmatism is the simplest explanation for linear aligned evolved centres containing both strongly alkaline and subalkaline rock types, representing the direct magmatic manifestation of an incipient mantle plume, however, the absence of a clear ageprogression observed in many aligned igneous provinces like the Damaraland shows that their tectono-magmatic origin is often much more complicated.If the Damaraland intrusions formed by successive emplacement of melts derived from the plume following the NE plate movement of Africa over the plume, decreasing ages with a difference of up to ~12 Ma would be expected between the NE and the SW complexes (Okorusu and Cape Cross; ~350 km apart) assuming a plate velocity of 3 cm/year (Gaina et al., 2013), which is, however, inconsistent with our observation.
In addition to the age pattern, asymmetric development of the rifted Paraná and Etendeka provinces with large-volume magma extrusion in a ~1000 km wide, decentralized extension domain in South America, and more focused magmatism localized along major deformation structures in southwestern Africa, point towards a strong control of lithosphere architecture and deformation patterns on the development of the LIP (Salomon et al., 2017).The pressure effect of the lithosphere on the mantle plume melting (lithospheric lid effect) is also to be considered (Niu, 2021).A thick lithospheric lid over a plume can effectively control the depth and degree of mantle melting, and the presence of a thick lithosphere also favours magmatic underplating by acting as a physical barrier against ascending melts (Niu, 2021).
Our preferred model for the genesis of the Damaraland igneous province (Fig. 9) is that mantle was first fertilized by the incoming Tristan-Gough plume or by a previous subduction event, and after the high degree melts of the incipient plume head erupted during the main flood basalt event, continued rifting of the South Atlantic induced the reactivation of pre-existing lineaments, lithospheric extension, and associated asthenospheric mantle upwelling.This in turn caused low degree melting of the previously fertilized mantle portions to produce magmas that are channelled towards the surface through the reactivated fractures.After the dissipation of the broad plume head during the late evolution of the Paraná-Etendeka LIP, the remaining plume stem continued to deliver hot mantle material to the surface, leading to the formation of the age-progressive seamount tracks in the South Atlantic Ocean for ca.26 Myr after the Paraná-Etendeka event (Fig. 1A; Rohde et al., 2013).
Evidence for local extension in pull-apart basins along transform faults related to the South Atlantic opening exist (e.g., Holzförster et al., 1999;Raab et al., 2002) and are supportive of a magma origin by mantle melting in local upwelling zones.Lineaments already present in the Pan African basement (Fig. 1B) were plausibly reactivated during the Early Cretaceous and used as pathway for the emplacement of the complexes at shallower crustal level, which also explains the SW-NE alignment of centres in the Damaraland province.The spatial proximity of the Early Cretaceous Damaraland complexes to the 274.1 ± 7.8 Ma Ondurakorume alkaline complex and the 242 ± 6.5 Ma Kalkfeld alkaline complex (Rb-Sr isochron dating, Kramm et al., 2017), both of which cannot be genetically related to the Tristan-Gough or any other mantle plume known, provides further evidence that local reactivation of lineaments in the Damara lithosphere can plausibly generate alkaline magmas even without a direct melt supply coming from an active mantle plume (Kramm et al., 2017).
Modern examples of magmatism along active lithospheric fracture zones without evidence of plume interaction can be found in both continental and oceanic settings: the aligned volcanic centres of the Cameroon line are interpreted as being produced by progressive lithospheric thinning through basal erosion caused by an asthenospheric channel flow (Elsheikh et al., 2014), and the marine volcanic chain of the Comores-Madagascar line developed along a right-lateral transform boundary in the oceanic lithosphere (e.g., Famin et al., 2020).

Conclusions
The Early Cretaceous Damaraland centres represent a ~1.4 My-long, late-stage episode in the evolution of the Paraná-Etendeka LIP and were emplaced from 133.2 to 131.8 Ma.The temporal distribution and style of Damaraland magmatism are inconsistent with the classic Hawaii-type hotspot model.Instead, we propose that the Damaraland centres formed by the reactivation of heterogeneous mantle domains previously metasomatized/fertilized either by (1) subduction of crustal material during the formation of the Damara mobile belt in the Neoproterozoic, (2) the Tristan-Gough plume through basal erosion of the lithosphere, and/or (3) the intrusion of the lithospheric mantle in form of metasomatic veins, and not as a direct consequence of melting in the plume.Lithospheric extension and consequent mantle upwelling produced different types of melts, which then intruded into the Damaraland basement along preexisting fault zones reactivated to accommodate the rifting of the South Atlantic.Timing and location of emplacement are therefore mainly controlled by lithospheric processes rather than by direct effects of the mantle plume.While the predominantly alkaline complexes (e.g., Messum, Okenyenya, Etaneno, Okorusu) show evidence for magma evolution paths in near-closed systems, the centres that evolved towards the granite minimum to produce granites or rhyolites (e.g., Erongo, Spitzkoppe, Otjihorongo, Paresis) locally assimilated parts of the host crustal material.
the work reported in this paper.

Y
.Sun et al.

Fig. 1 .
Fig. 1.Overview maps showing (A) the location of the Damaraland province within the Paraná-Etendeka LIP, and locations of associated Early Cretaceous dyke swarms in the Paraná basin, cratonic domains, and the age-progressive Walvis ridge/Tristan-Gough hotspot track in the South Atlantic Ocean.Red circles refer to sites dated by 40 Ar/ 39 Ar (Rohde et al., 2013) and the present-day locations of the Tristan da Cunha and Gough islands.(B) Locations of sampling sites, major Precambrian lineament zones, and important Early Cretaceous volcanic outcrops in the Damaraland province.Each complex is marked with an individual symbol with symbols in red denoting exclusively subalkaline complexes, symbols in green representing bimodal complexes, and symbols in blue referring to exclusively alkaline complexes.(C) Locations of the alkaline complexes of the Kalkfeld group with the emplacement ages shown for the older generation of intrusions.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
In the alkaline complexes, zircon from Etaneno syenite and Okorusu syenite are characterized by high Zr/Hf values (70-110).While the zircon grains from the subalkaline complexes have mostly negative εHf(t) values from − 5 to − 12, zircons from bimodal intrusions tend to show positive to moderately negative εHf(t) of +6 to − 6. Brandberg is the only subalkaline intrusion with εHf(t) around 0. Zircon crystals from the alkaline complexes are generally characterized by high zircon Zr/Hf of 70-100 Y. Sun et al.

Fig. 3 .
Fig. 3. Uranium-Pb concordia diagrams of Damaraland zircon and baddeleyite dated by LA-ICP-MS, with concordia ages calculated in IsoplotR(Vermeesch, 2018).Error ellipses are plotted as 2σ, and reported MSWD and chi-square probabilities were calculated in IsoplotR for combined concordance and equivalence.

Fig. 4 .
Fig. 4. Comparison of single-crystal 206 Pb/ 238 U Damaraland zircon dates by CA-ID-TIMS (uncertainties stated and plotted as 2σ) with published weighted mean zircon/baddeleyite TIMS ages for silicic magmas and dyke rocks in the Paraná province, Brazil(Rocha et al., 2020;Rocha et al., 2023).Dates shown in grey are excluded from the calculation of the weighted mean age.

Fig. 5 .
Fig. 5.Total alkali versus silica diagram showing the compositional variability of igneous rocks in the Damaraland province, which follow two different melt evolution paths towards the phonolite minimum or the granite minimum.Literature data compiled fromBowden et al., 1990;Haapala et al., 2007;Harris et al., 1999;Miller, 2008;Mingram et al., 2000.
both.The εHf(t) values in zircons from Cape Cross granite and Brandberg granite are around 0 or only moderately negative (0 to − 5), suggesting little (if any) path contamination.In the Spitzkoppe, Otjihorongo and Erongo complexes, which exclusively evolved towards the granite minimum, the much lower εHf(t) between − 5 and − 12 most likely stems from considerable path contamination.The most likely contaminants are Panafrican Damara basement granites and Flysch series metasediments including abundant quartz-plagioclase-biotite schists, which are characterized by εNd (at 132 Ma) values from − 4 to − 19, and

Fig. 6 .
Fig. 6.Geochemistry plot showing the Hf isotope compositions and Zr/Hf ratios of Damaraland zircon grains (dated to 140-124 Ma) from subalkaline, bimodal, and alkaline complexes.Negative εHf(t) reflect either source and/or path contamination in the host lithology, and the Zr/Hf reflect differences in melt primary sources.

Fig. 7 .
Fig. 7. Rank-order plot of zircon/baddeleyite U-Pb TIMS and plagioclase 40 Ar/ 39 Ar ages for Early Cretaceous intrusive rocks and volcanic rocks from the Damaraland and Paraná provinces compiled from this study and literature data, showing the duration and temporal evolution of magmatism in the Paraná-Etendeka LIP.Error bars represent 2σ uncertainties.

Fig. 8 .
Fig. 8. Compilation of bulk-rock Sr-Nd isotope compositions at 132 Ma for the Damaraland complexes and Etendeka volcanic rocks (grey fields) showing the role of crustal assimilation in the genesis of the complexes.Alkaline complexes in blue are uncontaminated, whereas granitic complexes in red (e.g., Erongo) and bimodal complexes in green (e.g., Paresis, Okenyenya) can display significant crustal input.Fields with dashed outlines represent Sr-Nd isotopic compositions of possible crustal contaminants, including Panafrican basement granite and metasediments of the Damara mobile belt, and Pre-Damara basement gneisses of the Congo craton.Data compiled fromBühn and Trumbull, 2003;Ewart et al., 1998;Ewart et al., 2004;Frindt et al., 2004;Harris et al., 1999;Le Roex and Lanyon, 1998;Milner and Le Roex, 1996;Marsh et al., 2001;Milner et al., 1993;Mingram et al., 2000;Schmitt et al., 2000;Teklay et al., 2020;Trumbull et al., 2000;Trumbull et al., 2003;Trumbull et al., 2004;Zhou et al., 2020.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)