Abstract
A supercontinent is generally considered to reflect the assembly of all, or most, of the Earth's continental lithosphere. Previous studies have used geological, atmospheric and biogenic ‘geomarkers’ to supplement supercontinent identification. However, there is no formal definition of how much continental material is required to be assembled, or indeed which geomarkers need to be present. Pannotia is a hypothesized landmass that existed in the interval c. 0.65–0.54 Ga and was comprised of Gondwana, Laurentia, Baltica and possibly Siberia. Although Pannotia was considerably smaller than Pangaea (and also fleeting in its existence), the presence of geomarkers in the geological record support its identification as a supercontinent. Using 3D mantle convection models, we simulate the evolution of the mantle in response to the convergence leading to amalgamation of Rodinia and Pangaea. We then compare this supercontinent ‘fingerprint’ to Pannotian activity. For the first time, we show that Pannotian continental convergence could have generated a mantle signature in keeping with that of a simulated supercontinent. As a result, we posit that any formal identification of a supercontinent must take into consideration the thermal evolution of the mantle associated with convergence leading to continental amalgamation, rather than simply the size of the connected continental landmass.
The formation of supercontinent Pangaea (Fig. 1a) had a profound effect on the thermal evolution of the mantle. An increase in deep mantle upwelling and plume formation is thought to follow supercontinent formation (Gurnis 1988; Lowman and Jarvis 1993; Zhong and Gurnis 1993; Zhong et al. 2007; Li and Zhong 2009; Yoshida 2010; Heron and Lowman 2010; Yoshida 2013; Gamal El Dien et al. 2019), with large igneous provinces (LIPs) being the surface manifestations of such plumes (e.g. Yale and Carpenter 1998; Courtillot et al. 1999; Ernst et al. 2005; Ernst and Bleeker 2010; Sobolev et al. 2011). Additional ‘geomarkers’ (see Nance et al. 1986) have also been identified for supercontinent Pangaea, including global-scale orogenesis (Nance et al. 1988; Santosh 2010a; Condie 2011; Müller et al. 2013), crustal growth (Hawkesworth et al. 2010, 2016), the distribution of metamorphic belts (Brown 2007) and mineral deposits (Goldfarb et al. 2010), rapid climate swings (Hoffman et al. 1998; Strand 2012; Young 2012), major events in the evolution of life and the atmosphere (Lindsay and Brasier 2002; Santosh 2010b; Knoll 2013; Melezhik et al. 2013), biogeochemical cycles (Nance et al. 1986), and profound sea-level change (Worsley et al. 1984; Nance and Murphy 2013). Such geomarkers are interpreted as signals of supercontinent convergence leading to amalgamation and may be useful proxies for identifying the existence of previous supercontinents for which palaeocontinental reconstructions may be less refined (Nance and Murphy 2019).
Continental configurations at (a) 320 Ma (Matthews et al. 2016) for Pangaea formation; (b) 900 Ma (Merdith et al. 2017) for Rodinia formation; (c) 700 Ma (Merdith et al. 2017) for Pannotia amalgamation ; (d) 650 Ma (Merdith et al. 2017) for Pannotia formation and change in subduction pattern; (e) 600 Ma (Merdith et al. 2017) showing Pannotia formation and region of CIMP; (f) 550 Ma (Merdith et al. 2017) for Pannotia dispersal. S, Siberia; L, Laurentia; B, Baltica; A, Amazonia; CIMP, Central Iapetus Magmatic Province. In (c–e) convergence zones leading to the amalgamation of the Gondwanan portion of Pannotia (blue lines) correspond with the approximate location of subduction zones.
Although their ages of amalgamation and break-up are still being refined, a number of pre-Pangaean supercontinents have been proposed (based on the identification of geomarkers in the geological record). Rodinia (Fig. 1b), with Laurentia at its centre, was formed by 1.1–0.9 Ga global-scale orogenesis (Dalziel 1991; Hoffman 1991; Moores 1991; Torsvik 2003; Li et al. 2004, 2008). It was fully assembled by c. 900 Ma and is believed to have rifted apart in two separate episodes (c. 0.85–0.70 Ga and c. 0.62–0.54 Ga) (Li et al. 2008, 2013; Li and Evans 2011). Earlier supercontinents include Nuna/Columbia (Hoffman 1997; Rogers and Santosh 2002) whose existence is most recently dated at c. 1.6–1.4 Ga (e.g. Pehrsson et al. 2015), Kenorland (Williams et al. 1991) thought to have existed during the interval c. 2.7–2.5 Ga, and possibly Ur (c. 3.0 Ga), although the latter is better described as a supercraton (i.e. transient, late Archean landmasses which broke up to form cratons, e.g. Bleeker 2003).
Although there are geological markers for supercontinent assembly (e.g. Nance et al. 2014), no formal system exists for supercontinent identification, as a result of which there is no formal definition of a supercontinent (e.g. Meert 2012). Even a basic understanding of how much contiguous landmass is required is open to debate. For example, Meert (2012) proposed a figure of 75% of available continent crust but, as pointed out by Nance and Murphy (2019), this is both difficult to estimate and views supercontinents as objects rather than the consequence of a geodynamic process.
The status of Pannotia (c. 0.65–0.54 Ga; Gondwana–Laurentia–Baltica and possibly Siberia, Fig. 1c–f), is an example of the ongoing debate about what constitutes a supercontinent (see Murphy et al. 2020; Evans, in press, this volume). Given the uncertainties inherent in Neoproterozoic palaeogeographical reconstructions, Nance and Murphy (2019) argue that Pannotia's status as a supercontinent is apparent from a number of ‘supercontinent’ markers around the time of its proposed formation (e.g. global-scale orogeny, rapid continental growth, profound changes in the chemistry of the oceans and atmosphere, rapid and dramatic climate swings, as well as an explosion in biological activity; Nance et al. 1986; Hoffman 1991; Hoffman et al. 1998; Maruyama and Santosh 2008; Knoll 2013). However, there is general agreement that Pannotia's size and duration were significantly less than that of Pangaea, drawing controversy as to whether it can be classified as a unique supercontinent or was simply part of Rodinia after initial break-up (e.g. Mitchell et al. 2012; Li et al. 2019).
An area of Pannotian dynamics that has not been fully explored is the response of the mantle to the continental convergence leading to its amalgamation. Late Neoproterozoic continental convergence and protracted subduction leading to the assembly of the Gondwanan portion of Pannotia is strongly indicated by widespread ‘Pan-African’ collisional orogenic activity between 0.62 and 0.53 Ga (e.g. Li et al. 2008; Merdith et al. 2017). Investigating the effect of continental convergence and assembly on Late Neoproterozoic–Early Paleozoic (0.6–0.5 Ga) mantle convection patterns, could allow for supercontinent identification through its impact on mantle dynamics (e.g. Heron 2019). Furthermore, as the legacy of such an impact would affect Paleozoic mantle convection patterns, it would also need to be factored into models purporting to explain the amalgamation of Pangaea.
The Central Iapetus Magmatic Province (CIMP, Fig. 1e, f) is one potential manifestation of changes in deep mantle circulation that may have been an outcome of continental convergence (in a similar way to the Central Atlantic Magmatic Province and the convergence of Pangaea). CIMP comprises several LIPs and predominantly consists of mafic dyke swarms, basalts, bimodal basalt–rhyolite complexes, lamprophyres and carbonatites (Kamo et al. 1995; Higgins and van Breemen 1998; Cawood et al. 2001; Puffer 2002; Keppie et al. 2006; Pisarevsky et al. 2008; Ernst and Bell 2010; Thomas et al. 2012; Ernst 2014). CIMP initiated along a rift–rift–rift triple junction between Laurentia, Baltica and Amazonia (Fig. 1e, f) with pulses of magmatism (peaking at c. 550 Ma) widely attributed to plumes that herald the opening of the Iapetus Ocean and Tornquist Sea (Cawood et al. 2001; Pisarevsky et al. 2008; Thomas et al. 2012).
In this study, we use 3D thermo-chemical global mantle convection experiments that take into account the palaeo-subduction history from Rodinia to Pangaea (Matthews et al. 2016; Merdith et al. 2017). We analyse the thermal response of the mantle driven by the evolution of surface tectonics to identify any mantle signatures generated by convergence leading to Pangaea and Rodinia assembly, and then apply these to Pannotia. The models are checked against geological markers such as plume and hotspot locations to gauge the influence of the convergence of continental material in influencing mantle dynamics.
Modelling global mantle convection
We model global mantle convection using the ASPECT code (Kronbichler et al. 2012; Heister et al. 2017). ASPECT solves the following set of equations for compressible convection in the Earth's mantle, describing the mass, momentum and energy balance (taking into account adiabatic heating, shear heating and radiogenic heat production), and the transport of chemical composition:
The equations are solved for velocity u, pressure p, temperature T and chemical composition C (e.g. harzburgite or recycled ocean crust). η is the viscosity,
Our models begin from an initial condition that has an average mantle temperature gradient and composition that is appropriate for the present-day mantle (Fig. 2d) and in keeping with previous studies (e.g. Dannberg and Gassmöller 2018; Zhang and Li 2018). The initial mantle temperature is laterally homogeneous, following an adiabatic profile with thermal boundary layers at the top and bottom (Dannberg and Gassmöller 2018). We fix the surface temperature to 273 K and core–mantle boundary temperature to 3700 K (Fig. 2d). Although 3700 K is a higher temperature than some previous studies (e.g. Zhang and Li 2018), there is uncertainty in the exact value for the present day (e.g. Nomura et al. 2014). There is greater uncertainty in historical core-mantle temperature, with previous studies indicating a higher basal temperature in the Hadean (e.g. 4000 K, O'Neill and Zhang 2018). Here, by fixing the core–mantle boundary temperature throughout the study (from 1 Ga to 0 Ma), we do not capture the evolution of the interior average mantle temperature. However, the change in temperature between 1 Ga and the present day will not be significant given the uncertainty. As a result, a choice of 3700 K is appropriate for all geological times in the simulation.
The model setup. Global mantle convection models start from a homogenous initial condition (a) and are driven by plate velocities imposed at the surface (b) to produce subduction and (return-flow) plumes (c). The depth profiles of temperature (d), viscosity (e), density (f) and thermal expansivity (g) are given for the initial condition of a thermo-chemical model. The excess temperature is the difference between the temperature and adiabatic profiles (e.g. the non-adiabatic temperature).
Viscosity in the simulations (Fig. 2e) is depth- and temperature-dependent (as outlined in Heister et al. 2017) and based on a published viscosity model incorporating constraints from mineral physics, geoid deformation and seismic tomography (Steinberger and Calderwood 2006). Our model produces viscosity variations of four orders of magnitude over a temperature range appropriate for the Earth (ΔT = 3700 K) (Fig. 2d). However, this viscosity range is manually limited to reduce the computation time and make model runs feasible. In the Earth, lateral viscosity contrasts are likely to be larger than in our models. As in previous studies (Dannberg and Gassmöller 2018), we assume an average mantle composition of 82% harzburgite and 18% recycled oceanic crust (Xu et al. 2008) and compute the material properties ρ, α and Cp using PerpleX (Connolly 2009) and a mineral physics database (Stixrude and Lithgow-Bertelloni 2011). Accordingly, all material properties include the effects of phase transitions, with thermal expansivity and specific heat taking into account the corresponding latent heat release and consumption (Nakagawa et al. 2009).
Although we do not prescribe continental and oceanic material present at the surface, the model uses temperature-dependent viscosity which form the outlines of plates bound by the ridges and subduction zones produced from the applied plate reconstruction velocities. In that respect, we can generate areas of reduced heat flux due to cold material thickening in response to plate tectonic processes. As a result, our models can act to produce areas of the surface that could be described as a ‘continent’ and ‘ocean’ (e.g. Fig. 3).
Output of Model PangT for the present day (0 Ma). (a) A view into the mantle beneath the northern Atlantic and Europe for Model PangT. Excess mantle temperature anomaly contoured surfaces are given for warm (300 K, red) and cold ( − 500 K, blue) regions. (b) Temperature at 1 km depth from the surface for the model, highlighting the present-day ridges and continents (continental outline shown by white lines, Matthews et al. 2016).
For this study, we used ASPECT version 2.0.0-pre (commit 572f967) but compared models with a more recent version of the code. This was done to ensure planned methodological changes to a new version of the global mantle convection setup used here did not impact the conclusions drawn from our (older) simulations. All model parameters are given in Table 1, and the initial radial profiles of temperature- and depth-dependent material properties are plotted in Figure 2.
Model material parameters
The boundary condition at the surface is provided by plate reconstructions (Fig. 2b) and is therefore kinematic (e.g. it does not change in accordance with the force balance of the underlying mantle) (Fig. 2c). In this study we use two different plate reconstructions which produce (using GPlates 2.0) surface velocities at 1 myr intervals at a resolution of 1° × 1°. For our set of Pangaea models, we use the Matthews et al. (2016) plate reconstructions, which extend back in time 410 million years. The Matthews et al. (2016) global model includes the Paleozoic model of Domeier and Torsvik (2014) and the Mesozoic–Cenozoic model of Müller et al. (2016). For our set of Rodinia/Pannotia models, we use Merdith et al. (2017) plate reconstructions that cover the interval 1000 to 520 Ma (480 myr in total). Although these reconstructions do not yield a coherent Pannotia supercontinent, we utilize the surface velocity fields during the convergence that culminated in the Pan-African collisional orogenies and the assembly of Gondwana. At present, a plate reconstruction history between 0.52 and 0.41 Ga is not available, therefore we could not use a single simulation from 1 Ga to the present day.
We apply a 190 myr ‘spin-up’ period to the start of all models to reduce the impact of the initial condition on mantle dynamics. In starting our models from a spherically symmetrical setup (e.g. Fig. 2), there is a period where the mantle overturns – generating a pulse in core–mantle boundary heat flux. By applying a ‘spin-up’ phase of 190 myr, we ensure no signal from the initialization of the model affects our results. The time 190 myr was chosen as it is long enough that the mantle is well mixed but not so long that the models become ‘expensive’ to run on a supercomputer (i.e. slow). In model development, 50, 100 and 190 myr spin-up phases were tested – the overall trends of the heat flux in the models were similar in all cases ensuring the results were robust. The simulations last more than 600 myr in total, but our analysis begins after 190 myr when the mantle is well mixed. Further details of the model setup can be found in Dannberg and Gassmöller (2018).
In the following sections we describe the mantle dynamics that occur for two separate models: Model PangT, which uses 410 myr of surface velocities from the Matthews et al. (2016) plate reconstruction; and Model NeoproT, which applies surface velocities from the Merdith et al. (2017) plate reconstruction for 1000 to 520 Ma. Both models are purely thermal simulations, meaning there are no thermo-chemical piles present. We first look at the behaviour of the mantle in response to the formation and dispersal of the supercontinent Pangaea (e.g. Model PangT), and then apply these responses to Pannotia to see whether this putative supercontinent can be identified from the evolution of the Neoproterozoic mantle (Model NeoproT).
The response of the mantle to supercontinent Pangaea
Studies analysing the timing of continent–continent collisions and age of rifting sequences show that the land masses of Gondwana (the African, Antarctic, India, Australian and South American plates) and Laurasia (the Eurasian and North American plates) formed the supercontinent Pangaea near the equator at c. 320 Ma (Smith et al. 1981; Hoffman 1991; Scotese 2001). The break-up of Pangaea is thought to have occurred in a number of stages with North America separating from the land mass at c. 180 Ma (starting the opening of the Central Atlantic Ocean), the dispersal of the Antarctic–Australian, Indian and South American continents occurring between 140 and 100 Ma (Smith et al. 1981; Scotese 2001) and, finally, Australia separating from Antarctica during the Paleocene (Veevers and McElhinny 1976).
The surface velocities from Matthews et al. (2016) span 410 myr (410–0 Ma) and therefore encapsulate the formation and breakup of Pangaea. Although there is no continental material specified at the surface, our models apply appropriate subduction zone positions over time dictated by the plate reconstruction history. Figure 3 shows the output for Model PangT at 0 Ma (present day), with plumes rising from the core–mantle boundary to the surface alongside cold downwelling slabs interacting with the mantle flow (Fig. 3a). The temperature close to the surface shows the presence of ridges in appropriate positions in the mid-Atlantic, east Pacific, and Indian oceans for the present day (Fig. 3b).
Figure 4 shows the positioning of plumes over time for Model PangT. Plumes are calculated as areas with an excess temperature greater than 400 K at 700 km depth (shown in orange on Fig. 4). The hottest regions of plume activity are shown in red (>550 K). In our model, there are more plumes forming in regions where continental material would be located during supercontinent formation (7 red plumes in Fig. 4d) than there are during dispersal (3–4 in Fig. 4a, b). At 40 Ma in Model PangT (Fig. 4b), plumes form in the location of Iceland and the High Arctic, consistent with geological evidence (e.g. White and McKenzie 1989; Jones et al. 2002; Rickers et al. 2013). Likewise, during continental dispersal, a number of plumes are positioned along the mid-Atlantic ridge. There is also a distinct change in mantle dynamics from Pangaea formation to dispersal with plumes appearing to be mostly ‘hot’ (>550 K) (Fig. 4).
Evolution of plume positions for Model PangT for the past 340 million years. (a) present day, (b) 40 Ma, (c) 220 Ma and (d) 340 Ma. Figures show continental positions as given by Matthews et al. (2016) and the plume positions at 700 km depth represented by excess temperature values of >550 K (red) and >400 K (orange). Present-day hotspot locations from French and Romanowicz (2015) and Steinberger (2000) are shown in yellow circles in (a).
At 0 Ma (Fig. 4a), there is some agreement (and some disagreement) between the location of model plumes and the present-day hotspot database (Steinberger 2000; French and Romanowicz 2015). In particular, plumes form near the hotspots of Iceland, Hawaii, Samoa, Louisville, San Fernandez, Tristan, the Azores, Cameroon, and the East African Rift System. In testing the robustness of the simulation, we formally analysed the positioning of the hot thermal anomalies. For half the database Model PangT produces a plume within 1500 km of a present-day hotspot. This is a performance that is better than 80% of 10 000 simulations that randomly position the same number of plumes as Model Pang T (45) and then compare these to the 29 hotspots in the database.
By analysing the core–mantle boundary heat flux (Fig. 5a), we find a signature in the response of the mantle to the formation of the supercontinent Pangaea: the basal heat flux (e.g. core–mantle boundary) increases from 16 to 18 terawatt (TW) while Pangaea is contiguous (e.g. 320–180 Ma, Fig. 1a) and remains relatively stable (c. 18 TW) during dispersal (c. 180 Ma). It should be noted that the actual values given in our simulations for heat flux are strongly related to assumptions regarding lower mantle rheology. For example, changing the core–mantle temperature or thermal conductivity in the model would strongly influence the actual heat flux values – but not the overall trend.
Thermal response of the mantle to supercontinent formation and dispersal. (a) and (c) show the basal heat flux for Model PangT and Model NeoproT, respectively. Panels (b) and (d) show the evolution of plume coverage over time for the two models. Plume coverage is defined by the number of 1° × 1° boxes that have an excess temperature value >400 K at 700 km depth (e.g. a hot thermal anomaly). For example, at 400 Ma in (b), Model PangT has approximately 150 1° × 1° boxes that feature a ‘plume’ (which represents 0.5% of the potential coverage available). Times when a supercontinent may have been in existence are given by the shaded regions (Smith et al. 1981; Hoffman 1991; Scotese 2001; Li et al. 2008, 2013; Li and Evans 2011; Nance and Murphy 2019). The significant difference between the basal heat flux and plume coverage values for the end of Model NeoproT (520 Ma) as compared to the beginning of Model PangT (410 Ma) is an artefact of the two models using separate (non-continuous) plate reconstructions and beginning from different initial conditions.
Our models infer that convergence leading to Pangaea amalgamation may have been sufficient by 400 Ma to impact core–mantle boundary heat flux, as shown by the increase in Figure 5a. This timing corresponds to the 100 Ma record of convergence leading to the collision of promontories of Gondwana with Laurussia, a major phase in the amalgamation of Pangaea (e.g. Kroner and Romer 2013; Arenas et al. 2014; Wu et al. 2020). Since our global mantle convection Model PangT performs well against plume location tests (Fig. 4), we can apply the model's mantle response to Pangaea (Fig. 5a) to Neoproterozoic time to investigate whether Rodinia and/or Pannotia display a similar signature.
Pannotia's mantle signature?
Figures 5c, d, 6 and 7 show the response of the mantle to imposed surface velocities spanning 1000 to 520 Ma (Merdith et al. 2017). Based on the criteria specified in the previous section, we can identify the assembly and dispersal of Rodinia from the core–mantle boundary heat flux. Figure 5c shows the increase in heat flux during the assembly and lifespan of Rodinia – an increase which stops during its break-up and dispersal. A similar pattern of increasing basal heat flux is produced in our model between c. 660 and 570 Ma (Fig. 5c), which corresponds to continental convergence leading to the Pan-African collisions during which Pannotia is thought to have assembled. This increase stops post-570 Ma, with our model producing an increase in plume coverage (Fig. 5d). This timing corresponds with peak activity of CIMP. Some of the plumes form under the site of the putative supercontinent Pannotia at 550 Ma (Figs 6 & 7a).
Evolution of plume positions for Model NeoproT for the interval 900–550 Ma. (a) 550 Ma, (b) 640 Ma, (c) 740 Ma and (d) 900 Ma. Figures show continental and convergence locations as given by Merdith et al. (2017) and the plume positions at 700 km depth represented by excess temperature values of >550 K (red) and >400 K (orange). Dashed purple lines in (a) outline region displayed in Figure 7. Marker A is used as a reference for plume location in Figure 7a.
Output of Model NeoproT post-Pannotia (550 Ma). (a) A view into the mantle in the area outlined in Figure 6a. Excess mantle temperature anomaly contoured surfaces are given for warm (300 K, red) and cold (− 500 K, blue) regions. Marker A is in reference to the plume location in Figure 6a. (b) Temperature at 1 km depth from the surface for the model, highlighting the position of ridges and continents (continental outline at 550 Ma from Merdith et al. 2017 shown in white).
The potential role of large low shear velocity provinces
Large structures of low shear-wave seismic velocities (large low shear velocity provinces or LLSVPs) observed in the lower mantle by seismic tomography studies and interpreted as regions of elevated temperatures and/or compositional anomaly, appear to play a defining role in the origin and positioning of mantle plumes (Steinberger 2000; Torsvik et al. 2006; Zhong et al. 2007; Deschamps et al. 2012; French and Romanowicz 2015; Cottaar and Lekic 2016). Although the composition and origin of LLSVPs are debated (e.g. Davies et al. 2012; Koelemeijer et al. 2017; Zaroli et al. 2017), there is growing consensus (Tan et al. 2002; Torsvik et al. 2006; Garnero et al. 2016) that plumes may form within or adjacent to these regions. In turn, the formation of plumes can influence supercontinent dispersal, leading to a further repositioning of subduction zones and a change in the pattern of mantle flow (van Hinsbergen et al. 2011). The relative importance of LLSVPs in controlling the positioning of mantle plumes is also controversial. The end-members of this debate are (i) that LLSVPs are spatially stable antipodal structures over time, with plumes generated from their edges (e.g. Torsvik et al. 2006, 2010), or that (ii) LLSVPs are entirely passive in mantle convection processes controlled by subduction flow (e.g. McNamara and Zhong 2005; Zhong et al. 2007; Davies et al. 2012; Heron et al. 2015). The uncertainty and ambiguity in whether these geophysical features signify a dense, chemically distinct layer in the deep mantle (e.g. a thermo-chemical pile) or represent purely thermal anomalies make modelling mantle convection difficult (for reviews see Garnero et al. 2016; McNamara 2019).
The previous two models (PangT and NeoproT) do not feature any thermo-chemical layer. However, in Model NeoproTC we apply a 150 km thick thermo-chemical layer to the initial condition of the model (Fig. 2g). Following a number of previous studies (e.g. Heister et al. 2017; Dannberg and Gassmöller 2018), our model LLSVPs have the composition of mid-ocean ridge basalt and are therefore chemically distinct from the surrounding mantle. Here, the thermo-chemical piles have a corresponding buoyancy number
Figure 8 shows the evolution of the thermo-chemical pile in Model NeoproTC from an initial uniform layer. The yellow colours outline the shape of the model LLSVP at 550 Ma with material located mainly in the present-day positions of the western Pacific Ocean, the Atlantic Ocean and Europe. Model plumes, given by the red markers, are generated mainly from the edges of the thermo-chemical piles rather than their interior (Fig. 8a). However, some plumes form at a distance from the thermo-chemical pile.
Output of Model NeoproTC post-Pannotia (550 Ma). (a) Excess temperature of Model NeoproTC near the core-mantle boundary (2800 km depth) showing the outline of the thermo-chemical pile (yellow) alongside model plume positions (red and orange) and paleo-continent locations (black outline, Merdith et al. 2017). The white contour corresponds to the present-day positions of large low shear velocity provinces. This low velocity contour is generated from runs where 9 out of 14 seismic tomography models agree that there are slow anomalies (Shephard et al. 2017). Numbers 1 to 7 indicate model deep mantle plumes that originate close to the edge of present-day LLSVPs (white). Projection is centred on Greenwich Meridian. (b) Basal heat flux and (c) evolution of plume coverage for Model NeoproTC (purple) and Model NeoproT (green). Due to the presence of a thermo-chemical pile, there is a larger plume coverage in Model NeoproTC than NeoproT. Numbers 1 to 7 indicate model deep mantle plumes that originate close to the edge of present-day LLSVPs (white).
In comparison to the purely thermal model (NeoproT, Fig. 6a), the thermo-chemical simulation produces fewer hot plumes (red) under the oceanic lithosphere, although plumes still develop beneath similar continental locations (Fig. 8a). This simulation suggests that the overall pattern of plume formation is not significantly affected by the presence of a thermo-chemical layer (e.g. Hassan et al. 2015). Furthermore, the presence of a thermo-chemical layer acts as an insulator to the core, reducing the heat flow through the boundary (Fig. 8b).
In analysing the criteria for identifying a supercontinent, the basal heat flux in the presence of a thermo-chemical layer masks the dispersal of Rodinia and the formation of Pannotia. In Model NeoproTC, the core–mantle boundary heat flux continues to increase from 750 to 650 Ma (Fig. 8b), but in Model NeoproT it stabilized over the same time interval (Fig. 5c). However, the impact of the convergence leading to Pannotia assembly can still be identified at 650 Ma through the increased rate of heat leaving the core into the mantle. In terms of plume coverage, the evolution of Model NeoproTC (Fig. 8c) is similar to that of NeoproT (Fig. 5d): while Rodinia is being assembled there is an increase in plume coverage that stops following Rodinia break-up, but increases once more during the continental convergence leading to Pannotia. Due to the presence of a thermo-chemical pile, there is a larger plume coverage in Model NeoproTC than NeoproT.
A number of previous studies have proposed that LLSVPs may have been relatively fixed in their location on timescales greater than that of supercontinent formation and dispersal (Torsvik et al. 2006, 2010). Although spatially fixed thermo-chemical piles are not modelled in this study, we show the outline of present-day LLSVPs in Figure 8a to highlight where (potentially) a thermo-chemical pile would have been located with respect to Pannotia's continental configuration. The break-up of Pannotia could have been facilitated by mantle plumes forming at the edges of fixed LLSVPs, since Baltica, Siberia and southwestern Gondwana are situated above the margins of the present-day seismic anomalies. However, a number of plumes in both Model NeoproTC and NeoproT do form at what would be the edge of a present-day LLSVP location (e.g. markers 1–7, Fig. 8a), although no thermo-chemical piles are situated in the present-day positions.
Discussion
We developed 3D global mantle convection models that feature surface velocities linked to plate reconstructions to identify a mantle signature of supercontinent convergence, assembly and dispersal. By analysing the response of the mantle to the convergence, assembly and dispersal of Pangaea (410–0 Ma), we were able to identify such a signature, which we then compared with that of Rodinia and with convergence leading to assembly of hypothesized supercontinent Pannotia (620–540 Ma). For Pangaea and Rodinia, the models produce an increase in heat flux across the core–mantle boundary during the lifespan of the supercontinent (Fig. 5), and a stabilization in basal heat flux following break-up. For convergence leading to Pannotia assembly, the models produced a similar thermal response at the base of the mantle over similar timescales (Fig. 5).
These supercontinent ‘markers’ of core–mantle boundary heat flux and sub-continental plume formation (i.e. deep mantle convection patterns) are difficult to quantify in a geological context. However, there is a potential link between geomagnetic variations and mantle convection (e.g. Davies 1999; Buffett 2002). Simulations of the geodynamo suggest that transitions from periods of rapid polarity reversals to periods of prolonged stability (e.g. geomagnetic superchrons) may have been triggered by changes in core–mantle boundary heat flow either globally or in equatorial regions – a response that has been previously applied to the supercontinent cycle (e.g. Olson et al. 2010, 2015; Zhang and Zhong 2011; Biggin et al. 2012; Choblet et al. 2016). Our results are in agreement with previous studies that have applied similar methodologies. Zhang and Zhong (2011) used 3D numerical simulations and 450 myr of plate reconstruction history to analyse core–mantle boundary heat flux during supercontinent formation and dispersal. Despite generating very similar trends to this study for Pangaean heat flux with and without thermo-chemical piles, Zhang and Zhong (2011) highlighted the fluctuating heat flux between 100 Ma and the present day as being most important. The trend of increasing heat flux between 350 and 250 Ma (Fig. 5) that we present here is also shown in Zhang and Zhong (2011).
Our study marks a first step in the classification of periods of continental convergence and amalgamation over geological time through identifying mantle thermal responses to the convergence, assembly and dispersal stages of supercontinents. We do not model continental material which, in itself, has an impact on mantle dynamics through modifying the wavelength of mantle flow (Grigné et al. 2007; Phillips et al. 2009; Rolf et al. 2012) and the initiation and polarity of subduction (e.g. Crameri and Tackley 2014), alongside topographic and sea-level influences. However, we do prescribe surface velocities from plate reconstructions (Matthews et al. 2016; Merdith et al. 2017) that may help to artificially capture some aspects of these mantle dynamics. As a result, it would be possible to compare times of enhanced basal heat flux with the timing and frequency of orogenesis and potential sea-level variations, which are characteristics of supercontinent assembly (e.g. Nance and Murphy 2019). A study of the impact of our mantle convection models on crustal growth, climate, the evolution of life and biogeochemical cycles (i.e. other potential factors in supercontinent identification) would also be possible but is beyond the scope of this study. Although the models we present here do not take into consideration crustal growth (e.g. Hawkesworth et al. 2010; Rozel et al. 2017; Jain et al. 2019a, b), we do analyse the amount of plume coverage arriving from the deep mantle (e.g. Fig. 5). Such upwellings may help in understanding mantle circulation patterns over time when compared to the LIP record (Ernst et al. 2005; Ernst and Bleeker 2010).
The CIMP initiated along a rift–rift–rift triple junction between Laurentia, Baltica and Gondwana (Amazonia) and most reconstructions for the early stages of the opening of the Iapetus Ocean and Tornquist Sea invoke a plume centered on this triple junction (Higgins and van Breemen 1998; Cawood et al. 2001; Pisarevsky et al. 2008). In our numerical models a plume (and/or series of plumes) develops in the region where CIMP is thought to have originated (e.g. marker A in Figs 6a & 7a for Model NeoproT, and markers 1–3 in Fig. 8a for Model NeoproTC). Although we do not produce a suite of models that systematically test the potential density and rheology of the thermo-chemical layer at the base of the mantle, it is important to note that ‘CIMP’ plumes occur regardless of whether a model ‘LLSVP’ features in the models or not. This outcome indicates that the plate reconstruction history may play a more dominant role in mantle dynamics than the composition of the deep mantle (i.e. top-down over bottom-up), as previously proposed by McNamara and Zhong (2005) and Davies et al. (2012). The CIMP identification in our study can be used as an indication that the models are performing well. Future work to further analyse model LIPs during the timeframe of Pannotia will help us understand how robust such numerical models are in producing geological events.
The impact of continental insulation on the thermal evolution of our models cannot be assessed because the models do not include continental material. Nevertheless, some studies indicate that continental insulation is a potential driving mechanism for the supercontinent cycle. Given that Pangaea was relatively stable and largely stationary for >100 million years, Anderson (1982) proposed that continental insulation could have had a dramatic effect on the underlying mantle. Over this length of time, the excess heat trapped by the supercontinent would cause uplift (through thermal expansion), partial melting of the mantle and, ultimately, the dispersal of continental material (Anderson 1982; Worsley et al. 1984; Nance et al. 1986, 1988; Gurnis 1988). In addition, a large geoid high similar to that over present-day Africa would be generated sub-supercontinent as a result of the thermal expansion caused by the continental insulation. Subsequent numerical (e.g. Lowman and Jarvis 1993, 1999; Zhong and Gurnis 1993; Yoshida et al. 1999; Lenardic et al. 2005, 2011; Phillips and Bunge 2005; Coltice et al. 2007, 2009; Phillips et al. 2009; Phillips and Coltice 2010; Yoshida 2010; Rolf et al. 2012) and geochemical (Brandl et al. 2013; Brown and Johnson 2018) studies suggest an important role for continental insulation in producing increased sub-continental mantle temperatures.
Although the models presented here do not take into consideration continental insulation, there remains a thermal response to the formation of a supercontinent that generates elevated sub-continental temperatures (e.g. Fig. 5). However, we find that subduction history is the most important factor in driving mantle convection patterns; an observation that follows a number of previous studies that suggest a diminished role of continental insulation in global mantle dynamics (Lowman and Jarvis 1995, 1996; Heron and Lowman 2010, 2011, 2014; Yoshida 2013; Heron et al. 2015).
The results presented here represent a first step in understanding global mantle convection during the late Neoproterozoic. A comparison between different numerical methodologies would be beneficial to test the robustness of our findings. A better understanding of the role of the initial thermal condition in the generation of a mantle ‘signal’ in response to a supercontinent is also needed. Furthermore, we only use one interpretation of what the plate reconstruction history could have been during the Neoproterozoic. For instance, Scotese and Elling (2017) present reconstruction for the time period with a contiguous Pannotia. Given the difficulty in generating plate reconstructions at any given time period (let alone a billion years ago), we acknowledge the uncertainty surrounding the plate velocities used in our mantle convection models. Further work in this field must be used to contrast and compare the impact of different plate reconstructions on mantle convection models.
The availability of a plate reconstruction history between 0.52 and 0.41 Ga could be used to produce a single simulation covering a billion years, rather than the two separate models for 1 to 0.52 Ga (Merdith et al. 2017) and 0.41 Ga to the present day (Matthews et al. 2016) that we implement. Future work that incorporates a full plate reconstruction history from 1 Ga to the present day would allow for a better understanding of the transition between supercontinent dispersal and formation, in particular for Pangaea. In the models presented here, there is a significant difference between the basal heat flux and plume coverage values for the end of Model NeoproT (520 Ma, Fig. 5c, d) as compared to the beginning of Model PangT (410 Ma, Fig. 5a, b). This difference is an artefact of the two models using separate (non-continuous) plate reconstructions and beginning from different initial conditions.
Our models indicate that the convergence, assembly and dispersal of Pannotia had a mantle response that we can measure. However, it is unlikely that the continental extent of Pannotia would fulfil the criterion that a supercontinent should consist of at least 75% of the preserved continental crust prior to initial breakup (Meert 2012). To test further the robustness of our models, it would be appropriate to compare the thermal evolution of the mantle from numerical simulations that do not feature prescribed plate surface velocities (e.g. Phillips et al. 2009; Phillips and Coltice 2010; Rolf and Tackley 2011; Rolf et al. 2012, 2014, 2018). By allowing the surface to evolve dynamically, rather than kinematically (e.g. McNamara and Zhong 2005; Heron and Lowman 2010, 2011, 2014; Davies et al. 2012; Heron et al. 2015; Li and Zhong 2017), the mantle interior is not ‘forced’ into particular mantle convection patterns. Modelling with dynamic boundary conditions could provide further information on the level of continental amalgamation needed to generate a ‘significant’ mantle response, and whether or not a ‘small’ supercontinent can generate a global geodynamics signature (as shown here).
Conclusions
We present the first 3D global mantle convection models aimed at identifying a thermal response to the formation of Pannotia. By analysing the mantle signature of the amalgamation of Pangaea and Rodinia (Figs 4 & 6), we find that an increase in core–mantle boundary heat flux and plume formation may be markers of a supercontinent (Fig. 5). Our results suggest that the convergence, assembly and dispersal of Pannotia had a thermal legacy that needs to be applied to models of Pangaea amalgamation. Such a legacy is not factored into models of tectonic evolution in which the transition from Rodinia to Pangaea represents one supercontinent cycle (Murphy et al. 2020). This thermal signature occurs in our numerical models during the time of Pannotia, indicating that such a continental amalgamation may have the ability to impact mantle dynamics (Fig. 5). This finding is supported by the development of plumes in the spatial and temporal location of the Central Iapetus Magmatic Province (CIMP). However, our findings are inconclusive on the potential role of large low shear velocity provinces (LLSVPs) in the development of plumes and subsequent break-up of supercontinents (Fig. 8). We did not find a significant difference in the plume locations and mantle dynamics of numerical simulations with and without thermo-chemical piles.
The results here represent an important contribution to a growing body of research that is posing questions on the fundamentals of global mantle convection: how do we define a supercontinent? In this initial study, we show that the formation of Pannotia generates an identifiable change in mantle flow, in line with processes related to the (less controversial) supercontinents of Pangaea and Rodinia. Here, we posit that identifying the mantle response to a supercontinent is key to formalizing a definition. This thermal signature occurs in our models during the time of continental convergence leading to Pannotia amalgamation, indicating that such a continental convergence may have had the ability to impact mantle dynamics (Fig. 5).
Acknowledgements
We thank reviewers Sergei Pisarevsky and Tobias Rolf for their insightful comments which led to significant improvements in the manuscript. This paper is part of UNESCO IGCP Project 648: Supercontinent Cycles and Global Geodynamics. Computations were performed on the Niagara supercomputer at the SciNet HPC Consortium (Loken et al. 2010).
Funding
PJH is grateful for funding from a Marie Skłodowska Curie Actions Individual Fellowship, proposal 749664. JBM acknowledges the continuing support of N.S.E.R.C. Canada (A06223) and RDN acknowledges support provided by Project CZ.02.2.69/0.0/0.0/16_015/0002362 of the Czech Operational Programme Research, Development and Education. SciNet is funded by the Canada Foundation for Innovation; the Government of Ontario; Ontario Research Fund-Research Excellence; and the University of Toronto. We thank the Computational Infrastructure for Geodynamics (geodynamics.org) which is funded by the National Science Foundation under award EAR-0949446 and EAR-1550901 for supporting the development of ASPECT.
Author contributions
PJH: data curation (lead), formal analysis (equal), funding acquisition (lead), investigation (lead), methodology (lead), project administration (lead), software (lead), validation (lead), visualization (lead), writing – original draft (equal); JBM: conceptualization (equal), formal analysis (equal), investigation (equal), writing – original draft (equal); RDN: conceptualization (equal), formal analysis (equal), investigation (equal), writing – original draft (equal); RNP: formal analysis (equal), funding acquisition (equal), resources (lead), writing – original draft (lead).
Data availability statement
Data and input files related to the numerical modelling of the paper will be made available at https://github.com/heronphi/Heron_Pannotia. A number of figures in this manuscript were generated using Generic Mapping Tools (Wessel et al. 2013) and applying Scientific Colour Maps (Crameri 2018).
- © 2020 The Author(s). Published by The Geological Society of London
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/).