Refine
Document Type
- Article (22)
- Postprint (3)
- Doctoral Thesis (1)
- Habilitation Thesis (1)
- Other (1)
- Review (1)
Language
- English (29)
Is part of the Bibliography
- yes (29)
Keywords
- subduction zone (4)
- geodynamics (3)
- models (3)
- numerical model (3)
- numerical modeling (3)
- plume (3)
- Pyrenees (2)
- architecture (2)
- continental lithosphere (2)
- extension (2)
Continental rift systems form by propagation of isolated rift segments that interact, and eventually evolve into continuous zones of deformation. This process impacts many aspects of rifting including rift morphology at breakup, and eventual ocean-ridge segmentation. Yet, rift segment growth and interaction remain enigmatic. Here we present geological data from the poorly documented Ririba rift (South Ethiopia) that reveals how two major sectors of the East African rift, the Kenyan and Ethiopian rifts, interact. We show that the Ririba rift formed from the southward propagation of the Ethiopian rift during the Pliocene but this propagation was short-lived and aborted close to the Pliocene-Pleistocene boundary. Seismicity data support the abandonment of laterally offset, overlapping tips of the Ethiopian and Kenyan rifts. Integration with new numerical models indicates that rift abandonment resulted from progressive focusing of the tectonic and magmatic activity into an oblique, throughgoing rift zone of near pure extension directly connecting the rift sectors.
Breakup Without Borders
(2019)
Relative plate motions during continental rifting result from the interplay of local with far-field forces. Here we study the dynamics of rifting and breakup using large-scale numerical simulations of mantle convection with self-consistent evolution of plate boundaries. We show that continental separation follows a characteristic evolution with four distinctive phases: (1) an initial slow rifting phase with low divergence velocities and maximum tensional stresses, (2) a synrift speed-up phase featuring an abrupt increase of extension rate with a simultaneous drop of tensional stress, (3) the breakup phase with inception of fast sea-floor spreading, and (4) a deceleration phase occurring in most but not all models where extensional velocities decrease. We find that the speed-up during rifting is compensated by subduction acceleration or subduction initiation even in distant localities. Our study illustrates new links between local rift dynamics, plate motions, and subduction kinematics during times of continental separation.
Inherited rheological structures in the lithosphere are expected to have large impact on the architecture of continental rifts. The Turkana depression in the East African Rift connects the Main Ethiopian Rift to the north with the Kenya rift in the south. This region is characterized by a NW-SE trending band of thinned crust inherited from a Mesozoic rifting event, which is cutting the present-day N-S rift trend at high angle. In striking contrast to the narrow rifts in Ethiopia and Kenya, extension in the Turkana region is accommodated in subparallel deformation domains that are laterally distributed over several hundred kilometers. We present both analog experiments and numerical models that reproduce the along-axis transition from narrow rifting in Ethiopia and Kenya to a distributed deformation within the Turkana depression. Similarly to natural observations, our models show that the Ethiopian and Kenyan rifts bend away from each other within the Turkana region, thus forming a right-lateral step over and avoiding a direct link to form a continuous N-S depression. The models reveal five potential types of rift linkage across the preexisting basin: three types where rifts bend away from the inherited structure connecting via a (1) wide or (2) narrow rift or by (3) forming a rotating microplate, (4) a type where rifts bend towards it, and (5) straight rift linkage. The fact that linkage type 1 is realized in the Turkana region provides new insights on the rheological configuration of the Mesozoic rift system at the onset of the recent rift episode. Plain Language Summary The Turkana depression in the Kenya/Ethiopia borderland is most famous for its several million years old human fossils, but it also holds a rich geological history of continental separation. The Turkana region is a lowland located between the East African and Ethiopian domes because its crust and mantle have been stretched in a continent-wide rift event during Cretaceous times about 140-120 Ma ago. This thin lithosphere exerted paramount control on the dynamics of East African rifting in this area, which commenced around 15 Ma ago and persists until today. Combining analog "sandbox" experiments with numerical geodynamic modeling, we find that inherited rift structures explain the dramatic change in rift style from deep, narrow rift basins north and south of the Turkana area to wide, distributed deformation within the Turkana depression. The failed Cretaceous rift is also responsible for the eastward bend of the Ethiopian rift and the westward bend of the Kenyan rift when entering the Turkana depression, which generated the characteristic hook-shaped form of present-day Lake Turkana. Combing two independent modeling techniques-analog and numerical experiments-is a very promising approach allowing to draw robust conclusions about the processes that shape the surface of our planet.
Complex, time-dependent, and asymmetric rift geometries are observed throughout the East African Rift System (EARS) and are well documented, for instance, in the Kenya Rift. To unravel asymmetric rifting processes in this region, we conduct 2D geodynamic models. We use the finite element software ASPECT employing visco-plastic rheologies, mesh-refinement, distributed random noise seeding, and a free surface. In contrast to many previous numerical modeling studies that aimed at understanding final rifted margin symmetry, we explicitly focus on initial rifting stages to assess geodynamic controls on strain localization and fault evolution. We thereby link to geological and geophysical observations from the Southern and Central Kenya Rift. Our models suggest a three-stage early rift evolution that dynamically bridges previously inferred fault-configuration phases of the eastern EARS branch: (1) accommodation of initial strain localization by a single border fault and flexure of the hanging-wall crust, (2) faulting in the hanging-wall and increasing upper-crustal faulting in the rift-basin center, and (3) loss of pronounced early stage asymmetry prior to basinward localization of deformation. This evolution may provide a template for understanding early extensional faulting in other branches of the East African Rift and in asymmetric rifts worldwide. By modifying the initial random noise distribution that approximates small-scale tectonic inheritance, we show that a spectrum of first-order fault configurations with variable symmetry can be produced in models with an otherwise identical setup. This approach sheds new light on along-strike rift variability controls in active asymmetric rifts and proximal rifted margins.
Linking deep seismic profiles with regional-scale gravity inversion is a powerful tool to deduce the architecture of rifted margins and their structural evolution. Here we map upper and lower crustal thicknesses of the northern South China Sea (SCS) margin in order to investigate the occurrence of depth-dependent crustal extension from the proximal to the distal margin. By comparing upper and lower crustal stretching factors, we find that the northern margin of the SCS is segmented in three parts: (1) sedimentary basins where upper crust is stretched more than lower crust, (2) distal margin where lower crust is stretched more than upper crust, (3) mostly proximal margin regions where the two layers have similar stretching factors. Our results suggest that sedimentary basins and distal margin prominently feature depth-dependent extension, however accommodated by different processes. While differential thinning within sedimentary basins appears to be governed by lateral pressure variations inducing lower crustal flow, we suggest the distal margin to be affected by a combination of mantle flow-induced lower crustal shearing and sequential fault activity during crustal hyper-extension.
Observations of rift and rifted margin architecture suggest that significant spatial and temporal structural heterogeneity develops during the multiphase evolution of continental rifting. Inheritance is often invoked to explain this heterogeneity, such as preexisting anisotropies in rock composition, rheology, and deformation. Here, we use high-resolution 3-D thermal-mechanical numerical models of continental extension to demonstrate that rift-parallel heterogeneity may develop solely through fault network evolution during the transition from distributed to localized deformation. In our models, the initial phase of distributed normal faulting is seeded through randomized initial strength perturbations in an otherwise laterally homogeneous lithosphere extending at a constant rate. Continued extension localizes deformation onto lithosphere-scale faults, which are laterally offset by tens of km and discontinuous along-strike. These results demonstrate that rift- and margin-parallel heterogeneity of large-scale fault patterns may in-part be a natural byproduct of fault network coalescence.
Observations of rift and rifted margin architecture suggest that significant spatial and temporal structural heterogeneity develops during the multiphase evolution of continental rifting. Inheritance is often invoked to explain this heterogeneity, such as preexisting anisotropies in rock composition, rheology, and deformation. Here, we use high-resolution 3-D thermal-mechanical numerical models of continental extension to demonstrate that rift-parallel heterogeneity may develop solely through fault network evolution during the transition from distributed to localized deformation. In our models, the initial phase of distributed normal faulting is seeded through randomized initial strength perturbations in an otherwise laterally homogeneous lithosphere extending at a constant rate. Continued extension localizes deformation onto lithosphere-scale faults, which are laterally offset by tens of km and discontinuous along-strike. These results demonstrate that rift- and margin-parallel heterogeneity of large-scale fault patterns may in-part be a natural byproduct of fault network coalescence.
Continental rifting is responsible for the generation of major sedimentary basins, both during rift inception and during the formation of rifted continental margins. Geophysical and field studies revealed that rifts feature complex networks of normal faults but the factors controlling fault network properties and their evolution are still matter of debate. Here, we employ high-resolution 2D geodynamic models (ASPECT) including two-way coupling to a surface processes (SP) code (FastScape) to conduct 12 models of major rift types that are exposed to various degrees of erosion and sedimentation. We further present a novel quantitative fault analysis toolbox (Fatbox), which allows us to isolate fault growth patterns, the number of faults, and their length and displacement throughout rift history. Our analysis reveals that rift fault networks may evolve through five major phases: (a) distributed deformation and coalescence, (b) fault system growth, (c) fault system decline and basinward localization, (d) rift migration, and (e) breakup. These phases can be correlated to distinct rifted margin domains. Models of asymmetric rifting suggest rift migration is facilitated through both ductile and brittle deformation within a weak exhumation channel that rotates subhorizontally and remains active at low angles. In sedimentation-starved settings, this channel satisfies the conditions for serpentinization. We find that SP are not only able to enhance strain localization and to increase fault longevity but that they also reduce the total length of the fault system, prolong rift phases and delay continental breakup.
Flexural strike-slip basins
(2021)
Strike-slip faults are classically associated with pull-apart basins where continental crust is thinned between two laterally offset fault segments. We propose a subsidence mechanism to explain the formation of a new type of basin where no substantial segment offset or synstrike-slip thinning is observed. Such "flexural strike-slip basins" form due to a sediment load creating accommodation space by bending the lithosphere. We use a two-way coupling between the geodynamic code ASPECT and surface-processes code FastScape to show that flexural strike-slip basins emerge if sediment is deposited on thin lithosphere close to a strike slip fault. These conditions were met at the Andaman Basin Central fault (Andaman Sea, Indian Ocean), where seismic reflection data provide evidence of a laterally extensive flexural basin with a depocenter located parallel to the strike-slip fault trace.
With the Late Cretaceous onset of Africa-Iberia-Europe convergence Central Europe experienced a pulse of intraplate shortening lasting some 15-20 Myr. This deformation event documents area-wide deviatoric compression of Europe and has been interpreted as a far-field response to Africa-Iberia-Europe convergence. However, the factors that governed the compression of Europe and conditioned the transient character of the deformation event have remained unclear. Based on mechanical considerations, numerical simulations, and geological reconstructions, we examine how the dynamics of intraplate deformation were governed by the formation of a convergent plate boundary fault between Iberia and Europe. During the Late Cretaceous, plate convergence was accommodated by the inversion of a young hyperextended rift system separating Iberia from Europe. Our analysis shows that the strength of the lithosphere beneath this rift was initially sufficient to transmit large compressive stresses far into Europe, though the lithosphere beneath the rift was thinned and thermally weakened. Continued convergence forced the formation of the plate boundary fault between Iberia and Europe. The fault evolved progressively and constituted a lithospheric-scale structure at the southern margin of Europe that weakened rheologically. This development caused a decrease in mechanical coupling between Iberia and Europe and a reduction of compressional far field stresses, which eventually terminated intraplate deformation in Central Europe. Taken together, our findings suggest that the Late Cretaceous intraplate deformation event records a high force transient that relates to the earliest strength evolution of a lithospheric-scale plate boundary fault.
We evaluate the spatial and temporal evolution of Earth's long-wavelength surface dynamic topography since the Jurassic using a series of high-resolution global mantle convection models. These models are Earth-like in terms of convective vigour, thermal structure, surface heat-flux and the geographic distribution of heterogeneity. The models generate a degree-2-dominated spectrum of dynamic topography with negative amplitudes above subducted slabs (i.e. circum-Pacific regions and southern Eurasia) and positive amplitudes elsewhere (i.e. Africa, north-western Eurasia and the central Pacific). Model predictions are compared with published observations and subsidence patterns from well data, both globally and for the Australian and southern African regions. We find that our models reproduce the long-wavelength component of these observations, although observed smaller-scale variations are not reproduced. We subsequently define "geodynamic rules" for how different surface tectonic settings are affected by mantle processes: (i) locations in the vicinity of a subduction zone show large negative dynamic topography amplitudes; (ii) regions far away from convergent margins feature long-term positive dynamic topography; and (iii) rapid variations in dynamic support occur along the margins of overriding plates (e.g. the western US) and at points located on a plate that rapidly approaches a subduction zone (e.g. India and the Arabia Peninsula). Our models provide a predictive quantitative framework linking mantle convection with plate tectonics and sedimentary basin evolution, thus improving our understanding of how subduction and mantle convection affect the spatio-temporal evolution of basin architecture.
Rock deformation at depths in the Earth’s crust is often localized in high temperature shear zones occurring at different scales in a variety of lithologies. The presence of material heterogeneities is known to trigger shear zone development, but the mechanisms controlling initiation and evolution of localization are not fully understood. To investigate the effect of loading conditions on shear zone nucleation along heterogeneities, we performed torsion experiments under constant twist rate (CTR) and constant torque (CT) conditions in a Paterson-type deformation apparatus. The sample assemblage consisted of cylindrical Carrara marble specimens containing a thin plate of Solnhofen limestone perpendicular to the cylinder’s longitudinal axis. Under experimental conditions (900 °C, 400 MPa confining pressure), samples were plastically deformed and limestone is about 9 times weaker than marble, acting as a weak inclusion in a strong matrix. CTR experiments were performed at maximum bulk shear strain rates of ~ 2*10-4s-1, yielding peak shear stresses of ~ 20 MPa. CT tests were conducted at shear stresses of ~ 20 MPa resulting in bulk shear strain rates of 1-4*10-4s-1. Experiments were terminated at maximum bulk shear strains of ~ 0.3 and 1.0.Strain was localized within the Carrara marble in front of the inclusion in an area of strongly deformed grains and intense grain size reduction. Locally, evidences for coexisting brittle deformation are also observed regardless of the imposed loading conditions. The local shear strain at the inclusion tipis up to 30 times higher than the strain in the adjacent host rock, rapidly dropping to 5times higher at larger distance from the inclusion. At both investigated bulk strains, the evolution of microstructural and textural parameters is independent of loading conditions. Ourresults suggest that loading conditions do not significantly affect material heterogeneity-induced strain localization during its nucleation and transient stages.
The lithosphere is often assumed to reside in a thermal steady-state when quantitatively describing the temperature distribution in continental interiors and sedimentary basins, but also at active plate boundaries. Here, we investigate the applicability limit of this assumption at slowly deforming continental rifts. To this aim, we assess the tectonic thermal imprint in numerical experiments that cover a range of realistic rift configurations. For each model scenario, the deviation from thermal equilibrium is evaluated. This is done by comparing the transient temperature field of every model to a corresponding steady-state model with an identical structural configuration. We find that the validity of the thermal steady-state assumption strongly depends on rift type, divergence velocity, sampling location, and depth within the rift. Maximum differences between transient and steady-state models occur in narrow rifts, at the rift sides, and if the extension rate exceeds 0.5-2 mm/a. Wide rifts, however, reside close to thermal steady-state even for high extension velocities. The transient imprint of rifting appears to be overall negligible for shallow isotherms with a temperature less than 100 degrees C. Contrarily, a steady-state treatment of deep crustal isotherms leads to an underestimation of crustal temperatures, especially for narrow rift settings. Thus, not only relatively fast rifts like the Gulf of Corinth, Red Sea, and Main Ethiopian Rift, but even slow rifts like the Kenya Rift, Rhine Graben, and Rio Grande Rift must be expected to feature a pronounced transient component in the temperature field and to therefore violate the thermal steady-state assumption for deeper crustal isotherms.
Submarine landslides can generate local tsunamis posing a hazard to human lives and coastal facilities. Two major related problems are: (i) quantitative estimation of tsunami hazard and (ii) early detection of the most dangerous landslides. This thesis focuses on both those issues by providing numerical modeling of landslide-induced tsunamis and by suggesting and justifying a new method for fast detection of tsunamigenic landslides by means of tiltmeters. Due to the proximity to the Sunda subduction zone, Indonesian coasts are prone to earthquake, but also landslide tsunamis. The aim of the GITEWS-project (German-Indonesian Tsunami Early Warning System) is to provide fast and reliable tsunami warnings, but also to deepen the knowledge about tsunami hazards. New bathymetric data at the Sunda Arc provide the opportunity to evaluate the hazard potential of landslide tsunamis for the adjacent Indonesian islands. I present nine large mass movements in proximity to Sumatra, Java, Sumbawa and Sumba, whereof the largest event displaced 20 km³ of sediments. Using numerical modeling, I compute the generated tsunami of each event, its propagation and runup at the coast. Moreover, I investigate the age of the largest slope failures by relating them to the Great 1977 Sumba earthquake. Continental slopes off northwest Europe are well known for their history of huge underwater landslides. The current geological situation west of Spitsbergen is comparable to the continental margin off Norway after the last glaciation, when the large tsunamigenic Storegga slide took place. The influence of Arctic warming on the stability of the Svalbard glacial margin is discussed. Based on new geophysical data, I present four possible landslide scenarios and compute the generated tsunamis. Waves of 6 m height would be capable of reaching northwest Europe threatening coastal areas. I present a novel technique to detect large submarine landslides using an array of tiltmeters, as a possible tool in future tsunami early warning systems. The dislocation of a large amount of sediment during a landslide produces a permanent elastic response of the earth. I analyze this response with a mathematical model and calculate the theoretical tilt signal. Applications to the hypothetical Spitsbergen event and the historical Storegga slide show tilt signals exceeding 1000 nrad. The amplitude of landslide tsunamis is controlled by the product of slide volume and maximal velocity (slide tsunamigenic potential). I introduce an inversion routine that provides slide location and tsunamigenic potential, based on tiltmeter measurements. The accuracy of the inversion and of the estimated tsunami height near the coast depends on the noise level of tiltmeter measurements, the distance of tiltmeters from the slide, and the slide tsunamigenic potential. Finally, I estimate the applicability scope of this method by employing it to known landslide events worldwide.
Continental rift systems open up unique possibilities to study the geodynamic system of our planet: geodynamic localization processes are imprinted in the morphology of the rift by governing the time-dependent activity of faults, the topographic evolution of the rift or by controlling whether a rift is symmetric or asymmetric. Since lithospheric necking localizes strain towards the rift centre, deformation structures of previous rift phases are often well preserved and passive margins, the end product of continental rifting, retain key information about the tectonic history from rift inception to continental rupture.
Current understanding of continental rift evolution is based on combining observations from active rifts with data collected at rifted margins. Connecting these isolated data sets is often accomplished in a conceptual way and leaves room for subjective interpretation. Geodynamic forward models, however, have the potential to link individual data sets in a quantitative manner, using additional constraints from rock mechanics and rheology, which allows to transcend previous conceptual models of rift evolution. By quantifying geodynamic processes within continental rifts, numerical modelling allows key insight to tectonic processes that operate also in other plate boundary settings, such as mid ocean ridges, collisional mountain chains or subduction zones.
In this thesis, I combine numerical, plate-tectonic, analytical, and analogue modelling approaches, whereas numerical thermomechanical modelling constitutes the primary tool. This method advanced rapidly during the last two decades owing to dedicated software development and the availability of massively parallel computer facilities. Nevertheless, only recently the geodynamical modelling community was able to capture 3D lithospheric-scale rift dynamics from onset of extension to final continental rupture.
The first chapter of this thesis provides a broad introduction to continental rifting, a summary of the applied rift modelling methods and a short overview of previews studies. The following chapters, which constitute the main part of this thesis feature studies on plate boundary dynamics in two and three dimension followed by global scale analyses (Fig. 1).
Chapter II focuses on 2D geodynamic modelling of rifted margin formation. It highlights the formation of wide areas of hyperextended crustal slivers via rift migration as a key process that affected many rifted margins worldwide. This chapter also contains a study of rift velocity evolution, showing that rift strength loss and extension velocity are linked through a dynamic feed-back. This process results in abrupt accelerations of the involved plates during rifting illustrating for the first time that rift dynamics plays a role in changing global-scale plate motions. Since rift velocity affects key processes like faulting, melting and lower crustal flow, this study also implies that the slow-fast velocity evolution should be imprinted in rifted margin structures.
Chapter III relies on 3D Cartesian rift models in order to investigate various aspects of rift obliquity. Oblique rifting occurs if the extension direction is not orthogonal to the rift trend. Using 3D lithospheric-scale models from rift initialisation to breakup I could isolate a characteristic evolution of dominant fault orientations. Further work in Chapter III addresses the impact of rift obliquity on the strength of the rift system. We illustrate that oblique rifting is mechanically preferred over orthogonal rifting, because the brittle yielding requires a lower tectonic force. This mechanism elucidates rift competition during South Atlantic rifting, where the more oblique Equatorial Atlantic Rift proceeded to breakup while the simultaneously active but less oblique West African rift system became a failed rift. Finally this Chapter also investigates the impact of a previous rift phase on current tectonic activity in the linkage area of the Kenyan with Ethiopian rift. We show that the along strike changes in rift style are not caused by changes in crustal rheology. Instead the rift linkage pattern in this area can be explained when accounting for the thinned crust and lithosphere of a Mesozoic rift event.
Chapter IV investigates rifting from the global perspective. A first study extends the oblique rift topic of the previous chapter to global scale by investigating the frequency of oblique rifting during the last 230 million years. We find that approximately 70% of all ocean-forming rift segments involved an oblique component of extension where obliquities exceed 20°. This highlights the relevance of 3D approaches in modelling, surveying, and interpretation of many rifted margins. In a final study, we propose a link between continental rift activity, diffuse CO2 degassing and Mesozoic/Cenozoic climate changes. We used recent CO2 flux measurements in continental rifts to estimate worldwide rift-related CO2 release, which we based on the global extent of rifts through time. The first-order correlation to paleo-atmospheric CO2 proxy data suggests that rifts constitute a major element of the global carbon cycle.
Over the past few decades, azimuthal seismic anisotropy measurements have been widely used proxy to study past and present-day deformation of the lithosphere and to characterize convection in the mantle. Beneath continental regions, distinguishing between shallow and deep sources of anisotropy remains difficult due to poor depth constraints of measurements and a lack of regional-scale geodynamic modeling. Here, we constrain the sources of seismic anisotropy beneath Madagascar where a complex pattern cannot be explained by a single process such as absolute plate motion, global mantle flow, or geology. We test the hypotheses that either Edge-Driven Convection (EDC) or mantle flow derived from mantle wind interactions with lithospheric topography is the dominant source of anisotropy beneath Madagascar. We, therefore, simulate two sets of mantle convection models using regional-scale 3-D computational modeling. We then calculate Lattice Preferred Orientation that develops along pathlines of the mantle flow models and use them to calculate synthetic splitting parameters. Comparison of predicted with observed seismic anisotropy shows a good fit in northern and southern Madagascar for the EDC model, but the mantle wind case only fits well in northern Madagascar. This result suggests the dominant control of the measured anisotropy may be from EDC, but the role of localized fossil anisotropy in narrow shear zones cannot be ruled out in southern Madagascar. Our results suggest that the asthenosphere beneath northern and southern Madagascar is dominated by dislocation creep. Dislocation creep rheology may be dominant in the upper asthenosphere beneath other regions of continental lithosphere.
Movements of tectonic plates often induce oblique deformation at divergent plate boundaries. This is in striking contrast with traditional conceptual models of rifting and rifted margin formation, which often assume 2-D deformation where the rift velocity is oriented perpendicular to the plate boundary. Here we quantify the validity of this assumption by analysing the kinematics of major continent-scale rift systems in a global plate tectonic reconstruction from the onset of Pangea breakup until the present day. We evaluate rift obliquity by joint examination of relative extension velocity and local rift trend using the script-based plate reconstruction software pyGPlates. Our results show that the global mean rift obliquity since 230 Ma amounts to 34 degrees with a standard deviation of 24 degrees, using the convention that the angle of obliquity is spanned by extension direction and rift trend normal. We find that more than similar to 70 % of all rift segments exceeded an obliquity of 20 degrees demonstrating that oblique rifting should be considered the rule, not the exception. In many cases, rift obliquity and extension velocity increase during rift evolution (e.g. Australia-Antarctica, Gulf of California, South Atlantic, India-Antarctica), which suggests an underlying geodynamic correlation via obliquity-dependent rift strength. Oblique rifting produces 3-D stress and strain fields that cannot be accounted for in simplified 2-D plane strain analysis. We therefore highlight the importance of 3-D approaches in modelling, surveying, and interpretation of most rift segments on Earth where oblique rifting is the dominant mode of deformation.
Initiation of subduction following the impingement of a hot buoyant mantle plume is one of the few scenarios that allow breaking the lithosphere and recycling a stagnant lid without requiring any preexisting weak zones. Here, we investigate factors controlling the number and shape of retreating subducting slabs formed by plume-lithosphere interaction. Using 3-D thermomechanical models we show that the deformation regime, which defines formation of single-slab or multi-slab subduction, depends on several parameters such as age of oceanic lithosphere, thickness of the crust and large-scale lithospheric extension rate. Our model results indicate that on present-day Earth multi-slab plume-induced subduction is initiated only if the oceanic lithosphere is relatively young (<30-40 Myr, but >10 Myr), and the crust has a typical thickness of 8 km. In turn, development of single-slab subduction is facilitated by older lithosphere and pre-imposed extensional stresses. In early Earth, plume-lithosphere interaction could have led to formation of either episodic short-lived circular subduction when the oceanic lithosphere was young or to multi-slab subduction when the lithosphere was old.