Mesoarchean magmatism is widespread on the eastern margin of the Kaapvaal Craton, but its origin is still poorly understood, mainly because geochemical data is rare. To shed some light on the source of this Mesoarchean magmatism and to relate different Mesoarchean volcanic sequences to each other, we provide major and trace element data as well as Hf-Nd isotope compositions of amphibolites sampled close to the Kubuta Ranch in south-central Eswatini. These amphibolites, so far, were of unknown correlation to any volcanic sequence in Eswatini or South Africa. Hence, the aim of our study is to characterize the mantle source composition of these volcanic rocks and, furthermore, to constrain their genetic relation to other volcanic sequences in Eswatini and South Africa. Our findings reveal that, based on coherent trace element patterns and similar Nd isotope characteristics, the Kubuta volcanic rocks can be genetically linked to the ca. 3.0 Ga Usushwana Igneous Complex in West-Central Eswatini and the ca. 2.9 Ga Hlagothi Complex located in the KwaZulu-Natal province. In contrast, the coeval ca. 3.0 Ga Nsuze and ca. 2.9 Ga Mozaan Groups (Pongola Supergroup) of south-central Eswatini and northern KwaZulu-Natal province have slightly enriched compositions compared to the newly sampled Kubuta volcanic rocks. Our results suggest that the Nsuze and Mozaan Groups were sourced from a primitive mantle reservoir, whereas the Usushwana, Hlagothi, and Kubuta mafic rocks were derived by melting of a more depleted mantle source comparable to that of modern depleted MORB. Furthermore, our assimilation-fractional crystallization (AFC) calculations and Nd isotope constraints reveal that some samples were contaminated by the crust and that the crustal contaminants possibly represent felsic rocks related to the ca. 3.5 Ga crust-forming event in the Ancient Gneiss Complex. Alternatively, melting of a metasomatized mantle or plume-lithospheric mantle interaction may also produce the trace element and isotopic compositions observed in the samples. From a synthesis of our geochemical observations and age data from the literature, we propose a refined petrogenetic model, for a continental flood basalt setting in a Mesoarchean large igneous province on the eastern Kaapvaal Craton. Our petrogenetic model envisages two magma pulses sourced from a primitive mantle reservoir that led to the formation of the Nsuze (first) and Mozaan (second) lavas. Conductive heating of ambient depleted mantle by the mantle plumes caused partial melting that led to the formation of the Usushwana Igneous Complex associated with the first magmatic event (Nsuze) and the Hlagothi Igneous Complex associated with the second magmatic event (Mozaan). However, due to lacking age data of sufficient resolution, it is not possible to affiliate the Kubuta lavas to either the first or the second magmatic event.
The Archean eon covers much of Earth’s history and was a period of major crust forming events (e.g., [1–3]) that culminated in the stabilization of the ancient cratons. One of these major crust-forming events occurred in the Mesoarchean (ca. 3.0 to 2.9 Ga) on the eastern Kaapvaal Craton in southern Africa producing large volumes of continental flood basalts in an intracratonic rift basin (e.g., [4–7]). These Mesoarchean basalts have well-preserved magmatic textures and preserve their primary geochemical composition, which can be used to place constraints on the genetic relationships to other units of the Kaapvaal Craton.
The correlation of lavas of the eastern Pongola basin and western Witwatersrand basin in terms of age and trace element composition led Gumsley et al.  to the conclusion of the existence of a so far unknown Mesoarchean large igneous province, which might extend further West to the Dominion Group . Mainly based on precise age determinations and field observations, these authors proposed that the large igneous province was formed by two magma pulses derived from a short-lived mantle plume. In their petrogenetic model, the first magma pulse formed the basaltic lavas of the ca. 3.0 Ga Nsuze Group (Pongola Supergroup) and the gabbros and cumulate rocks of the Usushwana Igneous Complex (Figure 1; Gumsley et al. ). After some time of subsidence and sedimentation, the second pulse of magma formed the ca. 2.9 Ga lavas of the Mozaan Group (Pongola Supergroup; Figure 1) and the mafic sills of the Hlagothi Complex .
However, there is some disagreement on the nature of the parental magmas that formed the Nsuze volcanic rocks. While Hegner et al.  proposed an ultramafic komatiitic parental magma for the Nsuze volcanic rocks, Armstrong et al.  suggested a primitive basaltic parental magma. Likewise, the proposed crustal contamination processes that affected these parental magmas vary from bulk assimilation  and assimilation-fractional crystallization  to magma mixing of basaltic and felsic melts . However, all aforementioned studies agree that a contamination by felsic crustal rocks from the Ancient Gneiss Complex (AGC) of Eswatini occurred. Yet, these studies do not specify, if the crustal contaminant to the Nsuze volcanic rocks belongs to the 3.2 Ga or 3.5 Ga crust-forming event within the AGC.
In this study, we sampled amphibolites from the Kubuta area of southwestern Eswatini. The affiliation of the Kubuta samples is unknown, but from a combination of major and trace element, as well as Lu-Hf and Sm-Nd isotope data, we will show that they could either be related to volcanic rocks of the 2.98-2.84 Ga Pongola Supergroup (Figure 1) [10–14] or the ca. 2.99 to 2.87 Ga Usushwana Igneous Complex (Figure 1) [4, 5] and the ca. 2.86 Ga Hlagothi Complex . Alternatively, the Kubuta volcanic rocks could be related to the ca. 2.66 Ga White Mfolozi Dike Swarm (WMDS)  further south in the KwaZulu-Natal province.
Using our extensive data set, we aim to better constrain the composition of the mantle sources and the contamination history of the aforementioned sequences to test the reliability of the petrogenetic model of Gumsley et al.  that was mostly based on geochronology. In combination with literature data, we propose a refined petrogenetic model for the Nsuze, Mozaan, Usushwana, Hlagothi, and Kubuta volcanic rocks that also takes geochemical evidence into account.
2. Geological Overview and Sample Petrography
2.1. Geological Overview
The eastern Kaapvaal Craton, located in Eswatini and eastern South Africa, records a continuous geological history from the early Archean to the Neoarchean recording the stabilization of continental lithosphere in this region. The first stabilized continental crust is represented by the 3.66-3.20 Ga Ancient Gneiss Complex (AGC) [19–22] that predominantly contains complexly deformed Paleo- to Mesoarchean grey gneisses that are rarely interlayered with amphibolites (e.g., [23–25]). The oldest component of the AGC is the ca. 3.66-3.45 Ga Ngwane Gneiss (Figure 1) that underwent ductile metamorphism under amphibolite to granulite conditions [20, 26]. The Ngwane Gneiss is intruded by the genetically linked 3.48-3.43 Ga Tsawela Gneiss [26–28] as well as by several younger gneisses and granitoids like the 3.23-3.22 Ga Usutu Intrusive Suite in northern and central Eswatini  and the 3.28-3.24 Ga Nhlangano and Mahamba Gneisses in south-central Eswatini (Figure 1) . The ca. 3.1 Ga Pigg’s Peak and Mpuluzi Batholiths [31, 32] represent the last major intrusions exposed in the eastern Kaapvaal Craton [22, 30, 33]. These intrusions separate the AGC in the Southeast from the Barberton Greenstone Belt in the Northwest. In its most southerly occurrence, the Mpuluzi Batholith forms the basement to the ca. 2.98-2.84 Ga Pongola Supergroup (Figure 1) [10–14].
The Pongola Supergroup was deposited in an area of about 24000 km2 extending from south-central Eswatini into northern and central KwaZulu-Natal in eastern South Africa and belongs to one of the oldest cratonic cover sequences in the world . Recently, new age determinations lead to the finding that the Pongola Supergroup also extends further West into the 2.96 Ga Dominion Group revealing that a craton-wide igneous province was active between ca. 2.99 and 2.96 Ga . The succession was deposited in a volcanic rift setting  between contrasting crustal blocks of the AGC with generally older granitoids ranging in age from ca. 3.6 to 3.2 Ga to the Northeast  and ca. 3.3 to 3.2 Ga to the Southwest  of the Pongola basin. In the northern and central areas of the basin, the metamorphic grade of the volcanic and sedimentary rocks of the Pongola Supergroup ranges from largely greenschist up to granulite facies in some areas .
The Pongola Supergroup is divided into the lower, mostly volcanic, ca. 2.98 Ga Nsuze Group and the upper, predominantly sedimentary, ca. 2.86 Ga Mozaan Group (Figure 1) [37, 38]. The volcanic rocks of the Nsuze Group range from metabasalts to metarhyolites  and were deposited in an intracontinental environment  under subaerial to locally subaqueous conditions [9, 40]. The Mozaan Group rests unconformably on the Nsuze Group and mostly contains sandstone, shale, mudstone, and locally banded iron formations and stromatolites indicating deposition in a shallow marine environment [37, 41–44]. Towards the top of the succession, volcanic rocks become more dominant again .
The Nsuze and Mozaan Groups were intruded by the Usushwana Igneous Complex (UIC; Figure 1), the Thole Complex, and the Hlagothi Complex. The UIC has an emplacement age of Ma . However, recent U-Pb baddeleyite age dating by Gumsley et al.  on gabbroic samples from the Piet Retief Suite (now re-named to Mkhondo Suite) that is part of the UIC yielded significantly older ages up to 2989 Ma falling in the age range of Nsuze volcanic rocks. Based on their new ages, Gumsley et al.  interpreted the Usushwana Igneous Complex as being a feeder dike system for Nsuze lavas. The Ma Hlagothi Complex and the related Ma Thole Complex consist of several gabbroic sills that are found in the KwaZulu-Natal province and are interpreted to be the feeder dikes for the Mozaan flood basalts . Based on the similarity of ages and composition, Gumsley et al.  proposed that the UIC and the Thole and Hlagothi Complexes, together with the Mozaan volcanic rocks form part of a Mesoarchean large igneous province located in the southeastern part of the Kaapvaal Craton and possibly extending into the eastern Witwatersrand block.
Later intrusions of ca. 2.85-2.70 Ga granitoids  were triggered by extensional deformation that also caused local deformation of the Pongola Supergroup . Following this event, ca. 2.70 to 2.66 Ga extensional dolerite dike swarms, sometimes referred to White Mfolozi Dike Swarm, crosscut the eastern Kaapvaal Craton terminating its Archean crust formation history (e.g., [15, 46–48]).
2.2. Sampling and Sample Petrography
Seven samples have been collected from the area of the Kubuta Ranch in south-central Eswatini (Figure 2). These samples are referred to as ‘Kubuta volcanic rocks’ in the following. Samples KS-SW 31 to 35 have been collected east of the Sibowe shear zone at the MR 25 (Figure 2) in a geologic succession that was mapped by the Geological Survey of Eswatini as “amphibolite of uncertain correlation” . In contrast, samples KS-SW 45 and 47 have been collected west of the Sibowe shear zone as boulders in the river bed of the Matsanjeni River (Figure 2). Hence, one of the aims of this study is to constrain the relationship of these rocks in terms of mantle source and contamination history to other volcanic sequences in this region, i.e., the 2.98-2.84 Ga Pongola Supergroup (Figure 1) [10–14], the ca. 2.99 to 2.87 Ga Usushwana Igneous Complex (Figure 1) [4, 5], the ca. 2.86 Ga Hlagothi Complex , or the ca. 2.66 Ga White Mfolozi Dike Swarm (WMDS)  further south in the KwaZulu-Natal province.
All samples appear relatively fresh with preserved primary magmatic mineral textures. However, when taking a closer look at the samples in thin section, the rocks can be identified to have undergone greenschist-facies metamorphism expressed by variably chloritized mineral assemblages and formation of epidote. Sample KS-SW 31 has a microgabbroic texture indicating a moderate cooling rate and containing minor relict olivine and abundant relict clinopyroxene and plagioclase (Figures 3(a) and 3(b)). Samples KS-SW 33 to 35 exhibit an intergranular texture with usually bigger needle-shaped amphibole crystals and a medium-grained matrix of clinopyroxene, orthopyroxene, plagioclase, and epidote. Sample KS-SW-32 is a cumulate containing orthopyroxene of up to 5 mm in a diameter. Samples KS-SW 33 and 35 have radiating and branching fibrous plagioclase crystals preserved (Figures 3(c) and 3(d)). Carbonatization along mineral cracks is rare in most samples, except in KS-SW 32, which shows relatively strong carbonatization and chloritization. However, the original magmatic texture of this sample is still preserved.
In contrast, samples KS-SW 45 and 47 show a porphyric texture having either 1-2 mm amphibole prisms and clinopyroxene crystals floating in a fine-grained plagioclase-bearing matrix (KS-SW 45; Figure 3(e)) or plagioclase and clinopyroxene crystals floating in a fine-grained orthopyroxene-bearing matrix (KS-SW 47). These two samples show almost no carbonatization. Furthermore, apatite can be found in plagioclase phenocrysts and different opaque phases occur as accessory components within the matrix.
3. Analytical Techniques
The rock samples were cut using a diamond blade saw to remove alteration rims, and thin sections were prepared. The remaining sample material was crushed in a steel jaw crusher and ground in an agate mill. Whole rock major and some trace element compositions including the contents of H2O and CO2 were analyzed by X-ray fluorescence (XRF) on Li2B4O7-flux fusion discs using a PANalytical Axios X-ray spectrometer at the ElMiE Lab at the German Research Centre for Geosciences (GFZ) in Potsdam, Germany (Table 1). Loss on ignition (LOI) was determined by weight difference after fusion. Reproducibility was determined on three certified reference materials (crm) and is within the analytical precision, which is better than 2% for main elements and better than 10% for trace elements.
The trace element compositions of all samples were analyzed using a Thermo Finnigan Element XR SF-ICP-MS at Freie Universität Berlin following the protocols of Schneider et al. . Sample powders were digested in Parr® pressure vessels in concentrated HF-HNO3 for 72 hours and were subsequently dried down. The evaporation step was repeated with 2 ml concentrated HNO3. After evaporation to dryness, samples were diluted in 3% HNO3 for trace element measurements. The analytical procedure for trace element analysis was carried out following Schneider et al. . Relative total uncertainties (RSD) for all measured elements usually range from 4% to 18%. The accuracy of our measurements can be inferred from multiple analyses of the BHVO-2 reference material and was better than 5% for most elements and between 6% and 12% for Ga, Sr, Y, and Nb (Table 2). Additionally, the total procedural uncertainty can be inferred from sample replicates that were analyzed in four analytical sessions and were better than 5% for the REEs, Ba, Th, U, and Co, better than 11% for Ga, Y, Ta, Sc, V, and Cu, and 25% for Mo and 40% for Pb (Table 2).
Hafnium and Nd isotope compositions and Lu, Hf, Sm, and Nd concentrations by isotope dilution were analyzed (Table 3) using a Thermo Scientific Neptune multicollector ICP-MS at the University of Cologne (Germany), operated in static mode. The analyses were obtained by isotope dilution using mixed 176Lu-180Hf and 149Sm-150Nd tracers following the protocols of Münker et al.  and Weyer et al. . Chemical separation of Lu, Hf, Sm, and Nd followed the procedures outlined in Pin and Zalduegui  and Münker et al. .
Values of 176Hf/177Hf and 143Nd/144Nd were corrected for mass bias using the exponential law and assuming and , respectively. Repeated analyses of the Cologne AMES Hf solution standard yielded an external reproducibility of ±41 ppm (2 s.d., ) for 176Hf/177Hf, whereas repeated analysis of the AMES Nd and JNdi-1 standards yielded external reproducibilities of 21 ppm (2 s.d., ) and 70 ppm (2 s.d., ) for 143Nd/144Nd, respectively, over the course of the analytical session. Reported values are given relative to for JNdi (i.e., for La Jolla Nd ) and for the Münster AMES Hf standard, which is isotopically indistinguishable from the JMC-475 standard for Hf . All procedural blanks handled along with the sample batches during the course of this study were better than <3 pg for Lu, <37 pg for Hf, <22 pg for Sm, and<132 pg for Nd and therefore were insignificant. CHUR parameters used to calculate epsilon values are from Bouvier et al.  and an age of 2980 Ma was estimated (see discussion in Section 6.2).
4.1. Major Elements
The SiO2 and MgO contents in the samples range from 44.4 to 53.2 wt.% and from 7.17 to 22.9 wt.%, respectively (Table 1). Magnesium numbers (Mg#) range from 54 to 78 corresponding to the variability in MgO and FeOtotal concentrations. Compared to literature data, our samples have much higher MgO contents than the Nsuze Group volcanic rocks (Pongola Supergroup; Figures 4(a)–4(d)) [9–11, 46] and volcanic rocks from the WMDS . However, our Kubuta samples agree in MgO content with the volcanic rocks of the Mozaan Group (Pongola Supergroup) [10, 57], the Hlagothi Complex , and the Usushwana Igneous Complex . The two samples of the latter unit from the study of Nhleko  were originally grouped as Mozaan volcanic rocks. However, Nhleko  supposed in his thesis that they might rather belong to the Usushwana Igneous Complex based on their distinct geochemical composition. Adopting this interpretation, we plot these two samples as UIC in the diagrams of this study. Unfortunately, there is currently no more literature data from the UIC to test a possible overlap in element abundances of our Kubuta samples with the UIC; only chondrite-normalized REE data are published for this unit .
In Fe2O3 vs. MgO space, a positive correlation of samples of the Nsuze and Mozaan suites below 5 wt.% MgO (Figure 4(a)) suggest fractionation of pyroxene and titanomagnetite. In addition, fractionation of plagioclase and pyroxene can be identified at MgO contents below 5 wt.% by the positive correlations of Al2O3 and CaO vs. MgO (Figures 4(b)–4(d)), respectively. Furthermore, the positive correlation of Eu anomalies with CaO contents indicates plagioclase fractionation. Fractionation of apatite is indicated by the negative correlation of P2O5 above a SiO2 content of 60 to 65 wt.% (Figure 4(e)), and fractionation of zircon is indicated by the negative correlation of Zr and SiO2 above a SiO2 concentration of around 70 wt.% (Figure 4(f)). However, another trend is visible in Figure 4(f) which cannot be explained by zircon fractionation, where Zr increases with increasing SiO2, possibly representing the involvement of Zr-rich partial melts of the continental crust. Compared to the fractionation trends observed from literature data, our samples mostly plot at the high-Mg side of these trends, resembling the high-Mg samples from the Usushwana Igneous Complex and Hlagothi Complex.
4.2. Trace Elements
The trace element concentrations of the samples are summarized in Table 2. The primitive mantle-normalized trace element patterns for Kubuta samples show an enrichment of light rare earth elements (LREE) and flat to slightly depleted heavy REEs (Figures 5(a) and 5(b)). Consequently, most samples exhibit (La/Yb)CN (i.e., chondrite-normalized) values greater than one, ranging from 1.90 to 6.15, except for sample KS-SW 47 that has (La/Yb)CN of around one. Four samples have (Gd/Yb)CN values around one and in three samples the (Gd/Yb)CN values range from 1.28 to 2.21. Pronounced negative anomalies are observed for Nb , Ta , and Ti. Values of Th/U range from 3.3 to 6.0 in agreement with Archean upper crustal values, i.e., (e.g., ). Furthermore, negative Eu anomalies are present for four samples ranging from , indicating fractionation of plagioclase, whereas sample KS-SW 31 shows a small positive of 1.15. All other samples lack Eu anomalies. Furthermore, all samples show elevated concentrations of fluid mobile elements such as Rb, Ba, and Pb.
When compared to the literature data for the Nsuze Group [9–11, 46] and Mozaan Group [10, 57] volcanic rocks, our samples are generally less enriched in incompatible trace elements, yet exhibit similar (La/Yb)PMN and (Gd/Yb)PMN values and elemental anomalies (Figure 5(a)). When compared to samples from the Usushwana Igneous Complex  and Hlagothi Complex , our samples match well their elemental abundances and exhibit similar (La/Yb)PMN and (Gd/Yb)PMN values and elemental anomalies (Figures 5(a) and 5(b)). However, compared to volcanic rocks from the White Mfolozi Dike Swarm , the Kubuta samples have higher (La/Yb)PMN and (Gd/Yb)PMN values and pronounced negative Nb and Ta anomalies that are not observed in the volcanic rocks of the WMDS.
4.3. Hafnium-Neodymium Isotope Compositions
The Lu-Hf and Sm-Nd isotopic compositions and Lu, Hf, Sm, and Nd concentrations of our sample suite are summarized in Table 3. For calculation of initial εHf and εNd values, an age of 2980 Ma was applied (; see discussion in Section 6.2). Most of the Kubuta samples plot on the terrestrial array in εHf(t) vs. εNd(t) space after Vervoort et al.  that we use as a measure for a possible decoupling of initial Hf and Nd isotope compositions. However, sample KS-SW 31 has much higher εHf(t) (+26.6) at its εNd(t) value (+6.3) and plots off the array (see in Table 3). Thus, this sample was not included in the plot in Figure 6, considering it as disturbed. All other samples have εHf(2980) and εNd(2980) values ranging from -4.8 to +3.5 and from -2.5 to +1.1, respectively.
Literature data for rocks of the Nsuze Group and Usushwana Igneous Complex yielded εNd(t) values ranging from -1.8 to -7.8 [5, 11] and -3.7 to -1.6, respectively. These values are in good agreement with three of our samples yielding negative εNd(t) (and negative εHf(t)) values. However, the other three Kubuta samples have chondritic to slightly elevated εNd(t) (and εHf(t)) values at comparable neodymium concentrations (Figure 6). Unfortunately, whole-rock Lu-Hf isotope compositions were not available in the literature for further comparison. Likewise, no Hf-Nd isotope data were available for the Mozaan volcanic rocks, the Hlagothi Complex, and the White Mfolozi Dike Swarm to compare with our Kubuta samples.
5.1. Effects of Alteration and Metamorphism
Great care must be taken when interpreting the geochemistry of Archean volcanic rocks, because these rocks underwent metamorphism that likely caused mobility of some elements which lead to alteration of the element inventory of these rocks. In Archean volcanic rocks, elements typically known as relatively mobile during metamorphism include Rb, K, Na, Sr, Ba, Ca, Fe, P, and Pb [62, 63]. Thus, these elements should be used with caution to constrain magmatic processes. In contrast, the effects of metamorphic alteration on high field strength elements (HFSE), REEs, Al, Cr, and Ni are usually minor, and therefore, these elements are mostly considered immobile on the whole-rock scale [63–66].
To distinguish the elemental disturbance of Archean greenstones, Polat et al.  developed five requirements that need to be fulfilled in order to consider the greenstones to be largely undisturbed in their trace element composition. Most of our investigated greenschist-facies metabasalts meet all five criteria of Polat et al.  in that they (1) generally have volatile (i.e., H2O plus CO2) contents lower than 6 wt.% (Table 1); (2) show coherent trace elements patterns of immobile elements normalized to PM (Figure 5); (3) display correlations of REEs, Th, Ti, and Nb abundances with Zr content; (4) show no correlation of their Th-Nb-La interelement ratios with volatile (i.e., H2O plus CO2) abundance; and (5) show no correlation of element abundance with the degree of chloritization or carbonatization. However, sample KS-SW 32 has a high volatile content (i.e., H2O plus CO2) of 7.9 wt.%, resulting from its relatively strong chloritization and carbonatization. This sample also has much higher MgO (22.9 wt.%), lower Na2O (0.01 wt.%) and K2O (0.03 wt.%) contents than less altered samples, likely indicating mobility of Mg, Na, and K during greenschist-facies metamorphism. Consequently, we consider the major element composition of sample KS-SW 32 to be disturbed. In contrast, chloritization and carbonatization appear not to have disturbed the REE, HFSE, and the Hf-Nd isotope systematics of this sample. All other samples do not seem to have disturbed major and trace elemental compositions.
As stated above, sample KS-SW 31 plots off the modern terrestrial array  in εHf(t) vs. εNd(t) space towards much higher εHf(t) at a given εNd(t) compared to the other Kubuta samples and compared to literature Nd isotope data for the Pongola and Usushwana suites [5, 11]. These high εHf(2980 Ma) at εNd(2980 Ma) could be a result of incorrect age assignment: This sample might belong to a different, much younger, volcanic suite, possibly to the Karoo flood basalts, leading to overestimated radiogenic εHf(t) at εNd(t) values. Other arguments in favor of this hypothesis are, for instance, the different petrographic textures (microgabbro) leading to the positive Eu anomaly of this rock compared to the other samples. Consequently, this sample is excluded from our further discussion.
5.2. Classification and Age Constraints
Because of possible disturbance of the major element composition of our samples, we used the discrimination diagram of Winchester and Floyd  that was updated by Pearce  for all magma types applying the immobile trace element ratios Zr/Ti and Nb/Y, instead of the total alkalis vs. silica (TAS) diagram after Le Maitre et al. . In the Zr/Ti vs. Nb/Y diagram, the former ratio reflects the degree of differentiation, and therefore is comparable to SiO2, and the latter tracks enrichment, and, hence, is comparable to the alkali content. With this approach, the Kubuta lavas can be classified as subalkaline () basalts and minor basaltic andesites overlapping with the Hlagothi and Usushwana Igneous Complexes, the White Mfolozi Dyke Swarm, and most of the Nszue and Mozaan lavas (Figure 7).
Based on our isotope analysis, we obtained a Sm-Nd whole-rock regression line corresponding to an errorchron of Ma (; Figure 8(a)) and a Lu-Hf whole-rock regression line that corresponds to an errorchron of Ma (; Figure 8(b)) for five samples of the Kubuta sample suite. However, the uncertainties on the calculated ages are quite large which may result from the relatively small spread in 147Sm/144Nd and 176Lu/177Hf ratios in our samples. Alternatively, the large uncertainties on the regression lines may indicate that our samples might not be cogenetic or be contaminated by crustal material. However, the data still give a first-order approximation on the geological unit to which the Kubuta metabasalts could be assigned to. Hence, the most likely units that could correlate with the Kubuta samples in terms of age and composition are the mafic/ultramafic Nsuze and Mozaan Groups of the Pongola Supergroup, as well as the mafic/ultramafic Usushwana Igneous Complex (UIC), the mafic/ultramafic Hlagothi Complex, and the mafic volcanic rocks of the White Mfolozi Dike Swarm.
The age of the Nsuze Group is well constrained by zircon U-Pb dating of evolved volcanic rocks which range from Ma to Ma in age [10–12]. In contrast, a Sm-Nd whole-rock errorchron ‘age’ obtained from basaltic to rhyolitic samples yielded a slightly younger, but within uncertainty overlapping, age of Ma . The age of the overlying Mozaan Group ranges from Ma  to Ma , respectively. The mostly intrusive Usushwana Igneous Complex (UIC) yielded Rb-Sr whole-rock age data of Ma  and Ma , postdating the Nsuze volcanic rocks. However, recent U-Pb baddelyite dating by Gumsley et al.  yielded slightly older ages ranging from Ma to Ma that indicate a synchronous evolution of the UIC and the Nsuze Group. Additionally, the ages for the Mozaan Group and UIC overlap within error with the ca. 2866 Ma mafic to ultramafic sills of the Hlagothi Complex  that intrude the basal units of the Nsuze Group. The WMDS was dated to Ma , also overlapping with our Sm-Nd and Lu-Hf whole-rock errorchron dates.
In summary, based on our Sm-Nd and Lu-Hf errorchron ages, the Kubuta volcanic rocks could be assigned to all the aforementioned volcanic units. However, in terms of trace elemental patterns combined with the Nd isotope composition, the Kubuta volcanic rocks share more similarities with the lavas of the Nsuze and Mozaan Groups [9–11, 57] as well as with the Usushwana [c.f., 5] and Hlagothi Igneous Complexes [c.f., 6] rather than with the lavas of the White Mfolozi Dike Swarm [c.f., 15] (Figures 5(a), 5(b), and 6(a)).
Hence, we exclude the possibility that the Kubuta volcanic rocks are related to the White Mfolozi Dike Swarm. Consequently, based on the good match of major, trace element, and Nd isotope abundances, we assign the Kubuta samples to be part of the proposed large igneous province of Gumsley et al. , consisting of mafic and ultramafic rocks from the Mozaan Group, the Hlagothi Complex, and the Usushwana Igneous Complex. However, the role of these units in the proposed large igneous province and their relation to the Nsuze lavas needs further investigation.
5.3. Origin of the Enriched Geochemical Composition of the Lavas
5.3.1. Crustal Contamination Mechanisms
Crustal contamination of Nsuze volcanic rocks has been suggested by several authors based on the elemental composition and negative εNd(t) values of the lavas [5, 9, 10, 46]. As potential contaminants to the Nsuze volcanic rocks, granitoids from the Ancient Gneiss Complex (AGC) were favored [5, 9]. The proposed assimilation processes vary from bulk assimilation of ca. 3.5 Ga old granitoid crust  to assimilation-fractional crystallization (AFC) of more granitic AGC crust . In addition, Nhleko  proposed magma mixing of basaltic melts and felsic melts derived from partial melting of the continental crust to explain the diverse compositions of Nsuze volcanic rocks.
The arguments in favor of an involvement of continental crust in the petrogenesis of the volcanic rocks from the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta sample suites are mostly their enrichment of incompatible elements, especially La and Th, in conjunction with negative anomalies for Nb and Ta (Figures 4 and 5). Additionally, consistent primitive mantle normalized trace element patterns and mostly negative εNd(t) values for the aforementioned volcanic rocks suggest a similar contamination history for all suits. In Th/Yb vs. Nb/Yb space, all samples plot well above the modern MORB-OIB array towards higher Th/Yb and Nb/Yb indicating substantial contamination of the parental magmas resulting from Archean crustal input and/or partial melting and fractional crystallization (e.g., ) (Figure 9). However, one type of crustal contaminant alone cannot account for the diversity in compositions of the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta samples (e.g., Figure 9).
Possible crustal contaminants that could have modified the composition of the aforementioned ca. 2.9 to 2.8 Ga volcanic rocks are granitoid gneisses from the AGC either belonging to the ca. 3.5 Ga (e.g., [5, 9, 10]) or the ca. 3.2 Ga crust-forming events (for a review on the ca. 3.5 and 3.2 Ga crust-forming events in the eastern Kaapvaal Craton, please refer to ). Contamination most likely occurred at the upper crustal levels (e.g., [5, 9, 10]). Arguments in favor of AFC are fractionation of pyroxene and plagioclase as indicated by the fractionation trends in Figure 4 and by the depletions in Ti and Eu illustrated in Figure 5, as well as by the typical crustal contamination trend displayed in εNd(t) vs. La/Sm space (Figure 6(b)).
In contrast, the increase in Zr concentration with increasing SiO2 for Nsuze rhyolites, illustrated in Figure 4(f), is likely a result of the increased alkalinity of these rocks which hampers zircon fractionation (e.g., ), suggesting that these rocks formed from crustal melts (e.g., ). Hence, magma mixing of mafic melts with felsic partial melts of older crust could explain the evolved compositions of the volcanic rocks of the Nsuze and Mozaan Groups, too, as previously suggested by Nhleko .
In the following, we will discuss the mechanisms of the two different contamination processes, i.e., magma mixing and assimilation-fractional crystallization.
(1) The Role of Magma Mixing. Partial melting of AGC crust and subsequent magma mixing with mafic underplated magma was suggested by Nhleko  to have formed the Nsuze volcanic rocks. One argument for partial melting of AGC crust is based on one rhyolite from the study of Hegner et al.  which yields a strongly negative εNd(t) value of -7.8. Partial melting of crustal material would shift the resulting magma composition towards higher La/Yb and Zr concentrations and towards slightly lower Th/Yb and Nb/Yb than the AGC crustal rocks. This is observed for rhyolites of the Nsuze and two andesites of the Mozaan and Hlagothi suites, respectively (Figure 9). Hence, it is possible that at least some andesites and basaltic andesites of the Nsuze and Mozaan Groups might have been formed by magma mixing with subsequent fractional crystallization.
To explore the role of magma mixing, we calculated mixing paths using the Microsoft Excel calculation spreadsheet of Ersoy and Helvacı . As a mafic mixing end member, we chose basalt AG-F of the Hlagothi Complex , because it has the lowest Nb/Yb at low Th/Yb, low La/Yb, and medium MgO (8.36 wt.%) and therefore serves as the best approximation of the composition of the parental melt. As felsic mixing end member, we used rhyolite NZ-39 of the study of Nhleko , because it has high SiO2 (72.2 wt.%), very low MgO (0.05 wt.%), and the highest Zr concentration (601 ppm), and consequently most likely represents a partial melt of the crust. The results of magma mixing are illustrated in Figure 10. With this approach, most of the Nsuze and Mozaan volcanic rocks having can be reproduced by 10 to 50% mixing of a rhyolitic melt with a basaltic melt (Figures 10(a) and 10(b)). However, many basalts, basaltic andesites, and andesites of the Nsuze and Mozaan Groups plot towards lower Th/Yb at similar Nb/Yb (Figure 10(a)). These lavas cannot be reproduced by magma mixing. Furthermore, magma mixing would have caused more heterogeneity in Nd isotope composition as a function of magma composition, i.e., the more felsic the magma, the lower the εNd(t) should be. However, this is not reflected by the relatively uniform negative εNd(t) values in most Nsuze, Usushwana, and three Kubuta volcanic rocks ranging between -2 and -4 (Figure 6). Consequently, a process other than simple magma mixing must have played a role in the formation of the Nsuze and Mozaan, but also the Hlagothi and Kubuta volcanic rocks, i.e., assimilation-fractional crystallization.
(2) The Role of Assimilation-Fractional Crystallization. To constrain the possible crustal contaminants to the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta volcanic rocks, we calculated AFC paths using the Microsoft Excel calculation spreadsheet of Ersoy and Helvacı . For modeling the contamination history of Nsuze and Mozaan volcanic rocks, a komatiite from the study of Schneider et al.  was used as analogue for the parental komatiitic magma, having 26.9 wt.% MgO (sample ZA 33). Our approach follows the conclusion of Hegner et al.  who proposed a komatiitic parental magma for Nsuze volcanic rocks. For all other volcanic suites including the Kubuta samples of this study, we tried different mafic endmembers and chose a modeled 30% batch melt of incompatible element-depleted mantle (after Salters and Stracke ; equations after Shaw ; partition coefficients for mantle melting are listed in Table S 1), because it gave the best fit for our modelling. As possible contaminants, different crustal rocks from the AGC were tested, ranging in age from ca. 3.5 to 3.2 Ga (Figure 11).
The crustal compositions which yield AFC paths that best describe the evolved compositions of the Nsuze and Mozaan volcanic rocks are two granitoid gneisses of the 3.2 Ga crust-forming event, of which one has a tonalitic composition and was sampled in the vicinity of the Kubuta samples of this study (AGC 473, 3241 Ma) . The other one (AGC 449, 3229 Ma)  has a granodioritic to granitic composition and was sampled from the Piggs Peak inlier, the northernmost outcrop of AGC crust, close to the Phophonyane shear zone. Although this latter sample is not from the vicinity of exposed Nsuze or Mozaan volcanic rocks, we take its composition as an approximation of the composition of the crustal contaminant, because the two granitoid gneisses belong to the same crust-forming event at ca. 3.2 Ga, which is considered to have been widespread in the AGC (e.g., ). Based on Hf and Nd whole-rock isotope data, the ca. 3.2 Ga old AGC gneisses are interpreted to have formed by melting of mafic lower crust and variable assimilation of older ca. 3.5 Ga felsic crust, resulting in variable enrichments of incompatible elements in these rocks (e.g., ).
The chosen potential crustal contaminants are distinct in their concentrations of incompatible elements. Tonalite AGC 449 has high Th/Yb, high Nb/Yb, a low Zr concentration, and high La/Yb, whereas granodiorite AGC 473 has lower Th/Yb at elevated Nb/Yb, a higher Zr concentration, and high La/Yb . These two granitoid gneisses sample the heterogeneity of ca. 3.2 Ga AGC crust and represent two distinct end members in our modeling (Figure 11(a)). With this approach, two different AFC trends for Nsuze and Mozaan andesites and basaltic andesites in Th/Yb vs. Nb/Yb space can be reproduced. The upper AFC trend in Figure 11(a) (red line) is described by a maximum of 18% AFC of AGC 449, whereas the lower trend in Figure 11(a) (black line) can be reproduced by maximum 12% AFC for AGC 473. Again, the rhyolites of the Nsuze Group likely represent partial melts of continental crust. In contrast, in Zr vs. La/Yb space (Figure 11(c)), the AFC modeling does not reproduce all Nsuze and Mozaan lavas indicating that the chosen gneisses do not fully match the composition of the actual crustal contaminant. Additionally, the calculated percentages of AFC for Zr vs. La/Yb are higher than for Th/Yb vs. Nb/Yb which could be due to an overestimation possibly caused by small amounts of magma mixing, because Zr is more sensitive to magma mixing than Th.
On the other hand, our modeling for the Hlagothi, Usushwana, and Kubuta units shows that parental magmas derived from a depleted mantle source likely underwent AFC processes involving a crustal contaminant that has a low concentration of incompatible over more compatible elements and a high Zr concentration. In our modeling, this crustal contaminant is represented by the 3534 Ma tonalitic gneiss sample AGC 200 from the Ngwane Gneiss Unit (Figure 11(b)) . In Th/Yb vs. Nb/Yb space, our calculations show that the Hlagothi rocks and our Kubuta samples can be reproduced by maximum 10% AFC and subsequent fractional crystallization of the parental magma shifting the composition along a parallel line to the MORB-OIB array to higher Th/Yb and Nb/Yb (Figure 11(b)). In addition, in Zr vs. La/Yb space, the Hlagothi and Kubuta rocks can be reproduced by maximum 12% AFC and subsequent fractional crystallization as indicated by the arrow in Figure 11(d).
In summary, our calculations show that the variable compositions of the volcanic rocks from the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta suites are largely the result of assimilation-fractional crystallization of either a komatiitic parental magma (Nsuze and Mozaan) or a parental magma similar to depleted mantle melt (Hlagothi, Usushwana, and Kubuta). However, magma mixing between a fractionated melt that already underwent assimilation on the one hand and a felsic crustal melt on the other hand can also account for compositions of some of the Nsuze and Mozaan Groups lavas. Furthermore, our modeling shows that the AGC crustal rocks from both, the ca. 3.5 Ga and the ca. 3.2 Ga crust-forming events are suitable contaminants to the aforementioned volcanic sequences.
This observation, however, is inconsistent with the εNd(t) values of Nsuze, Usushwana, and Kubuta rocks that generally are about -2 to -3 (Figure 6). We calculated that at the time of formation of these volcanic units, i.e., at ca. 3.0 Ga, the ca. 3.2 AGC crust developed only to εNd(3.0) values of about -3 to -4, whereas the ca. 3.5 Ga AGC crust had already developed to εNd(3.0) values of up to -8 (Figure 6). Considering the low calculated amounts of AFC ranging from 10% to 18%, crustal assimilation of the younger generation TTGs seems to be unlikely. Hence, taking the isotope compositions into account, it is more likely that ca. 3.5 Ga old AGC crust was assimilated, that had comparable elemental compositions to the ca. 3.2 Ga old crustal rocks which were used for modeling.
5.3.2. Melting of Metasomatized Lithospheric Mantle
An alternative explanation for the enriched trace elemental composition in the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta lavas could be a lithospheric mantle source beneath the Eastern Kaapvaal Craton that had previously been overprinted metasomatically, similar to what has been proposed for the lithospheric mantle beneath the Pilbara Craton (e.g., ). As a consequence, partial melting of metasomatized lithospheric mantle or plume-interaction with such a metasomatized mantle would produce similar geochemical characteristics to crustal contamination, i.e., high Th/Nb and elevated LREE (e.g., [85, 86]), as observed in the aforementioned volcanic sequences.
The addition of fluids that caused metasomatism of the lithospheric mantle prior to magma formation could have occurred in two ways, either by ancient subduction processes during crust formation (e.g., [22, 29, 87]) or by a process analogous to sagduction such as proposed for the Pilbara Craton (e.g., [85, 86, 88, 89]). These two possibilities are controversially discussed and determining which of the two possibilities is more plausible for a possible metasomatic overprint of the lithospheric mantle beneath the Kaapvaal Craton is beyond the scope of this paper. Nevertheless, we cannot rule out the probability that partial melting of a metasomatized lithospheric mantle could have produced the observed enrichment in the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta volcanic rocks.
5.4. Mantle Source Composition
Another way to test the relationship of our Kubuta samples to the Nsuze Group or the supposed large igneous province consisting of Mozaan, Hlagothi, and Usushwana rocks (e.g., ) is to constrain the possible mantle sources for these units. The coherence in elemental anomalies of primitive mantle-normalized trace element patterns of our Kubuta samples compared to the Nsuze, Mozaan, Hlagothi, and Usushwana rocks (Figures 5(a) and 5(b)) suggests that all these volcanic rocks could be genetically related. We note that in Th/Yb vs. Nb/Yb space (Figure 9), the Kubuta samples, together with the samples from the other units, plot consistently above the modern MORB-OIB array towards higher Th/Yb. As discussed in the previous chapter, the high Th/Nb and LREE enrichment could either result from contamination of the parental magmas of these rocks with continental crust or, alternatively, from the presence of metasomatized lithospheric mantle beneath the underlying continental basement.
In previous studies on the volcanic rocks of the Nsuze Group, different types of parental magmas for the evolved rocks were suggested. For example, Armstrong et al.  proposed a basaltic parental magma for Nsuze volcanic rocks derived from a primitive mantle source comparable to a garnet-lherzolite, whereas Nhleko  suggested that the basaltic parental magma was derived from a mixture of depleted mantle and enriched mantle melts. In contrast, Hegner et al.  proposed that the Nsuze and Usushwana volcanic rocks were derived from an ultramafic komatiitic parental magma sourced from a depleted mantle reservoir.
In the previous section, we showed that the Kubuta, Nsuze, Mozaan, Hlagothi, and Usushwana volcanic rocks lie on assimilation-fractional crystallization trends involving different crustal components and at least two mantle source compositions. As a result of our modeling, the parental melts of the Nsuze and Mozaan volcanic rocks appear to be derived from komatiitic melts, similar to Barberton komatiites. These results support the hypotheses of Armstrong et al.  and Hegner et al.  for the origin of Nsuze volcanic rocks.
In contrast, all other samples from the Hlagothi, Usushwana, and Kubuta units do not appear to be derived from primitive mantle melts. Instead, these more mafic to ultramafic rocks seem to originate from a more depleted mantle source (Figure 11(a)). This feature indicates that despite the coherence of primitive mantle-normalized trace element patterns, the Kubuta samples together with the volcanic rocks from the Hlagothi and Usushwana units derived from a mantle reservoir distinct from that of the volcanic Nsuze and Mozaan Group rocks. This observation contradicts the conclusions of Hegner et al.  who proposed a komatiitic parental magma for Usushwana volcanic rocks, but supports the possible existence of a common mantle reservoir for the Hlagothi, Usushwana, and Kubuta units, as proposed by Gumsley et al. . However, in our modeling we use a modern depleted mid-ocean ridge basalt (MORB) mantle source  which might not be representative of the depleted mantle in the Archean. However, taken as a rough approximation, our results require an Archean mantle comprising reservoirs that were at least as depleted as modern MORB mantle. These could either be represented by ambient depleted mantle or, alternatively, depleted subcontinental lithospheric mantle. Furthermore, different mantle reservoirs for Hlagothi, Usushwana, and Kubuta volcanic rocks on the one hand, and Nsuze and Mozaan Group rocks on the other hand, would argue for a more complex formation history of the large igneous province than previously proposed by Gumsley et al.  based on age dating.
5.5. Petrogenetic Model for Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta Volcanic Rocks
Petrogenetic models for Nsuze, Mozaan, Hlagothi, and Usushwana volcanic rocks have been proposed earlier considering zircon and baddeleyite age dating, as well as common geochemical characteristics. These petrogenetic models invoke a plume origin for the aforementioned volcanic rocks with at least two pulses of magma generation (e.g., [4–6, 46]). In these models, the first pulse of magma produced the parental magmas of the Nsuze lavas and the Usushwana Igneous Complex, whereas a later mantle plume event produced the parental magmas of the Mozaan lavas and the Hlagothi Complex [4, 6]. Based on similar ages for the Usushwana and Nsuze suites and the more mafic composition of the Usushwana unit, Hegner et al.  and Gumsley et al.  interpreted the Usushwana dikes as being feeder dikes for the Nsuze volcanic rocks. Likewise, Gumsley et al.  suggested that the similar trace element patterns of the synchronous Hlagothi and Mozaan units together with the more mafic composition of Hlagothi rocks are indicative of the Hlagothi Complex being the feeding system for the Mozaan volcanism. Based on their observations, Gumsley et al.  and Gumsley et al.  proposed that the sills of the Hlagothi Complex and the grabbros of the Usushwana Igneous Complex together with the extrusive Mozaan lavas form part of a large igneous province in the southeastern part of the Kaapvaal Craton. This large igneous province is proposed to extend into the western Witwatersrand basin, based on synchronous volcanic sequences, dikes, sills, and layered complexes that are found on the eastern Witwatersrand block [4, 6]. A relation of synchronous volcanic rocks of the Mozaan Group and Witwatersrand Supergroup has also been suggested by Mukasa et al.  and recently confirmed by Paprika et al.  by a stratigraphic correlation based on zircon U-Pb ages with the Dominion group.
In the present study, we investigated the geochemistry and probable relation of Nsuze, Mozaan, Hlagothi, Usushwana, and the newly sampled Kubuta volcanic rocks in more detail. At first glance, similar trace element and Nd isotope compositions (where available) suggest a common mantle source and a common crustal contamination process for the aforementioned volcanic sequences. Our assimilation-fractional crystallization modeling of Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta volcanic rocks showed that the lavas from the Nsuze and Mozaan Groups were derived from a common parental melt similar to komatiite from a weakly depleted mantle source. In contrast, the volcanic rocks from the Hlagothi, Usushwana, and Kubuta suites were derived from a distinct, much more depleted mantle source.
Our findings challenge the hypothesis of Hegner et al. , Gumsley et al. , and Gumsley et al.  that the Usushwana Igneous Complex and Hlagothi Complex would represent feeder systems for Nsuze and Mozaan lavas, respectively. Instead, our results suggest that the Usushwana IC and the Hlagothi Complex were derived from a more incompatible element-depleted mantle source than the extrusive rocks of the Nsuze and Mozaan Groups. However, a synchronous formation of the Usushwana IC and Nsuze lavas on the one hand and the Hlagothi Complex and Mozaan lavas on the other hand is supported by precise age determinations [4–6]. Hence, we propose a refined petrogenetic model for Nsuze, Usushwana, Mozaan, Hlagothi, and the newly sampled Kubuta volcanic rocks (Figure 12) combining our findings and the conclusions of Hegner et al. , Gumsley et al. , and Gumsley et al. 
We propose that the komatiitic parental magma to the Nsuze lavas was sourced from an ascending mantle plume leading to the first pulse of magma generation at around 2980 to 2970 Ma. The plume likely stagnated beneath the lithospheric mantle and mafic underplating of the komatiitic magma occurred between the lithospheric mantle and the felsic TTG crust that ultimately sourced the komatiitic parental magmas to the Nsuze lavas. Contemporaneously, the heat of the mantle plume likely caused partial melting of incompatible element-depleted ambient mantle or depleted subcontinental lithospheric mantle possibly by conductive heating. These partial melts represent the parental melts of the Usushwana IC. Subsequently, both, the komatiitic and depleted mantle melts filled magma chambers in the felsic AGC crust and underwent AFC processes that reduced the density of the magmas and the evolved melts could erupt (Figure 12(a)). After this first magmatic event, thermal subsidence led to the deposition of Mozaan Group sediments in a shallow marine environment, likely an epicontinental sea [4, 18]. In a second event at around 2860 Ma, magmatism was reactivated and plume-sourced komatiitic parental melts produced the parental magmas of the Mozaan lavas. Analogous to the first magmatic event, the hot komatiitic plume caused partial melting of the incompatible element-depleted ambient mantle or, alternatively, the depleted subcontinental lithospheric mantle which resulted in the production of the parental magmas of the Hlagothi Complex. Again, these parental melts filled magma chambers in the felsic AGC crust and underwent AFC prior to eruption (Figure 12(b)).
In our refined petrogenetic model, the Kubuta rocks represent differentiation products of depleted mantle melts similar to the lavas of the Usushwana Igneous Complex or Hlagothi Complex. Due to the poorly constrained age of the Kubuta samples, it is currently not possible to distinguish if these rocks formed in association with the first or second pulse of magma generation.
The newly sampled metavolcanic rocks from the vicinity of the Kubuta Ranch in southern Eswatini resemble the major and trace element compositions of the Hlagothi and Usushwana Igneous Complexes, having similar trace element characteristics that indicate combined mineral fractionation and crustal assimilation. Also, our results show that the Nsuze and Mozaan volcanic rocks of the Pongola Supergroup were formed from a komatiitic parental magma. When compared to the Nsuze and Mozaan lavas, the Kubuta, Hlagothi, and Usushwana suites are more mafic to ultramafic in composition. In addition, assimilation and partial melting of continental crust is indicated in Nsuze, Usushwana, and Kubuta volcanic rocks by generally negative εNd(t) values and negative Nb-Ta anomalies.
Assimilation-fractional crystallization modeling shows that the lavas of all volcanic suites (Nsuze, Mozaan, Usushwana, Hlagoti, and Kubuta) result from modification of their respective parental magmas by AFC processes involving Ancient Gneiss Complex crust of diverse trace elemental composition. The heterogeneity of AGC crust is represented in the crustal end members used for modeling that include samples of the 3.2 Ga and 3.5 Ga crust-forming events. However, based on the negative εNd(t) values of about -3 in Nsuze, Usushwana, and Kubuta lavas and the calculated low amounts of AFC, we consider it more likely that AGC crust of the 3.5 Ga crust-forming event was assimilated. However, melting of a metasomatized lithospheric mantle or plume-lithospheric mantle interaction may also produce magmas with the observed geochemical characteristics.
Based on our results and previous studies, we propose a refined petrogenetic model for the formation of Nsuze, Mozaan, Hlagothi, Usushwana, and the newly sampled Kubuta volcanic rocks. In this model, we follow the interpretation of Gumsley et al.  that the aforementioned volcanic suits formed a Mesoarchean large igneous province that probably extended into the eastern Witwatersrand block, including the Dominion group . However, unlike Gumsley et al. , we also assign the Nsuze lavas to belong to this large igneous province. Furthermore, our results show that the Nsuze, Mozaan, Hlagothi, Usushwana, and Kubuta lavas were not derived from a single mantle plume, as proposed by Gumsley et al. , but rather were derived from two distinct sources—a plume producing komatiitic magmas (Nsuze and Mozaan) on the one hand and partial melting of depleted mantle (Hlagothi, Usushwana, and Kubuta) triggered by the conductive heat of the mantle plume on the other hand.
To test our hypothesis of two different mantle sources and, additionally, to constrain the contamination history in more detail to be able to distinguish if 3.2 Ga or 3.5 Ga AGC crust was preferentially assimilated, more geochemical work on the Mozaan, Hlagothi, and Usushwana and Dominion volcanic rocks is needed in form of major and trace element as well as whole-rock Hf-Nd isotope analyses. Furthermore, precise age determinations of the Kubuta volcanic rocks are crucial to affiliate these rocks to either the first or second pulse of magma generation.
Nevertheless, our findings strengthen the geochemical fingerprint of the proposed large igneous province  and provide information to better constrain the petrogenesis of these Mesoarchean flood basalts that record an important crustal growth process in the eastern Kaapvaal Craton at ca. 2.9 Ga.
All data are available in the manuscript and supplementary material. Raw data can be provided on request from the corresponding author.
Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.
JEH, KPS, CM, and AK thank the Deutsche Forschungsgemeinschaft for grants HO4794/1-1, HO4794/1-2, Mu1406/8, Mu1406/19-1, KR590/94-1, and KR590/94-2. Axel Hofmann is thanked for logistic support in the field and sample shipping. We thank Philipp Gleißner, Yogita Kadlag, and Chunhui Li (Berlin) for support with the ICP-MS as well as Frank Wombacher and Almut Katzemich for lab support at Cologne. Monika Feth is thanked for laboratory assistance at Berlin. This is a contribution within the framework of the SPP1833 Project ‘Building a habitable Earth’ funded by the DFG. We acknowledge support by the Open Access Publication Fund of the Freie Universität Berlin.