We design a series of simulations in a 3-D Cartesian domain. We use a multi-material approach to simulate subduction of a composite visco-plastic plate into a viscous mantle, with no upper plate (i.e., “free subduction”). While numerous studies demonstrate that the upper plate can affect subduction dynamics significantly (e.g. van Hunen et al., 2004; Espurt et al., 2008), its absence allows us to analyse the dynamic forcing by the subducting plate which is widely agreed to be the main driver of subduction dynamics and diversity (e.g. Forsyth and Uyeda, 1975; Goes et al., 2017). Our approach neglects the thermal evolution of the slab and associated feedbacks on density and viscosity (e.g., Garel et al., 2014; Suchoy et al., 2021). This simplification enhances numerical efficiency, allowing us to perform a systematic 3-D study across a wide parameter space. Previous studies have demonstrated that the multi-material, mechanical approach used herein captures the first-order dynamics of subduction (e.g. Bellahsen et al., 2005; Capitanio et al., 2007; Stegman et al., 2010; Chen et al., 2022).

We use Fluidity, an adaptive, unstructured mesh, finite-element, control-volume computational modelling framework (e.g., Wilson, 2009; Davies et al., 2011; Kramer et al., 2012, 2021), to solve the conservation equations of mass and momentum for an incompressible fluid under the infinite Prandtl number and Boussinesq approximations: ∇ ⋅ u = 0 (1) ∇ ⋅ η ∇ ⃗ u + ∇ ⃗ u T − ∇ ⃗ p = − g Δ ρ Γ y ̂ (2) ∂ Γ ∂ t + u ⋅ ∇ ⃗ Γ = 0 (3) where u is the velocity, η is the dynamic viscosity, p is the pressure, g is the acceleration due to gravity, Δρ is the difference in density between different materials, Γ is the material volume fraction (Γ = 1 within a given material and Γ = 0 elsewhere) and y ̂ is a unit vector in the direction of gravity. We utilise an adaptive unstructured mesh of tetrahedral elements, with minimum and maximum element sizes of 3 and 300 km, respectively. This allows us to resolve fine-scale features where the gradients of velocity and viscosity are strong, whilst maintaining computational efficiency (see Davies et al., 2007, 2011, for further detail).

Our reference young and old “no-ridge” models (Table 1) follow two Cartesian model designs from Chen et al. (2022) (models W2400_young and W2400_ref, respectively). We use a domain of 4,000 × 2890 × 4,000 km³ (length×depth×width) to simulate half the subducting plate, assuming symmetry at the centre of the plate along the long axis (z = 0) (Figure 1).

We apply a free-surface boundary condition at the top of the domain and free-slip boundary conditions at all other boundaries. The initial half-plate is 2,200 km long (in the direction of subduction), 1,200 km wide (along strike) and 45 or 70 km thick, for young or old plates, respectively. We limit lateral flow from the mantle to the section below the plate by adding a side plate with an initial gap of 22 km away from the edge of the subducting plate (e.g., Holt and Becker, 2017). A 600 km gap was prescribed between the trailing edge of the subducting plate and the domain boundary at x = 0. The initial slab morphology follows an arc with a radius of 250 km along the top surface, to a dip of 77° (as in Garel et al., 2014; Suchoy et al., 2021), extending to an initial depth of 200 km.

The subducting plate is composed of 3 layers with a strong isoviscous core in the centre, and visco-plastic layers at the top and bottom. The core has a viscosity 100 times η _UM . The visco-plastic layers have this same initial viscosity, but follow a von Mises yield criterion (e.g., OzBench et al., 2008) so that: η P = τ I I 2 ⋅ ε ̇ I I , if τ I I < τ y τ y 2 ⋅ ε ̇ I I , if τ I I ≥ τ y (4) where the second invariant of the stress tensor τ I I = 2 η ⋅ ε ̇ I I , with ε ̇ I I the second invariant of the strain-rate tensor and τ _y the yield stress. This rheological layering captures the concentration of slab strength in the centre of the dense plate, as expected from temperature and strain-rate dependent mantle rheology (OzBench et al., 2008; Capitanio et al., 2009; Buffett and Becker, 2012). For the old reference plate, the core is 30 km thick and the top and bottom layers are each 20 km thick, while for the young plate each layer is 15 km thick. The density of the older plate is set to be higher than that of the young plate, such that the combined densities and thicknesses of the two background plates yield a range of plate buoyancies similar to that expected for 20 and 120 Myr lithosphere from plate cooling models (e.g., McKenzie et al., 2005). We applied the same density contrast to the full thickness of the plate, i.e., only considering joint buoyancy effects of the crust and mantle lithosphere. This results in the old plate being more negatively buoyant and effectively stronger than the young plate. The upper-lower mantle boundary (ULMB) is implemented as a viscosity jump at 660 km depth, with an upper mantle viscosity, η _UM , of 2 ⋅ 10²⁰ Pa ⋅ s and a lower mantle viscosity, η _LM , of 50 times η _UM . The side plate has a viscosity of 1,000 times η _UM . Material properties common to all simulations are summarised in Table 2.

To examine the effect of buoyant ridges (see Table 1), we decrease the density in a trench-perpendicular strip of the subducting plate (Figure 1B). We ran cases with a ridge positioned either at the centre of the plate (i.e., the symmetry edge of the half-plate, at z = 0) or positioned offset from the centre, in the centre of the half plate (initially centred at z = 600 km). The off-centre ridge cases effectively include two ridges (mirrored across the symmetry plane). The distance between these ridges is sufficiently large to allow us to investigate the effect of a single offset ridge. For off-centre ridges, we use a 200 km wide segment. For the central ridges, we use a 100 km wide segment of the plate, to represent half of a 200 km wide symmetrical ridge.

We consider 8 cases (Table 1). For each reference plate type (old and young), we ran a case without any ridges (“NR” or “no-ridge” models, Table 1), as well cases with a central ridge of high and low excess positive buoyancy (“HB” and “LB” models) relative to the buoyancy of the plate (i.e., counteracting the negative buoyancy of the downgoing plate). We also ran simulations for each plate type with a low positive buoyancy (“LB” cases) ridge at the side of the plate. Based on the topography of buoyant ridges on Earth today, we chose cases that yield a typical excess relief of 0.75 and 1.5 km (relative to background plate bathymetry) to define the parameters for our low and high buoyancy ridges, respectively. We calculated what difference in buoyancy for each ridge type relative to surrounding lithosphere would yield this excess relief assuming isostatic equilibrium (following Cloos, 1993; Gutscher et al., 2000; McKenzie et al., 2005, see supplementary material for further details). We then determined the corresponding density difference over the thickness of our modelled plates and reduced it accordingly. Since the modelled young plate is thinner than the old plate, this resulted in much greater density contrast between the ridge and the young plate, compared with the old plate.

When calculating subduction diagnostics, the plate-mantle interface was delineated as the iso-surface where the mantle material volume fraction is 0.5. Snapshots of slab morphologies for all 8 cases, following interaction with the ULMB, are displayed in Figure 2. The trench is defined as the front edge of the plate at 20 km depth (Figure 3). We measure the motion and velocity of the trench and the tip of the slab at 3 points along strike: 15 km from the symmetry boundary, 50 km from the edge of the plate and at a point half way along the trench (denoted as “symmetry,” “edge” and “middle,” respectively). The location of the edge and the middle was calculated throughout the simulations to account for temporal changes in trench shape (red diamonds in Figure 3).

We first discuss the behaviour of the two reference “no-ridge” models (NR_Old and NR_Young). The evolution of both is consistent with similar 2-D and 3-D models in previous studies (Funiciello et al., 2003; Schellart, 2005; Capitanio et al., 2007; Stegman et al., 2010) and includes two phases: initial sinking of the slab through the upper mantle and subsequent interaction of the slab with the lower mantle. We refer to the time of the first slab-ULMB interaction as t660.

In case NR_Old, the slab sinks through the upper mantle at an increasing rate, reaching the ULMB after ∼7 Myr (Figure 4A). Over the same period, the trench retreats at rates between 1 and 3 cm/yr, depending on the position along the trench. The highest cumulative trench retreat occurs at the middle of the half-plate (∼120 km in 7 Myr), with a lower value at the symmetry plane (∼100 km in 7 Myr) and the lowest value at the slab edge (∼70 km in 7 Myr; Figures 3A, 5A). As the slab approaches the ULMB, it flattens due to interaction with the more viscous lower mantle, which reduces the sinking velocity from an average rate of ∼6.5 cm/yr before t660 to ∼2.7 cm/yr after t660 (Figure 4A). Upon interaction with the ULMB, the trench stagnates and even advances slightly at the edge and at the symmetry boundary (Figures 3A, 5A). As a result, the trench develops a “W”-shape with a stagnation point, where trench motion tends to zero, at the centre (i.e., at the symmetry boundary) surrounded by faster retreating segments, with a maximum in retreat at the middle point of the half-plate (and its symmetrical equivalent).

In case NR_Young, the slab has less negative buoyancy and, accordingly, sinks at a slower rate of ∼2 cm/yr, reaching the ULMB after 24 Myr. Owing to its reduced thickness, the NR_Young plate is also weaker and, consequently, as the slab approaches the ULMB, it buckles (Figure 2E), as its sinking velocity reduces to ∼0.5 cm/yr (Figure 4E). Both before and after t660, the trench slowly advances at a steady rate of 0.1–0.5 cm/yr (Figure 5E). The rate of trench advance is similar everywhere along the trench, resulting in a reasonably straight “I”-shaped trench (Figure 3E).

In both the young and old reference models, the shape of the trench is determined by slab-mantle interactions (Schellart et al., 2007). The difference in trench retreat along strike for the old plate is due to return flow caused by slab rollback, which creates toroidal flow cells in the mantle near the edges of the plate. When the toroidal cell is smaller than the width of the half-plate, a stagnation point emerges at the centre of the trench (Schellart et al., 2007; Schellart, 2020). Our reference old plate is wide enough for such a central stagnation point to develop (Chen et al., 2022), thereby providing an excellent basis for evaluating the effect of ridge position. Narrower old plates develop a “C”-shaped trench instead, where the centre of the slab retreats more than its sides. The young plate sinks slower, which provides more time for deformation and bending at the trench. As a result, the young plate sinks almost vertically and drives very little trench motion (e.g., Capitanio et al., 2007). The slight curvature at the edge of the plate is the result of the weak toroidal cell caused by the displacement of mantle material due to plate sinking, and overall the trench stays reasonably straight (Chen et al., 2022).

Next, we examine cases with ridges at the centre of the plate (i.e., at the symmetry plane), with low and high positive buoyancy relative to the downgoing plate: cases HB_C and LB_C. In all cases, the addition of a buoyant ridge results in less trench retreat, or even trench advance, at the location where the ridge impacts the trench. Compared with the NR_Old reference case, trench retreat at the centre of the plate in case LB_C_Old is reduced from ∼100 km to ∼30 km by t660 (Figure 5B), while in case HB_C_Old, the centre of the trench advances ∼15 km by t660 (Figure 5A). In the NR_Young reference case, the centre of the trench is already advancing, and the addition of the ridge enhances this advance. For case LB_C_Young, the trench at the centre of the plate advances by ∼70 km by the reference t660 compared with ∼30 km in NR_Young case (Figure 5E). For the HB_C_Young case, the trench advances ∼110 km by the reference t660 (Figure 5D). In both young and old plate models, trench retreat at locations away from the ridge is increased or trench advance is decreased compared to the corresponding NR case.

The sinking velocity of the slab carrying the ridge is reduced locally due to the added positive buoyancy of the ridge. Particularly in the young plate models, the local plate sinking velocity at the ridge is significantly lower than in the rest of the plate. In model HB_C_Young, the ridge has not sunk at all by the reference t660 (Figure 4D). The differential sinking rate in this model is accommodated through thinning of the slab at the side of the ridge, which could facilitate slab tearing under certain rheological parameterisations. In comparison, the less buoyant ridge in model LB_C_Young has sunk ∼200 km, from its initial depth, by the reference t660 (Figure 4E). For older plate models, sinking velocity at the ridge deviates less from that in the rest of the plate. The part of the slab with the ridge has sunk ∼300 km for the HB_C_Old case and ∼400 km for LB_C_Old case (Figures 4A,B). Elsewhere along the slab, the sinking velocity remains similar to the reference case for both young and old plate models.

Previous 2-D models have demonstrated that, in free-subduction, higher slab pull and increased resistance to bending at the trench encourage trench retreat (Capitanio et al., 2007; Ribe, 2010). In our 3-D models, it is the difference between slab pull at the ridge and in the rest of the plate that drives the variable trench retreat, hindered at the ridge and enhanced elsewhere. This differential trench motion increases trench curvature, leading to the formation of a “W”-shaped trench in the young plate models and a more pronounced “W”-shape in the old plate, with more retreat at the edge of the plate than at the symmetry boundary (Figures 3B,C,F,G). Due to the low resistance to bending in young plates, such a “W”-shaped trench is difficult to generate by young plate subduction without along-strike buoyancy variations (Chen et al., 2022). The differential sinking between the ridge and rest of the plate also creates along-strike tension in the plate. The higher plate strength of the old plates facilitates stress transmission along strike. As a result, old plates can pull the ridge more effectively into the mantle, resulting in less along-strike variations in sinking and a smoother trench and slab shape in comparison to the corresponding young plates with ridges (Figures 2B,C,F,G and Figures 3B,C,F,G).

We evaluate the impact of ridge location by comparing results from models LB_C and LB_S. As for the models with a central ridge, subduction of a ridge on the side of the plate results in locally reduced trench retreat and, sometimes, induces a small amount of trench advance (Figures 3D,H). Although the side ridge cases effectively contain two ridges, there is no indication from the resulting trench shapes (Figures 3D,H) that there is an interference due to superposition of the influence of the two ridges, so these models can be used to study the effect of the subduction of a ridge offset from the slab centre.

The along-strike changes in trench motion compared to the reference cases are less for the cases with a side ridge than for the cases with a central ridge. For case LB_S_Old, the reduction in trench retreat by t660 is ∼50 km, compared to a reduction of ∼80 km for case LB_C_Old (Figures 5B,C). This is likely because the impact of the ridge in LB_S_Old is concentrated at the point where a ridge-free trench (case NR_Old) would retreat the most. Therefore, in this case, ridge subduction reduces the overall lateral deflection of the trench in comparison to the reference case (Figures 3A,D). For case LB_S_Young, the increase in trench advance relative to case NR_Young by the reference t660 is ∼20 km, compared with ∼35 km for case LB_C_Young (Figures 5E,F). All young plate cases with a ridge deform the trench from the straight reference shape of the NR_Young case to a curved trench with a cusp (Figures 3E,G,H).

As for cases with central ridges, the sinking velocity is reduced relative to the corresponding NR cases where the side ridge is subducted. The local reduction in sinking velocity is smaller for the LB_S cases than for the LB_C cases. The ridge in case LB_S_Old has sunk ∼400 km by the reference t660, compared with ∼390 km for case LB_C_Old (Figures 4B,C). Similarly, in LB_S_Young, the ridge has subducted ∼10 km deeper by the reference t660 than case LB_C_Young at the same time (Figures 4E,F). The side-ridge cases, that effectively contain two ridges, also demonstrate that a subducting buoyant ridge can lower the overall sinking velocity of the plate if the ridge adds enough positive buoyancy. The reduction in overall sinking velocity is small, resulting in a delay of ∼0.2 Myr in t660 in model LB_S_Old. In LB_C_Old, where the change in buoyancy is only half of that in the LB_S_Old case, the t660 delay is half as small, ∼0.1 Myr. The effect is somewhat larger for the LB_S_Young case, where t660 is delayed by ∼5 Myr compared with the reference NR_Young case (Figure 4F).

As in the cases with central ridges, a side ridge’s main effect is to locally reduce trench retreat and sinking velocity. These effects are both diminished when the ridge subducts along a part of the trench that is, otherwise strongly retreating. When the ridge impacts at the side, the composite trends of the plate’s tendency to retreat and the reduced trench retreat of the ridge leads to less along-strike variation in trench shape. The lower amount of along-strike bending in the LB_S cases may explain the slightly higher sinking velocities for the side ridges compared to the central ridges, as less potential energy is lost in plate deformation. Thus, our results imply that ridges that subduct along parts of the trench where trench motion is already hampered by other factors impact subduction more.

While 2-D models can elucidate some of the effects of buoyant feature subduction, they overestimate the buoyancy force exerted by the feature, due to the implicit assumption of infinite feature width. The full effects of the subduction of buoyant features of limited along-strike extent can therefore only be investigated using 3-D models.

All our ridge-subduction cases lead to trench deformation, with decreased trench retreat, trench stagnation or increased trench advance where the buoyant ridges subduct. Previous 3-D models that investigated the effect of buoyant feature subduction (Martinod et al., 2005; Mason et al., 2010; Flórez-Rodríguez et al., 2019) also found that trench retreat was impeded where the buoyant feature subducts, leading to differential retreat along strike. This observation was independent of size, shape and orientation and applied to both buoyant ridges and plateaus (Martinod et al., 2005, 2013). The resulting formation of trench cusps can be offset from the position of ridge impact at the trench if the feature subducts obliquely (Martinod et al., 2013). We show that significant trench deformation occurs in both old and young plates. Plate strength and the relative buoyancy of the ridges compared to the underlying plate determine how pronounced trench deformation is. That is, trenches where a weaker young plate subducts are deformed more locally than those where a stronger old plate subducts, for the same excess ridge buoyancy.

We also find that the position of buoyant ridges along the trench affects the response of the trench and resulting slab morphology. In older plates, ridges near intrinsic trench stagnation points accentuate along-strike variations in trench shape while ridges impacting the trench at other locations lead to less pronounced curvature. Schellart (2020), in 3-D models without buoyant features, found that the central trench stagnation point can, at a late stage of subduction zone evolution, facilitate local slab shallowing. Our results suggest that the subduction of a buoyant ridge at or near such a stagnation point would further enhance such behaviour. We note that trench stagnation points can also form in response to other factors, for example, at the edge of a subducting plate, due to interaction with a side plate or existence of a triple junction, and by interaction with a thick upper plate (Capitanio et al., 2011; Jadamec and Billen, 2012).

Although it is commonly assumed that the subduction of buoyant ridges leads to a decreased slab dip (e.g., Gutscher et al., 2000), our models display a range of behaviour (Figure 6). In old plate cases (Figures 6A–D), the slab is steepened where the ridge subducts, while the dip of the slab along the rest of the plate is slightly lower than in the reference case without a ridge. Although buoyant ridges act to inhibit sinking, the older subducting plates are sufficiently strong that the pull from the rest of the plate leads to a steeper slab where the buoyant feature enters the subduction zone. In this case, the lowest dips occur where the trench retreats most. In contrast, ridges on young plates develop shallower dip angles than the surrounding plate, and the shallowing is stronger for ridges with a higher positive buoyancy (Figures 6E–H). In these cases, the rest of the plate is weak, and is not negatively buoyant enough to pull down the segment with the ridge. Martinod et al. (2005) found that their slabs shallowed where the buoyant feature was subducting but steepened below. Flórez-Rodríguez et al. (2019) found yet another variation in their 3-D analogue models of subduction of buoyant ridges on relatively narrow plates, where the slab shallows in the centre of the ridge but steepens on its sides. In a spherical geometry, effective resistance to plate bending at the trench is higher than in equivalent Cartesian plates (Morra et al., 2006; Mahadevan et al., 2010; Chamolly and Ribe, 2021; Chen et al., 2022). Models of spherical plates display trench and slab deformation in response to along-strike variations in slab buoyancy comparable to Cartesian models (Morra et al., 2006), and are therefore likely to respond to subduction of buoyant ridge similarly. The differences in model behaviour illustrate that the effect of buoyant features on slab dip varies depending on plate strength and relative feature buoyancy, and a range of responses are possible.

Our results confirm that ridge subduction alone is insufficient to generate a flat slab, a finding that is, consistent with other 3-D free-subduction models (Martinod et al., 2005; Flórez-Rodríguez et al., 2019). Even in 2-D models (van Hunen et al., 2004; Gerya et al., 2009), where the buoyant features are essentially infinite in their along-strike dimension, the subduction of buoyant features alone does not generate a fully flattened slab. Instead, a number of studies have now confirmed that slab flattening requires additional mechanisms such as an upper plate that (by external forcing) overrides the lower plate faster than it can subduct (Van Hunen et al., 2002; Arcay et al., 2008; Currie and Beaumont, 2011; Martinod et al., 2013; Liu and Currie, 2016) or a thick upper plate that generates an upward suction force in the relatively high viscosity mantle wedge, which can pull a shallow-angle slab into a flat-slab geometry (Gerya et al., 2009; O’Driscoll et al., 2009; Manea et al., 2012; Rodríguez-González et al., 2012, 2014; Taramón et al., 2015; Schellart and Strak, 2021).

While upper plate forcing or thickness can result in a flat slab in the absence of buoyant feature subduction, these mechanisms have clear limitations (Hu et al., 2016). Overthrusting requires a highly forced upper plate velocity (3–5 cm/yr, e.g. van Hunen et al., 2004; Currie and Beaumont, 2011) and would usually apply to the entire upper plate, failing to localise flat subduction in bounded regions of several hundred kilometres, as is observed on Earth at present. A similar issue arises for a thickened upper plate which needs to be of a craton-like thickness (O’Driscoll et al., 2009; Manea et al., 2012; Taramón et al., 2015). While cratons comprise limited parts of the upper plate (e.g. Wyoming craton in North America; Amazonian craton in South America), they extend over scales of 1,000 km (e.g. Kusky et al., 2014, and references therein) and, therefore, can only explain very wide flat slabs. Conversely, modelling studies demonstrate that the subduction of buoyant features can aid flat slab formation, by locally reducing trench motion and slab dip (van Hunen et al., 2004; Gerya et al., 2009; Martinod et al., 2013). Thus the contribution of buoyant feature subduction can help to explain the limited lateral extent and temporal duration of flat slabs.

Together with published dynamic models, our results illustrate a number of different responses to buoyant feature subduction, both ridges and plateaus, are possible. These include three types: 1) a “trench deforming” response is an end-member in which subduction of a buoyant feature mainly deforms the trench with minor slab steepening or contortion. Our results demonstrate that this occurs when the positive buoyancy force of a buoyant ridge is low compared with overall plate strength. The upper plate in a “trench deforming” subduction style is typically thin or weak and therefore does not impede trench deformation; 2) a “slab deforming” response, in which subduction of the buoyant feature reduces slab dip locally whilst also deforming the trench. This happens when the buoyancy force associated with the feature is strong enough to locally decrease slab dip while any forcing by the upper plate does not play a significant role; 3) a “slab flattening” response, an end-member in which the slab is flattened at shallow depth over a distance of up to several hundred kilometres. This response occurs when the positive buoyancy of a feature is substantial, and it subducts under a thick, possibly advancing, upper plate. Observations of buoyant features subducting at the plate boundaries of the Pacific and Indian Oceans reflect these different expressions (see Figure 7 and Table 3), although many subduction zones do have additional complexities.

The best example of “trench-deforming” buoyant-feature subduction occurs along the Izu-Bonin-Mariana (IBM) system. The two main buoyant features interacting with the IBM subduction system and impeding retreat of the Mariana trench are both located near major cusps (i.e. the Caroline ridge near the Mariana-Yap junction and the Ogasawara ridge near Bonin-Mariana junction) and these cusps are inferred to have developed when these buoyant features began interacting with the trench (Miller et al., 2006a,c; Seton et al., 2012; Xia et al., 2020). The slab below the Ogasawara Plateau is slightly steeper than adjacent slab segments subducting north and south of it (Miller et al., 2006c; Hayes et al., 2018). The Caroline Island Ridge at the southern end of Mariana is located at a clear cusp where it has reduced convergence velocity significantly and blocked subduction (Miller et al., 2006a; Xia et al., 2021). The slab on the side of the Caroline Ridge, subducting below the southern Mariana trench (Miller et al., 2006a; Zhu et al., 2019) is steep, while the part of the slab along the Yap trench may have detached in response to impingement of the Caroline Ridge (Zhang and Zhang, 2020). Another example of “trench-deforming” buoyant-feature subduction is at the southern end of the Tonga trench, where the Hikurangi Plateau appears to have served as a pivot for the rotation and migration of the Tonga-Kermadec system (Schellart and Spakman, 2012; Seton et al., 2012). Although there is some shallowing of the Hikurangi-southern Kermedec slab above 100 km, this shallowing occurs mainly north of the plateau, while below 100 km depth, the slab is steeper than further north (Reyners et al., 2011; Hayes et al., 2018). These examples all display the expected behaviour for old plates with sufficient slab pull to carry buoyant features into the trench. Rosenbaum and Mo (2011) previously noted that this behaviour is typical of the western Pacific, where subduction of buoyant features often leads to locally hampered trench retreat and to slab steepening rather than slab shallowing.

Several other buoyant features in the Pacific and Indian Oceans are associated with “slab deforming” subduction, in which both a cusp in the trench and local shallowing of the slab dip is observed. Some cases correspond to smaller features on plates of young or intermediate age (as in our young plate models), and others to larger features on older plates (similar to the models of Flórez-Rodríguez et al., 2019). The most straightforward example of “slab deforming” buoyant-feature subduction may be the Cocos ridge, which coincides with a cusp and shallowing of the slab (Morell, 2015). Subduction of the Louisville Ridge below the Tonga trench introduces a minor cusp and some slab shallowing and contortion along the downdip projection of the obliquely subducting ridge (Bonnardot et al., 2007; Hayes et al., 2018). The Tonga-Kermadec trench is of moderate length, and might be expected to evolve to a “C” shaped trench (Schellart et al., 2007). However, the present-day shape of the Tonga-Kermadec trench is relatively straight, which may be due to early interaction between the middle of the trench, where it would have had the largest curvature, and the Louisville Ridge. This process would also account for the small size of the observed cusp where the ridge currently interacts with the Tonga-Kermadec trench. The cusp between the Japan and Kurile trench near Hokkaido is also associated with shallow slab below the island (Miller et al., 2006b; Kennett and Furumura, 2010; Hayes et al., 2018). Although there is no buoyant feature evident on the incoming plate, it has all the signatures of a recently subducted feature.

The interaction of the Ontong Java Plateau (OJP) with the South Pacific Trench is probably an extreme example of “slab deforming” style of interaction. The OJP separated from the Hikurangi plateau at ∼120 Myr and from the Manihiki plateau at ∼100 Myr (Taylor, 2006; Chandler et al., 2012; Seton et al., 2012) to reach the South Pacific trench at ∼20 Myr (Mann and Taira, 2004; Seton et al., 2012). The impinging of the OJP on the trench hampered trench retreat and deformed the trench which, until that time, had a typical concave shape similar to the present-day Japan trench (Mann and Taira, 2004; Seton et al., 2012). The large extent of the OJP made it too buoyant to subduct, perhaps acting similarly to the ridge in model HB_C_Young. As a result, subduction in the region ceased and eventually reversed polarity (Mann and Taira, 2004).

Slab deformation due to buoyant-feature subduction may expose parts of the slab to increased stresses, which result in slab tears (e.g. models by Mason et al., 2010). This may be the case for the subduction of the Emperor ridge at the cusp between the Kamchatka and Aleutian trench, which is associated with a low slab dip angle and a possible tear (Davaille and Lees, 2004). However, this cusp may predate subduction of the Emperor chain (Seton et al., 2012), and its formation may instead have been associated with the subduction of previous, arc-related, buoyant features (e.g., Vaes et al., 2019). Another example is the Palau-Kyushu ridge, located at the cusp between the Ryukyu and Nankai trenches, which has been proposed to contribute to shallow subduction along Nankai trench, possibly with a tear near the subducting ridge (Wu et al., 2016; Pownall et al., 2017; Xia et al., 2021). The low dip of the Nankai slab may also be affected by the fact that its northern end lies at the Boso triple junction, where trench motion is limited by interaction between the Nankai, Izu and Japan trenches. In the Indian ocean, subduction of the Investigator Fracture Zone and the extinct Wharton spreading ridge coincide with a cusp in the Sumatra trench, a local slab “kink” and possible slab tear (Hall and Spakman, 2015). The excess buoyancy is likely due to subduction of the Wharton ridge, as oceanic fracture zones like Investigator are not usually associated with changes to plate buoyancy (e.g., Jacob et al., 2014). Along the Java trench, where the Roo rise is starting to deform the trench (Kopp et al., 2006; Shulgin et al., 2011), it has been suggested that past subduction of similar buoyant features is responsible for a slab tear within the Java slab (Hall and Spakman, 2015).

Flat subduction is currently observed in several subduction zones in the Eastern Pacific, under Alaska, Mexico, and South America. In all cases, the strong westward motion of the Americas over the subducting Pacific, Cocos and Nazca plates may facilitate slab flattening (Van Hunen et al., 2002; Manea et al., 2017). In Alaska, the subduction of the Yakutat terrane occurs at a major cusp at the eastern end of the Aleutian trench and is associated with subduction of a narrow flat slab segment (Gulick et al., 2007). The Yakutat terrane subducts under a continental plate, although not of cratonic thickness, with a possible tear on the side (e.g., Daly et al., 2021). In South America, the Ecuador-Peru flat slab is associated with the subduction of the Carnegie and Nazca Ridges (which bound it from north and south) and the Chilean flat slab is associated with subduction of the Juan-Fernandez Ridge (Gutscher et al., 1999, 2000; Skinner and Clayton, 2013; Kumar et al., 2016; Bishop et al., 2017; Manea et al., 2017). The upper South American plate comprises several old platforms which may all have had thick roots, although not all persist to the present day (Feng et al., 2007). Although the Carnegie ridge interacts with the edge of the trench, where the trench motion is expected to be limited (Schellart et al., 2007), the associated Ecuador-Peru flat slab is located closer to the centre of the trench, as is the Chilean flat slab. Below the major cusp in the Andean trench, at the Bolivian Orocline, the slab is relatively steep. Although the cusp is in front of the subducting Iquique Ridge at present, it was formed before the arrival of the ridge to the trench, at 40–25 Ma (Allmendinger et al., 1997; Schepers et al., 2017). It may have formed above a previous flat slab segment which has been proposed to have existed between 40–20 Ma (Ramos and Folguera, 2009; Skinner and Clayton, 2013). A previously subducted buoyant feature may have formed the older flat slab and, aided by its location close to the centre of a long trench, have induced slab flattening. Capitanio et al. (2011) proposed that, in contrast to most trenches whose shape appears to be determined by subducting plate forcing, the recent ( < 20 Ma) rise of the Andes and corresponding shape of the South American trench may be controlled by variations in thickness of the upper plate. Finally, flat subduction below Mexico is generally thought to be too wide to be induced by the buoyant Tehuantepec ridge that subducts at its southern end (Gutscher et al., 2000; Skinner and Clayton, 2011). Factors that may contribute to the formation of the Mexican flat slab include the subduction of a spreading ridge near the northern end of the plateau, a continental overriding plate, and the flat slab location just north of the boundary between the North and South American plates, where shearing due to differential movement between the Americas may result in fast moving upper plate over slow and young subducting plate.