# Role of Mantle Drag on the Tectonics of Subduction Zones: Insights From Laboratory Models

Thomas Geffroy, B. Guillaume, M. Simoes, A. Replumaz, R. Lacassin, L.

Husson, J J Kermarrec

► To cite this version:

Thomas Geffroy, B. Guillaume, M. Simoes, A. Replumaz, R. Lacassin, et al. Role of Mantle Drag on the Tectonics of Subduction Zones: Insights From Laboratory Models. Tectonics, 2023, 42 (12), pp.e2023TC008018. 10.1029/2023tc008018. insu-04426143

HAL Id: insu-04426143

https://insu.hal.science/insu-04426143

Submitted on 30 Jan 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Role of Mantle Drag on the Tectonics of Subduction Zones: Insights From Laboratory Models

Thomas Geffroy, B. Guillaume, M. Simoes, A. Replumaz, R. Lacassin, L.

Husson, J J Kermarrec

► To cite this version:

Thomas Geffroy, B. Guillaume, M. Simoes, A. Replumaz, R. Lacassin, et al. Role of Mantle Drag on the Tectonics of Subduction Zones: Insights From Laboratory Models. Tectonics, 2023, 42 (12), pp.e2023TC008018. 10.1029/2023tc008018. insu-04426143

HAL Id: insu-04426143

https://insu.hal.science/insu-04426143v1

Submitted on 30 Jan 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH ARTICLE

10.1029/2023TC008018

**Key Points:**

* Our models show that mantle drag may exert a first-order control on subduction dynamics and upper plate tectonics

* Plate motion and overriding plate deformation linearly increases with the imposed mantle flow velocity in the models

* Mantle flow directed toward the trench favors upper plate trench-orthogonal shortening

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

T. Geffroy,
thomas.geffroy@univ-rennes.fr

**Citation:**

Geffroy, T., Guillaume, B., Simoes, M., Replumaz, A., Lacassin, R., Husson, L., & Kermarrec, J. J. (2023). Role of mantle drag on the tectonics of subduction zones: Insights from laboratory models. *Tectonics, 42*, e2023TC008018. https://doi.org/10.1029/2023TC008018

Received 19 JUL 2023
Accepted 18 NOV 2023

# Role of Mantle Drag on the Tectonics of Subduction Zones: Insights From Laboratory Models

<figure><img src="image_2.png" /></figure>

T. Geffroy¹ ID, B. Guillaume¹ ID, M. Simoes² ID, A. Replumaz³ ID, R. Lacassin² ID, L. Husson³ ID, and J. J. Kermarrec¹

¹CNRS, Géosciences Rennes, UMR 6118, Université de Rennes, Rennes, France, ²Institut de Physique du Globe de Paris, CNRS, Université Paris Cité, Paris, France, ³Institut des Sciences de la Terre (ISTerre), CNRS, USMB, IRD, IFSTTAR, Université Grenoble Alpes, Grenoble, France

**Abstract** Along convergent boundaries, the role played by mantle drag remains poorly understood despite its potential impact on subduction dynamics and in turn on the deformation regime of the overriding plate. In this study, we present 11 three-dimensional analog models of subduction including an overriding plate, in which mantle drag at the base of the lower or upper plate results from an imposed unidirectional horizontal mantle flow perpendicular to the trench, and in which the plate opposite to the flow is fixed. We varied the direction and the velocity of the imposed horizontal mantle flow between 0 and 10 cm/yr to quantify its impact on horizontal and vertical upper plate deformation, velocities of plates and subduction, and slab geometry. In our experiments, we show that a mantle flow lower than 5 cm/yr tends to laterally translate the slab rather than to generate internal deformation, resulting in limited differences in slab geometries between models. We also show that plate velocity correlates linearly with the imposed mantle flow velocity and associated mantle drag. The upper plate most often deforms by trench-orthogonal shortening, with shortening rates increasing linearly with mantle flow. Shortening rates are higher when mantle flow is directed toward the fixed upper plate and when the slab has not yet reached the upper-lower mantle discontinuity. Minimum trench-orthogonal shortening rates of $2.5 \times 10^{-15} \text{ s}^{-1}$ are required to thicken upper plates. This study suggests that mantle drag can exert first-order controls on the dynamics of subduction zones and associated tectonics.

**Plain Language Summary** The convective mantle and lithospheres interact to produce the motion of lithospheric plates at the Earth's surface. However, the time-evolution of their interactions remains to be fully characterized, particularly at subduction zones. Using three dimensional analog models, we test the role of mantle tractions on subduction dynamics by controlling the flow in the upper mantle. In our experiments, we show that the direction of the mantle flow has an impact on the deformation rates observed in the overriding plate, with faster deformation for a mantle flow directed toward the upper plate. The magnitude of the mantle flow also has an impact on the deformation of the overriding plate, with a linear increase of deformation rates with increasing mantle flow velocities. The arrival of the subducting plate at the upper-lower mantle discontinuity induces changes in the force equilibrium that leads to a decrease in plate velocity and upper plate deformation rates. From there, we propose that mantle drag is a key element controlling the time-evolution of subduction zones dynamics and the regime of deformation of the overriding plate.

## 1. Introduction

Tectonic plates and the underlying convective mantle interact with each other over geological time scales, leading to plate motion and deformation at plate boundaries, especially along subduction zones. It is generally accepted that subduction takes place owing to the negative buoyancy of the subducting lithosphere with respect to the surrounding mantle (Molnar & Gray, 1979). In turn, sinking of the slab displaces the mantle and generates a mantle flow that can ultimately drag the upper plate (Schellart & Moresi, 2013). In this case, the main force driving subduction and plate motion is the slab pull (Forsyth & Uyeda, 1975; Hager & O'Connell, 1979; Lithgow-Bertelloni & Richards, 1998). The subduction-driven convergence can lead to slab roll-back and trench retreat, resulting in trench-orthogonal stretching of the upper plate, as can be observed along many active subduction zones (e.g., Tonga, Ryukyu and Mariana; Lallemand et al., 2005). However, not all upper plates stretch above their subducting slabs. For instance, along the western margin of South America, the Andes formed by crustal shortening and thickening, indicating be taken into account. For instance, the absolute motion of both the subducting and overriding plates could exert a primary control on the regime of deformation of the upper plate
RESEARCH ARTICLE

10.1029/2023TC008018

**Key Points:**

* Our models show that mantle drag may exert a first-order control on subduction dynamics and upper plate tectonics

* Plate motion and overriding plate deformation linearly increases with the imposed mantle flow velocity in the models

* Mantle flow directed toward the trench favors upper plate trench-orthogonal shortening

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

T. Geffroy,
thomas.geffroy@univ-rennes.fr

**Citation:**

Geffroy, T., Guillaume, B., Simoes, M., Replumaz, A., Lacassin, R., Husson, L., & Kermarrec, J. J. (2023). Role of mantle drag on the tectonics of subduction zones: Insights from laboratory models. *Tectonics, 42*, e2023TC008018. https://doi.org/10.1029/2023TC008018

Received 19 JUL 2023
Accepted 18 NOV 2023

# Role of Mantle Drag on the Tectonics of Subduction Zones: Insights From Laboratory Models

<figure><img src="image_2.png" /></figure>

T. Geffroy¹ ID, B. Guillaume¹ ID, M. Simoes² ID, A. Replumaz³ ID, R. Lacassin² ID, L. Husson³ ID, and J. J. Kermarrec¹

¹CNRS, Géosciences Rennes, UMR 6118, Université de Rennes, Rennes, France, ²Institut de Physique du Globe de Paris, CNRS, Université Paris Cité, Paris, France, ³Institut des Sciences de la Terre (ISTerre), CNRS, USMB, IRD, IFSTTAR, Université Grenoble Alpes, Grenoble, France

**Abstract** Along convergent boundaries, the role played by mantle drag remains poorly understood despite its potential impact on subduction dynamics and in turn on the deformation regime of the overriding plate. In this study, we present 11 three-dimensional analog models of subduction including an overriding plate, in which mantle drag at the base of the lower or upper plate results from an imposed unidirectional horizontal mantle flow perpendicular to the trench, and in which the plate opposite to the flow is fixed. We varied the direction and the velocity of the imposed horizontal mantle flow between 0 and 10 cm/yr to quantify its impact on horizontal and vertical upper plate deformation, velocities of plates and subduction, and slab geometry. In our experiments, we show that a mantle flow lower than 5 cm/yr tends to laterally translate the slab rather than to generate internal deformation, resulting in limited differences in slab geometries between models. We also show that plate velocity correlates linearly with the imposed mantle flow velocity and associated mantle drag. The upper plate most often deforms by trench-orthogonal shortening, with shortening rates increasing linearly with mantle flow. Shortening rates are higher when mantle flow is directed toward the fixed upper plate and when the slab has not yet reached the upper-lower mantle discontinuity. Minimum trench-orthogonal shortening rates of $2.5 \times 10^{-15} \text{ s}^{-1}$ are required to thicken upper plates. This study suggests that mantle drag can exert first-order controls on the dynamics of subduction zones and associated tectonics.

**Plain Language Summary** The convective mantle and lithospheres interact to produce the motion of lithospheric plates at the Earth's surface. However, the time-evolution of their interactions remains to be fully characterized, particularly at subduction zones. Using three dimensional analog models, we test the role of mantle tractions on subduction dynamics by controlling the flow in the upper mantle. In our experiments, we show that the direction of the mantle flow has an impact on the deformation rates observed in the overriding plate, with faster deformation for a mantle flow directed toward the upper plate. The magnitude of the mantle flow also has an impact on the deformation of the overriding plate, with a linear increase of deformation rates with increasing mantle flow velocities. The arrival of the subducting plate at the upper-lower mantle discontinuity induces changes in the force equilibrium that leads to a decrease in plate velocity and upper plate deformation rates. From there, we propose that mantle drag is a key element controlling the time-evolution of subduction zones dynamics and the regime of deformation of the overriding plate.

## 1. Introduction

Tectonic plates and the underlying convective mantle interact with each other over geological time scales, leading to plate motion and deformation at plate boundaries, especially along subduction zones. It is generally accepted that subduction takes place owing to the negative buoyancy of the subducting lithosphere with respect to the surrounding mantle (Molnar & Gray, 1979). In turn, sinking of the slab displaces the mantle and generates a mantle flow that can ultimately drag the upper plate (Schellart & Moresi, 2013). In this case, the main force driving subduction and plate motion is the slab pull (Forsyth & Uyeda, 1975; Hager & O'Connell, 1979; Lithgow-Bertelloni & Richards, 1998). The subduction-driven convergence can lead to slab roll-back and trench retreat, resulting in trench-orthogonal stretching of the upper plate, as can be observed along many active subduction zones (e.g., Tonga, Ryukyu and Mariana; Lallemand et al., 2005). However, not all upper plates stretch above their subducting slabs. For instance, along the western margin of South America, the Andes formed by crustal shortening and thickening, indicating be taken into account. For instance, the absolute motion of both the subducting and overriding plates could exert a primary control on the regime of deformation of the upper plate
<figure><img src="image_3.png" /><figcaption>**Figure 1.** (a) Global map of the present-day mean of integrated mantle drags from both the upper and lower plates (annotated Σmd/2, colored dots), together with upper plate strain indicators (UPS) from Heuret and Lallemand (2005) going from highly extensive (−3) to highly compressive (+3) (colored bars); modified from Husson (2012). (b) Cross-section showing the geometry of the lower and upper plates in the case of the Tonga-Kermadec subduction zone (see location in panel (a)). Lithospheric thicknesses are extracted from the LITHO1.0 model (Pasyanos et al., 2014) and the slab geometry is from the tomographic model GAP-P4 (Obayashi et al., 2013). Black arrows indicate the direction of the integrated mantle drag beneath the upper and lower plates, and scale with its magnitude (after Husson (2012)). Red arrows give an indication of the direction and the magnitude of the horizontal trench displacement, using the GSRM-APM-1 reference frame (Kreemer, 2009). Panels (c) and (d) same as (b) for the Kuril and the South American subduction zones, respectively.</figcaption></figure>

(e.g., Cerpa et al., 2018; Guillaume et al., 2009, 2018; Schellart, 2005; Schepers et al., 2017). It should be noted, however, that previous studies based on analog or numerical models generally imposed a constant plate velocity, independent of the force equilibrium in the system, and therefore could potentially overestimate or underestimate the effect of absolute plate velocity. Other studies have argued for a combination of slab pull and pressure-driven asthenosphere flow to drive plate motion (Coltice et al., 2019; Höink et al., 2011; Semple & Lenardic, 2020). Coltice et al. (2019) suggested that for Earth-like models, around 35% on average of the surface may be dragged by the interior, therefore possibly forcing plate motion at subduction zones. Husson (2012) also proposed that the convective mantle could drag the base of the lithospheres, and that the resulting net trenchward force underneath both plates would ultimately control the upper plate deformation regime at subduction zones (Figure 1a). Upper plate stretching would preferentially occur where the net mean trenchward force is low or negative (e.g., Tonga-Kermadec subduction zone, Figure 1b). Instead, upper plate shortening would be found where the net mean trenchward force is highly positive. The latter situation could result either from a large trenchward mantle drag under the subducting plate and a smaller trenchward mantle drag under the upper plate (e.g., Kuril subduction zone, Figure 1c), or from large trenchward mantle drag on both sides of the subduction zone (e.g., South America subduction zone, Figure 1d).

Subduction zones are, in fact, not isolated systems, and the sub-lithospheric mantle drag may indeed not only result from the local mantle flow induced by slab sinking but may also be modulated by other processes at all spatial scales. The ascent of mantle plumes can induce flow in the mantle (e.g., Cande & Stegman, 2011; Obrebski et al., 2010) and have an impact on nearby subduction zones, as can be seen in the Tonga-Kermadec subduction zone with the presence of the Samoa plume close to the subduction (e.g., Chang et al., 2016). Interactions between close-by subduction zones may also modify mantle drag (Di Leo et al., 2012; Király et al., 2018; Peral et al., 2018). Some studies also proposed that an eastward directed global mantle flow could develop owing to the net westward displacement of plates at the surface and would be the cause for geometrical and kinematic differences between subduction zones verging to the west and those verging to the east (Crespi et al., 2007; Cuffaro & Doglioni, 2007; Ricard et al., 1991).

The role exerted by mantle drag on subduction dynamics and plate deformation remains poorly constrained, because of the large spatial and temporal scales of the involved processes (several hundreds km and several Myr)
<figure><img src="image_3.png" /></figure>

**Figure 1.** (a) Global map of the present-day mean of integrated mantle drags from both the upper and lower plates (annotated Σmd/2, colored dots), together with upper plate strain indicators (UPS) from Heuret and Lallemand (2005) going from highly extensive (−3) to highly compressive (+3) (colored bars); modified from Husson (2012). (b) Cross-section showing the geometry of the lower and upper plates in the case of the Tonga-Kermadec subduction zone (see location in panel (a)). Lithospheric thicknesses are extracted from the LITHO1.0 model (Pasyanos et al., 2014) and the slab geometry is from the tomographic model GAP-P4 (Obayashi et al., 2013). Black arrows indicate the direction of the integrated mantle drag beneath the upper and lower plates, and scale with its magnitude (after Husson (2012)). Red arrows give an indication of the direction and the magnitude of the horizontal trench displacement, using the GSRM-APM-1 reference frame (Kreemer, 2009). Panels (c) and (d) same as (b) for the Kuril and the South American subduction zones, respectively.

(e.g., Cerpa et al., 2018; Guillaume et al., 2009, 2018; Schellart, 2005; Schepers et al., 2017). It should be noted, however, that previous studies based on analog or numerical models generally imposed a constant plate velocity, independent of the force equilibrium in the system, and therefore could potentially overestimate or underestimate the effect of absolute plate velocity. Other studies have argued for a combination of slab pull and pressure-driven asthenosphere flow to drive plate motion (Coltice et al., 2019; Höink et al., 2011; Semple & Lenardic, 2020). Coltice et al. (2019) suggested that for Earth-like models, around 35% on average of the surface may be dragged by the interior, therefore possibly forcing plate motion at subduction zones. Husson (2012) also proposed that the convective mantle could drag the base of the lithospheres, and that the resulting net trenchward force underneath both plates would ultimately control the upper plate deformation regime at subduction zones (Figure 1a). Upper plate stretching would preferentially occur where the net mean trenchward force is low or negative (e.g., Tonga-Kermadec subduction zone, Figure 1b). Instead, upper plate shortening would be found where the net mean trenchward force is highly positive. The latter situation could result either from a large trenchward mantle drag under the subducting plate and a smaller trenchward mantle drag under the upper plate (e.g., Kuril subduction zone, Figure 1c), or from large trenchward mantle drag on both sides of the subduction zone (e.g., South America subduction zone, Figure 1d).

Subduction zones are, in fact, not isolated systems, and the sub-lithospheric mantle drag may indeed not only result from the local mantle flow induced by slab sinking but may also be modulated by other processes at all spatial scales. The ascent of mantle plumes can induce flow in the mantle (e.g., Cande & Stegman, 2011; Obrebski et al., 2010) and have an impact on nearby subduction zones, as can be seen in the Tonga-Kermadec subduction zone with the presence of the Samoa plume close to the subduction (e.g., Chang et al., 2016). Interactions between close-by subduction zones may also modify mantle drag (Di Leo et al., 2012; Király et al., 2018; Peral et al., 2018). Some studies also proposed that an eastward directed global mantle flow could develop owing to the net westward displacement of plates at the surface and would be the cause for geometrical and kinematic differences between subduction zones verging to the west and those verging to the east (Crespi et al., 2007; Cuffaro & Doglioni, 2007; Ricard et al., 1991).

The role exerted by mantle drag on subduction dynamics and plate deformation remains poorly constrained, because of the large spatial and temporal scales of the involved processes (several hundreds km and several Myr)
that can only be partially documented by geophysical and geological observations. Modeling efforts have been made to better understand retroactions between slab dynamics, mantle flow and upper plate deformation (e.g., Conrad & Behn, 2010; Conrad & Lithgow-Bertelloni, 2002; Espurt et al., 2008; Guillaume et al., 2010; Hertgen et al., 2020; Holt et al., 2015; Menant et al., 2016; Schellart, 2010). However, these models were generally not specifically designed to test the effect of mantle drag on plates and slab displacement. Plate displacement either resulted from subduction-driven mantle flow, or was forced by imposing constant velocities so that in this case it did not result from the equilibrium of forces within the system. The few 2D modeling studies that included a background mantle flow during subduction showed that slabs exhibit low dip when there is a forward flow (i.e., directed toward the upper plate), while steeper slabs are obtained when the flow is backward (i.e., directed toward the subducting plate) (Boutelier & Cruden, 2008; Ficini et al., 2017; Rodríguez-González et al., 2014). These studies also suggest that contrasted tectonic regimes can develop in the upper plate in correlation with varying slab geometries, a steeper slab being more likely associated with back-arc basins (Ficini et al., 2017). However, by conducting both 2-D models and 2.5-D models, one of the key components of mantle flow - namely the lateral (toroidal) flow - is neglected. Toroidal flow allows the mantle to move around the subduction zone, and can have important consequences on trench migration and upper plate deformation (Funiciello et al., 2003; Kincaid & Griffiths, 2003; Piromallo et al., 2006). Previous few 3D models have shown that the addition of flow in the mantle has little impact, on the geometry of the slab in the vertical plane, in the case of flow directed perpendicularly to the slab (Chertova et al., 2018; Guillaume et al., 2021), which is at odds with the results from 2D models. This highlights the necessity to include the toroidal component of mantle flow in modeling, as toroidal flow has been shown to be active around present-day subduction zones (e.g., Palano et al., 2017). Importantly, these various models did not include an overriding plate and when they did so, did not specifically investigate its internal deformation in light of subduction dynamics and boundary conditions. Recently, Habel et al. (2023) performed a series of 3D analog models of subduction zones, imposing a constant unidirectional mantle flow at 10 cm/yr and adding an overriding continental plate, with the idea of exploring the boundary conditions that lead to upper plate shortening and thickening. They found that these conditions are met and favored in the case of a mantle flow directed toward a fixed upper continental plate, while the absence of flow led to upper plate trench-orthogonal stretching. However, the systematics of these findings in terms of flow directions and velocities were not explored.

Taking into account the results and limitations of previous models, we considered a situation where the flow in the mantle is not simply the result of the slab's sinking but can be externally driven. Considering these boundary conditions, we aimed at constraining the dependence of upper plate strain regime and amount of deformation with mantle flow velocity and direction. We also aimed at evaluating if there is a systematic relationship between the applied mantle flow velocity and direction on one hand, and the geometry of the slab and plate kinematics on the other. To answer these questions, we built upon the study by Habel et al. (2023) by using a similar experimental set-up. Previous studies, based on models considering density-driven buoyancy variations determined from seismic investigations at the global scale, suggest that horizontal present-day upper mantle flow velocities could reach values up to 5 cm/yr (Becker et al., 2008; Bull et al., 2010). Additionally, global geodynamic models (Weismüller et al., 2015), and regional models considering Poiseuille flow in the asthenosphere resulting from plume activity (Colli et al., 2014; Stotz et al., 2023) propose that horizontal asthenosphere flow could reach velocities exceeding 10 cm/yr. In our study, we decided to systematically explore scaled velocities up to 10 cm/yr to account for possible uncertainties on present-day upper mantle horizontal velocity and variations of mantle flow in the past, by performing a total of 11 models. The plate opposite to the imposed flow was also fixed with respect to a lower mantle reference frame to avoid a global drift of the system that would have led to the interaction of the plates with the box boundaries. We monitored horizontal plate displacements, subduction velocity, slab geometry and upper plate vertical and horizontal deformation to evaluate how imposed boundary conditions controlled the motion and deformation of plates and slabs.

## 2. Materials and Methods

### 2.1. Experimental Setup

The experimental setup is adapted from Habel et al. (2023). We use a Plexiglas tank of 100 × 100 × 30 cm³ filled with glucose syrup. An intermediate Plexiglas panel with a thickness of 1 cm is fixed to the walls of the tank to simulate the discontinuity between the lower and upper mantles, considered here as impermeable. This panel
that can only be partially documented by geophysical and geological observations. Modeling efforts have been made to better understand retroactions between slab dynamics, mantle flow and upper plate deformation (e.g., Conrad & Behn, 2010; Conrad & Lithgow-Bertelloni, 2002; Espurt et al., 2008; Guillaume et al., 2010; Hertgen et al., 2020; Holt et al., 2015; Menant et al., 2016; Schellart, 2010). However, these models were generally not specifically designed to test the effect of mantle drag on plates and slab displacement. Plate displacement either resulted from subduction-driven mantle flow, or was forced by imposing constant velocities so that in this case it did not result from the equilibrium of forces within the system. The few 2D modeling studies that included a background mantle flow during subduction showed that slabs exhibit low dip when there is a forward flow (i.e., directed toward the upper plate), while steeper slabs are obtained when the flow is backward (i.e., directed toward the subducting plate) (Boutelier & Cruden, 2008; Ficini et al., 2017; Rodríguez-González et al., 2014). These studies also suggest that contrasted tectonic regimes can develop in the upper plate in correlation with varying slab geometries, a steeper slab being more likely associated with back-arc basins (Ficini et al., 2017). However, by conducting both 2-D models and 2.5-D models, one of the key components of mantle flow - namely the lateral (toroidal) flow - is neglected. Toroidal flow allows the mantle to move around the subduction zone, and can have important consequences on trench migration and upper plate deformation (Funiciello et al., 2003; Kincaid & Griffiths, 2003; Piromallo et al., 2006). Previous few 3D models have shown that the addition of flow in the mantle has little impact, on the geometry of the slab in the vertical plane, in the case of flow directed perpendicularly to the slab (Chertova et al., 2018; Guillaume et al., 2021), which is at odds with the results from 2D models. This highlights the necessity to include the toroidal component of mantle flow in modeling, as toroidal flow has been shown to be active around present-day subduction zones (e.g., Palano et al., 2017). Importantly, these various models did not include an overriding plate and when they did so, did not specifically investigate its internal deformation in light of subduction dynamics and boundary conditions. Recently, Habel et al. (2023) performed a series of 3D analog models of subduction zones, imposing a constant unidirectional mantle flow at 10 cm/yr and adding an overriding continental plate, with the idea of exploring the boundary conditions that lead to upper plate shortening and thickening. They found that these conditions are met and favored in the case of a mantle flow directed toward a fixed upper continental plate, while the absence of flow led to upper plate trench-orthogonal stretching. However, the systematics of these findings in terms of flow directions and velocities were not explored.

Taking into account the results and limitations of previous models, we considered a situation where the flow in the mantle is not simply the result of the slab's sinking but can be externally driven. Considering these boundary conditions, we aimed at constraining the dependence of upper plate strain regime and amount of deformation with mantle flow velocity and direction. We also aimed at evaluating if there is a systematic relationship between the applied mantle flow velocity and direction on one hand, and the geometry of the slab and plate kinematics on the other. To answer these questions, we built upon the study by Habel et al. (2023) by using a similar experimental set-up. Previous studies, based on models considering density-driven buoyancy variations determined from seismic investigations at the global scale, suggest that horizontal present-day upper mantle flow velocities could reach values up to 5 cm/yr (Becker et al., 2008; Bull et al., 2010). Additionally, global geodynamic models (Weismüller et al., 2015), and regional models considering Poiseuille flow in the asthenosphere resulting from plume activity (Colli et al., 2014; Stotz et al., 2023) propose that horizontal asthenosphere flow could reach velocities exceeding 10 cm/yr. In our study, we decided to systematically explore scaled velocities up to 10 cm/yr to account for possible uncertainties on present-day upper mantle horizontal velocity and variations of mantle flow in the past, by performing a total of 11 models. The plate opposite to the imposed flow was also fixed with respect to a lower mantle reference frame to avoid a global drift of the system that would have led to the interaction of the plates with the box boundaries. We monitored horizontal plate displacements, subduction velocity, slab geometry and upper plate vertical and horizontal deformation to evaluate how imposed boundary conditions controlled the motion and deformation of plates and slabs.

## 2. Materials and Methods

### 2.1. Experimental Setup

The experimental setup is adapted from Habel et al. (2023). We use a Plexiglas tank of 100 × 100 × 30 cm³ filled with glucose syrup. An intermediate Plexiglas panel with a thickness of 1 cm is fixed to the walls of the tank to simulate the discontinuity between the lower and upper mantles, considered here as impermeable. This panel
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Experimental setup seen from above (top) and from the side (bottom) for models with mantle flow directed toward the subducting plate ("backward flow," left), and toward the overriding plate ("forward flow," right). Dimensions are given in cm. In side views, the black horizontal line represents the lower boundary of the upper mantle.</figcaption></figure>

is 16 cm shorter than the width of the tank, allowing the syrup to flow up or down the sides of the setup. The syrup is set in motion by a moving piston under the intermediate Plexiglas panel, using a cord-pulley system connected to a step-motor. The movement of the piston is therefore opposite to that of the flow in the upper mantle (Figure 2). We vary the direction of the flow in the upper mantle relative to the vergence of the modeled subduction, testing two configurations: (a) a first one wherein the flow in the upper mantle is directed toward the upper plate (hereafter referred to as "forward flow"), with a free lower plate and a fixed upper plate (Figure 2, right) and (b) a second one in which the flow is directed toward the subducting plate (hereafter referred to as "backward flow"), with a free overriding plate and a fixed subducting plate (Figure 2, left). In both cases, we fixed one of the two plates to avoid a global drift of the system that would have ended with the interaction of the plates with the lateral boundaries of the tank. Here, the plates are placed far enough from the sides of the tank to avoid edge effects. Tested scaled velocities are 1.66, 3.33, 5, 7.5 and 10 cm/yr for each flow direction.

## 2.2. Materials

The Plexiglass tank is filled with glucose syrup (Cargill 01497), which is an analog material for the sub-lithospheric upper mantle. Glucose syrup is sensitive to temperature variations, as its viscosity follows an Arrhenius type function (Schellart, 2011). For the recorded laboratory temperatures, ranging from 15.4 to 19.1°C during the entire set of experiments, viscosity varied between 209 ± 5 and 453 ± 5 Pa·s, with an average value of 335 ± 70 Pa·s (Table 1). Glucose syrup is a Newtonian fluid and its viscosity does not vary with strain rate (Schellart, 2011). This is a simplification of the mantle rheology in nature as upper mantle is reported to deform under combined effects of diffusion (Newtonian) and dislocation (non-Newtonian) creep, their activation depending on stress, temperature, pressure, or grain size (e.g., Karato, 2010). Our approach therefore only considers diffusion creep in the sub-lithospheric upper mantle. The density of the syrup was measured to be 1,420 ± 5 kg/m³, and the temperature dependence of the density is considered insignificant given the range of laboratory temperatures.
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Experimental setup seen from above (top) and from the side (bottom) for models with mantle flow directed toward the subducting plate ("backward flow," left), and toward the overriding plate ("forward flow," right). Dimensions are given in cm. In side views, the black horizontal line represents the lower boundary of the upper mantle.</figcaption></figure>

is 16 cm shorter than the width of the tank, allowing the syrup to flow up or down the sides of the setup. The syrup is set in motion by a moving piston under the intermediate Plexiglas panel, using a cord-pulley system connected to a step-motor. The movement of the piston is therefore opposite to that of the flow in the upper mantle (Figure 2). We vary the direction of the flow in the upper mantle relative to the vergence of the modeled subduction, testing two configurations: (a) a first one wherein the flow in the upper mantle is directed toward the upper plate (hereafter referred to as "forward flow"), with a free lower plate and a fixed upper plate (Figure 2, right) and (b) a second one in which the flow is directed toward the subducting plate (hereafter referred to as "backward flow"), with a free overriding plate and a fixed subducting plate (Figure 2, left). In both cases, we fixed one of the two plates to avoid a global drift of the system that would have ended with the interaction of the plates with the lateral boundaries of the tank. Here, the plates are placed far enough from the sides of the tank to avoid edge effects. Tested scaled velocities are 1.66, 3.33, 5, 7.5 and 10 cm/yr for each flow direction.

## 2.2. Materials

The Plexiglass tank is filled with glucose syrup (Cargill 01497), which is an analog material for the sub-lithospheric upper mantle. Glucose syrup is sensitive to temperature variations, as its viscosity follows an Arrhenius type function (Schellart, 2011). For the recorded laboratory temperatures, ranging from 15.4 to 19.1°C during the entire set of experiments, viscosity varied between 209 ± 5 and 453 ± 5 Pa·s, with an average value of 335 ± 70 Pa·s (Table 1). Glucose syrup is a Newtonian fluid and its viscosity does not vary with strain rate (Schellart, 2011). This is a simplification of the mantle rheology in nature as upper mantle is reported to deform under combined effects of diffusion (Newtonian) and dislocation (non-Newtonian) creep, their activation depending on stress, temperature, pressure, or grain size (e.g., Karato, 2010). Our approach therefore only considers diffusion creep in the sub-lithospheric upper mantle. The density of the syrup was measured to be 1,420 ± 5 kg/m³, and the temperature dependence of the density is considered insignificant given the range of laboratory temperatures.
<table><caption>**Table 1**

Experimental Parameters, Boundary Conditions, and Scaling for the Eleven Models</caption><caption>Tableau 1</caption><thead><tr><th></th><th></th><th></th><th></th><th>SA06</th><th>SA09</th><th>SA10</th><th>SA12</th><th>SA13</th><th>SA14</th><th>SA19</th><th>SA20</th><th>SA21</th><th>SA22</th><th>SA24</th></tr></thead><tbody><tr><td rowspan="3">Lower plate</td><td>Thickness</td><td><i>h</i></td><td>mm</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td></tr><tr><td>Density</td><td>ρ</td><td>kg/m³</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td></tr><tr><td>Viscosity</td><td>η</td><td>Pa·s</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td></tr><tr><td rowspan="3">Upper plate</td><td>Thickness</td><td><i>h</i></td><td>mm</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td></tr><tr><td>Density</td><td>ρ</td><td>kg/m³</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td></tr><tr><td>Viscosity</td><td>η</td><td>Pa·s</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td></tr><tr><td rowspan="4">Sub-lithospheric mantle</td><td>Thickness</td><td><i>h</i></td><td>mm</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td></tr><tr><td>Density</td><td>ρ</td><td>kg/m³</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td></tr><tr><td>Viscosity</td><td>η</td><td>Pa·s</td><td>392</td><td>360</td><td>305</td><td>299</td><td>293</td><td>453</td><td>286</td><td>281</td><td>242</td><td>286</td><td>210</td></tr><tr><td>Temperature</td><td><i>T</i></td><td>°C</td><td>16.1</td><td>16.5</td><td>17.3</td><td>17.4</td><td>17.5</td><td>15.4</td><td>17.6</td><td>17.7</td><td>18.4</td><td>17.6</td><td>19.1</td></tr><tr><td rowspan="3">Scaling</td><td>Length</td><td>1 km</td><td>mm</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td></tr><tr><td>Time</td><td>1 Ma</td><td>s</td><td>164</td><td>150</td><td>127</td><td>125</td><td>122</td><td>189</td><td>119</td><td>117</td><td>101</td><td>119</td><td>88</td></tr><tr><td>Velocity</td><td>10 cm/yr</td><td>mm/min</td><td>5.6</td><td>6.1</td><td>7.1</td><td>7.3</td><td>7.4</td><td>4.8</td><td>7.6</td><td>7.7</td><td>9</td><td>7.6</td><td>10.4</td></tr><tr><td rowspan="2">Boundary conditions</td><td>Mantle velocity</td><td><i>v</i><sub>f</sub></td><td>cm/yr</td><td>—</td><td>1.66</td><td>3.33</td><td>5</td><td>10</td><td>7.5</td><td>-1.66</td><td>-3.33</td><td>-5</td><td>-7.5</td><td>-10</td></tr><tr><td>Fixed plate</td><td></td><td></td><td>—</td><td>Upper</td><td>Upper</td><td>Upper</td><td>Upper</td><td>Upper</td><td>Lower</td><td>Lower</td><td>Lower</td><td>Lower</td><td>Lower</td></tr></tbody></table>

Note. Uncertainties are of ±0.2 mm for thickness measurements, and of ±5 kg/m³ for densities. The imposed mantle flow velocity (boundary conditions) is reported for each one of the models, with scaled (and experimental) values.

The oceanic lithosphere is represented by a 1.4-cm thick plate made of PDMS silicone filled with iron powder. The density of the filled silicone was measured at 1,480 ± 5 kg/m³ and considered independent of temperature. The density difference between the glucose syrup and the silicone is therefore -60 kg/m³ and is considered to be fixed over time. The viscosity of the silicone is constant at around 99,800 ± 100 Pa·s for the temperature values recorded in the laboratory (Fernández-García, 2019) (Table 1). The silicone has a Newtonian rheology for strain rates lower than 10⁻² s⁻¹ (Rudolf et al., 2016), which is the case in our experiments, where strain rates do not exceed 10⁻³ s⁻¹. The viscosity ratio between the silicone plate and the glucose syrup in our models therefore ranges between 220 and 477, which is within the range of values for the ratio between slab and mantle in natural subduction zones estimated in between 50 and 500 (e.g., Funiciello et al., 2008; Loiselet et al., 2009; Wu et al., 2008).

The overriding continental lithosphere is represented by a 0.8-cm thick plate also made of PDMS silicone filled with iron powder. The viscosity of the silicone was measured to be 90,000 ± 100 Pa·s for the entire set of models and its density is 1,300 ± 5 kg/m³, corresponding to a density difference of +120 kg/m³ with respect to the glucose syrup (Table 1).

Both plates are decoupled by a thin layer (~1 mm) of a mix of 20% petrolatum and 80% paraffin oil deposited on top of the subducting plate, following the approach of Duarte et al. (2013). This material has a viscosity that is strain-rate and temperature dependent. For the estimated strain rates (10⁻³-10⁻² s⁻¹) and measured laboratory temperatures, it viscosity varied between around 200 and 2,000 Pa·s (after Duarte et al., 2013).

Given the rectangular shape of the plates, the subduction trench is initially linear at the onset of subduction. It then evolves freely according to the subduction dynamics.

## 2.3. Scaling

The experiments are properly scaled for gravity, length, density, viscosity, and velocity following the method of previous studies (e.g., Guillaume et al., 2021). The length scale ratio considered between laboratory and nature is 1.52 × 10⁻⁷, that is, 0.15 mm in the model corresponds to 1 km in nature. The initial thickness of the continental
<table><caption>**Table 1**
Experimental Parameters, Boundary Conditions, and Scaling for the Eleven Models</caption><caption>Tableau 1</caption><thead><tr><th></th><th></th><th></th><th></th><th>SA06</th><th>SA09</th><th>SA10</th><th>SA12</th><th>SA13</th><th>SA14</th><th>SA19</th><th>SA20</th><th>SA21</th><th>SA22</th><th>SA24</th></tr></thead><tbody><tr><td rowspan="3">Lower plate</td><td>Thickness</td><td><i>h</i></td><td>mm</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td><td>14</td></tr><tr><td>Density</td><td>ρ</td><td>kg/m<sup>3</sup></td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td><td>1,480</td></tr><tr><td>Viscosity</td><td>η</td><td>Pa·s</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td><td>99,800</td></tr><tr><td rowspan="3">Upper plate</td><td>Thickness</td><td><i>h</i></td><td>mm</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td><td>8</td></tr><tr><td>Density</td><td>ρ</td><td>kg/m<sup>3</sup></td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td><td>1,300</td></tr><tr><td>Viscosity</td><td>η</td><td>Pa·s</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td><td>90,000</td></tr><tr><td rowspan="4">Sub-lithospheric mantle</td><td>Thickness</td><td><i>h</i></td><td>mm</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td><td>100</td></tr><tr><td>Density</td><td>ρ</td><td>kg/m<sup>3</sup></td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td><td>1,420</td></tr><tr><td>Viscosity</td><td>η</td><td>Pa·s</td><td>392</td><td>360</td><td>305</td><td>299</td><td>293</td><td>453</td><td>286</td><td>281</td><td>242</td><td>286</td><td>210</td></tr><tr><td>Temperature</td><td><i>T</i></td><td>°C</td><td>16.1</td><td>16.5</td><td>17.3</td><td>17.4</td><td>17.5</td><td>15.4</td><td>17.6</td><td>17.7</td><td>18.4</td><td>17.6</td><td>19.1</td></tr><tr><td rowspan="3">Scaling</td><td>Length</td><td>1 km</td><td>mm</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td><td>0.15</td></tr><tr><td>Time</td><td>1 Ma</td><td>s</td><td>164</td><td>150</td><td>127</td><td>125</td><td>122</td><td>189</td><td>119</td><td>117</td><td>101</td><td>119</td><td>88</td></tr><tr><td>Velocity</td><td>10 cm/yr</td><td>mm/min</td><td>5.6</td><td>6.1</td><td>7.1</td><td>7.3</td><td>7.4</td><td>4.8</td><td>7.6</td><td>7.7</td><td>9</td><td>7.6</td><td>10.4</td></tr><tr><td rowspan="2">Boundary conditions</td><td>Mantle velocity</td><td><i>v</i><sub>f</sub></td><td>cm/yr</td><td>—</td><td>1.66</td><td>3.33</td><td>5</td><td>10</td><td>7.5</td><td>-1.66</td><td>-3.33</td><td>-5</td><td>-7.5</td><td>-10</td></tr><tr><td>Fixed plate</td><td></td><td></td><td>—</td><td>Upper</td><td>Upper</td><td>Upper</td><td>Upper</td><td>Upper</td><td>Lower</td><td>Lower</td><td>Lower</td><td>Lower</td><td>Lower</td></tr></tbody></table>

Note. Uncertainties are of ±0.2 mm for thickness measurements, and of ±5 kg/m³ for densities. The imposed mantle flow velocity (boundary conditions) is reported for each one of the models, with scaled (and experimental) values.

The oceanic lithosphere is represented by a 1.4-cm thick plate made of PDMS silicone filled with iron powder. The density of the filled silicone was measured at 1,480 ± 5 kg/m³ and considered independent of temperature. The density difference between the glucose syrup and the silicone is therefore -60 kg/m³ and is considered to be fixed over time. The viscosity of the silicone is constant at around 99,800 ± 100 Pa·s for the temperature values recorded in the laboratory (Fernández-García, 2019) (Table 1). The silicone has a Newtonian rheology for strain rates lower than 10⁻² s⁻¹ (Rudolf et al., 2016), which is the case in our experiments, where strain rates do not exceed 10⁻³ s⁻¹. The viscosity ratio between the silicone plate and the glucose syrup in our models therefore ranges between 220 and 477, which is within the range of values for the ratio between slab and mantle in natural subduction zones estimated in between 50 and 500 (e.g., Funiciello et al., 2008; Loiselet et al., 2009; Wu et al., 2008).

The overriding continental lithosphere is represented by a 0.8-cm thick plate also made of PDMS silicone filled with iron powder. The viscosity of the silicone was measured to be 90,000 ± 100 Pa·s for the entire set of models and its density is 1,300 ± 5 kg/m³, corresponding to a density difference of +120 kg/m³ with respect to the glucose syrup (Table 1).

Both plates are decoupled by a thin layer (~1 mm) of a mix of 20% petrolatum and 80% paraffin oil deposited on top of the subducting plate, following the approach of Duarte et al. (2013). This material has a viscosity that is strain-rate and temperature dependent. For the estimated strain rates (10⁻³-10⁻² s⁻¹) and measured laboratory temperatures, it viscosity varied between around 200 and 2,000 Pa·s (after Duarte et al., 2013).

Given the rectangular shape of the plates, the subduction trench is initially linear at the onset of subduction. It then evolves freely according to the subduction dynamics.

## 2.3. Scaling

The experiments are properly scaled for gravity, length, density, viscosity, and velocity following the method of previous studies (e.g., Guillaume et al., 2021). The length scale ratio considered between laboratory and nature is 1.52 × 10⁻⁷, that is, 0.15 mm in the model corresponds to 1 km in nature. The initial thickness of the continental
plate is 0.8 cm, which corresponds to ~53 km in nature (Table 1). We therefore consider that the upper plate is relatively “weak,” which allows us to better observe and quantify deformation in our experiments, as the amount of deformation scales with the plate thickness in our viscous models. The thickness of the oceanic plate is 1.4 cm, which is equivalent to ~92 km in nature. Considering the relationship between the age of the oceanic lithosphere t and its thickness h, established as $h = 2.32\sqrt{\kappa t}$ (Turcotte & Schubert, 2002), we obtain an age of 50 Ma for a thermal diffusivity $\kappa = 10^{-6}$ m²/s, which roughly corresponds to the median age of oceanic lithospheres at present-day subduction zones (after Müller et al. (2016)).

Scaling for time between models and nature is achieved following:

$$ \frac{t_n}{t_m} = \frac{(\Delta \rho g h)_m}{(\Delta \rho g h)_n} \frac{\eta_n}{\eta_m} \qquad (1) $$

The indexes n and m are associated to the parameters in nature or models, respectively, $\Delta \rho$ is the density contrast between the lithosphere and the sub-lithospheric mantle, and $\eta$ is the sub-lithospheric mantle viscosity. The effective viscosity of the Earth's mantle is thought to vary between $10^{19}$ and $10^{22}$ Pa·s depending on depth (e.g., King, 2016; Mitrovica & Forte, 2004). We consider here an average value for the entire sub-lithospheric upper mantle of $5 \times 10^{20}$ Pa·s. This means that 1 Myr in nature corresponds on average to a duration of $129 \pm 29$ s in the laboratory (Table 1). A velocity of 10 cm/yr in nature (maximum imposed velocity for the mantle flow) therefore corresponds on average to a velocity of $7.3 \pm 1.5$ mm/min in the laboratory.

Given that the viscosity of the sub-lithospheric mantle can vary significantly from one model to the other, in particular for the models conducted at the maximum and minimum laboratory temperatures, we choose to present time and velocities scaled to their natural equivalents for an easier inter-comparison of model results. Following this, for convenience, we also choose to describe our experiments using the nature-equivalent geodynamic terms (mantle, lithospheric plate) rather than their laboratory analogs (glucose syrup, silicone putties).

## 2.4. Forces at Work

The simplified subduction system is designed to impose a simple force balance, approximating that of plate tectonics (Figure 3). In the absence of external flow, the only driving force for subduction in the model is the slab pull ($F_B$), which depends on the density difference between the slab and the underlying mantle, and on the volume of slab sinking into the mantle. Activation of a background mantle flow generates an additional viscous force, namely the mantle drag (md) that applies at the base of the lithosphere, which contributes to plate displacement and favors subduction. When mantle flow is activated, we fix either the lower or the upper plate to a back wall. By doing this, we impose a no-motion boundary condition to the plate within the model box reference frame, generating a reaction force perpendicular to the wall ($F_r$). This reaction force can be assimilated to some extent to a ridge push force and/or to mantle drag in the opposite direction, even if not explicitly modeled (Habel et al., 2023). For models with backward flow (fixed subducting plate), this force constantly adjusts to balance the trenchward force that applies to the subducting plate, making $F_r$ a resisting force (Figure 3, blue). For models with forward flow (fixed upper plate), $F_r$ prevents the upper plate from moving away from the trench, increasing coupling between plates, and as such should also be considered as a resisting force (Figure 3, red).

Subduction is resisted by a set of forces. One of these forces is the slab bending force ($F_b$), which depends on the minimum radius of curvature of the bending slab, the viscosity of the slab, the thickness of the plate and the subduction velocity. Additional resistive forces include the viscous forces ($F_v$) exerted by the mantle on the slab, which have no analytical definition but depend on the vertical and horizontal velocities of the slab.

Friction at the interface between the two plates ($F_f$) also slows down subduction. Previous studies suggest that the subduction interface must be mechanically weak in the long term to match observed plate and convergence velocities (e.g., Behr et al., 2022; Duarte et al., 2014, 2015; Zhong & Gurnis, 1994).

When the slab approaches the bottom of the box, two additional forces contribute to slow down subduction. The first is the basal resistance force ($F_{BR}$), which prevents the slab from going further down. The second is the basal frictional force ($F_F$), which results from the resistance to the lateral motion of the slab at the interface between the plate and the bottom of the box. Subduction can therefore be separated into two stages: a first subduction stage until the slab reaches the bottom of the box, hereafter referred to as stage 1 (Figure 3, top), and a second subduction stage that takes place after the slab has reached the bottom of the box, hereafter referred to as stage 2 (Figure 3, bottom).
plate is 0.8 cm, which corresponds to ~53 km in nature (Table 1). We therefore consider that the upper plate is relatively “weak,” which allows us to better observe and quantify deformation in our experiments, as the amount of deformation scales with the plate thickness in our viscous models. The thickness of the oceanic plate is 1.4 cm, which is equivalent to ~92 km in nature. Considering the relationship between the age of the oceanic lithosphere t and its thickness h, established as $h = 2.32\sqrt{\kappa t}$ (Turcotte & Schubert, 2002), we obtain an age of 50 Ma for a thermal diffusivity $\kappa = 10^{-6}$ m²/s, which roughly corresponds to the median age of oceanic lithospheres at present-day subduction zones (after Müller et al. (2016)).

Scaling for time between models and nature is achieved following:

$$ \frac{t_n}{t_m} = \frac{(\Delta \rho g h)_m}{(\Delta \rho g h)_n} \frac{\eta_n}{\eta_m} \qquad (1) $$

The indexes n and m are associated to the parameters in nature or models, respectively, Δρ is the density contrast between the lithosphere and the sub-lithospheric mantle, and η is the sub-lithospheric mantle viscosity. The effective viscosity of the Earth's mantle is thought to vary between 10¹⁹ and 10²² Pa·s depending on depth (e.g., King, 2016; Mitrovica & Forte, 2004). We consider here an average value for the entire sub-lithospheric upper mantle of 5 × 10²⁰ Pa·s. This means that 1 Myr in nature corresponds on average to a duration of 129 ± 29 s in the laboratory (Table 1). A velocity of 10 cm/yr in nature (maximum imposed velocity for the mantle flow) therefore corresponds on average to a velocity of 7.3 ± 1.5 mm/min in the laboratory.

Given that the viscosity of the sub-lithospheric mantle can vary significantly from one model to the other, in particular for the models conducted at the maximum and minimum laboratory temperatures, we choose to present time and velocities scaled to their natural equivalents for an easier inter-comparison of model results. Following this, for convenience, we also choose to describe our experiments using the nature-equivalent geodynamic terms (mantle, lithospheric plate) rather than their laboratory analogs (glucose syrup, silicone putties).

## 2.4. Forces at Work

The simplified subduction system is designed to impose a simple force balance, approximating that of plate tectonics (Figure 3). In the absence of external flow, the only driving force for subduction in the model is the slab pull ($F_B$), which depends on the density difference between the slab and the underlying mantle, and on the volume of slab sinking into the mantle. Activation of a background mantle flow generates an additional viscous force, namely the mantle drag (md) that applies at the base of the lithosphere, which contributes to plate displacement and favors subduction. When mantle flow is activated, we fix either the lower or the upper plate to a back wall. By doing this, we impose a no-motion boundary condition to the plate within the model box reference frame, generating a reaction force perpendicular to the wall ($F_r$). This reaction force can be assimilated to some extent to a ridge push force and/or to mantle drag in the opposite direction, even if not explicitly modeled (Habel et al., 2023). For models with backward flow (fixed subducting plate), this force constantly adjusts to balance the trenchward force that applies to the subducting plate, making $F_r$ a resisting force (Figure 3, blue). For models with forward flow (fixed upper plate), $F_r$ prevents the upper plate from moving away from the trench, increasing coupling between plates, and as such should also be considered as a resisting force (Figure 3, red).

Subduction is resisted by a set of forces. One of these forces is the slab bending force ($F_b$), which depends on the minimum radius of curvature of the bending slab, the viscosity of the slab, the thickness of the plate and the subduction velocity. Additional resistive forces include the viscous forces ($F_v$) exerted by the mantle on the slab, which have no analytical definition but depend on the vertical and horizontal velocities of the slab.

Friction at the interface between the two plates ($F_f$) also slows down subduction. Previous studies suggest that the subduction interface must be mechanically weak in the long term to match observed plate and convergence velocities (e.g., Behr et al., 2022; Duarte et al., 2014, 2015; Zhong & Gurnis, 1994).

When the slab approaches the bottom of the box, two additional forces contribute to slow down subduction. The first is the basal resistance force ($F_{BR}$), which prevents the slab from going further down. The second is the basal frictional force ($F_F$), which results from the resistance to the lateral motion of the slab at the interface between the plate and the bottom of the box. Subduction can therefore be separated into two stages: a first subduction stage until the slab reaches the bottom of the box, hereafter referred to as stage 1 (Figure 3, top), and a second subduction stage that takes place after the slab has reached the bottom of the box, hereafter referred to as stage 2 (Figure 3, bottom).
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Forces acting on the modeled subduction system during the two stages of subduction. The balance of forces described here applies only to a vertical plane. The details of the viscous forces exerted by the mantle, in particular those related to the three-dimensional flow that affects it, are not represented here for simplicity.</figcaption></figure>

The balance of forces described here applies only to a vertical plane. While some of these forces (slab pull, slab bending, friction at the interface) vary little laterally in our models due to the quasi-cylindricity of subduction, this is not the case for all forces. In particular, the viscous forces exerted by the mantle can vary laterally due to the three-dimensional flow (poloidal and toroidal components) that affects it.

## 2.5. Experiments Description and Measurements

When active, mantle flow is started a few minutes before subduction initiation to ensure a steady flow. Subduction is initiated manually, by pushing the leading edge of the subducting plate into the mantle down to a depth of about 5–7 cm (equivalent to 330–460 km in nature). The conditions leading to subduction initiation are therefore not studied in this work. The overriding plate is set just after subduction initiation to control the distance between both plates, which is initially fixed at around 5 mm (equivalent to 33 km). The gap between the two plates only exists during the first minutes of subduction and gradually diminishes as the experiment evolves autonomously in response to both internal and external applied conditions.

We monitor the experiment by taking pictures every 3 min with DSLR cameras (Nikon D-3300) with a resolution of 6000 × 4000. One camera is located on one side of the experiment to observe the geometry of the slab at one edge, and one camera is placed on top of the model to follow the evolution of the subduction at the surface (Figure S2 in Supporting Information S1).

From the lateral view, we measure the time-evolution of slab horizontal velocity at upper mantle middepth ($v_m$) and the horizontal velocity of the slab tip ($v_{st(x)}$) with measurement uncertainty equivalent to 0.3–0.6 cm/yr. From the top view and along a central section in the direction orthogonal to the trench, we measure the time-evolution
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Forces acting on the modeled subduction system during the two stages of subduction. The balance of forces described here applies only to a vertical plane. The details of the viscous forces exerted by the mantle, in particular those related to the three-dimensional flow that affects it, are not represented here for simplicity.</figcaption></figure>

The balance of forces described here applies only to a vertical plane. While some of these forces (slab pull, slab bending, friction at the interface) vary little laterally in our models due to the quasi-cylindricity of subduction, this is not the case for all forces. In particular, the viscous forces exerted by the mantle can vary laterally due to the three-dimensional flow (poloidal and toroidal components) that affects it.

## 2.5. Experiments Description and Measurements

When active, mantle flow is started a few minutes before subduction initiation to ensure a steady flow. Subduction is initiated manually, by pushing the leading edge of the subducting plate into the mantle down to a depth of about 5–7 cm (equivalent to 330–460 km in nature). The conditions leading to subduction initiation are therefore not studied in this work. The overriding plate is set just after subduction initiation to control the distance between both plates, which is initially fixed at around 5 mm (equivalent to 33 km). The gap between the two plates only exists during the first minutes of subduction and gradually diminishes as the experiment evolves autonomously in response to both internal and external applied conditions.

We monitor the experiment by taking pictures every 3 min with DSLR cameras (Nikon D-3300) with a resolution of 6000 × 4000. One camera is located on one side of the experiment to observe the geometry of the slab at one edge, and one camera is placed on top of the model to follow the evolution of the subduction at the surface (Figure S2 in Supporting Information S1).

From the lateral view, we measure the time-evolution of slab horizontal velocity at upper mantle middepth ($v_m$) and the horizontal velocity of the slab tip ($v_{st(x)}$) with measurement uncertainty equivalent to 0.3–0.6 cm/yr. From the top view and along a central section in the direction orthogonal to the trench, we measure the time-evolution
of subduction velocity ($v_s$), trench velocity ($v_t$), velocity of the free plate $v_{fp}$ (lower plate for forward mantle flow or upper plate for backward mantle flow), and upper plate length changes by visually tracking the position of points corresponding to the plate's edges along the section, over successive frames. Velocities are reported within the model box reference frame for all models, that is, equivalent to a fixed lower mantle reference frame. Measurement uncertainties vary from 0.075 cm/yr to 0.15 cm/yr, and upper plate trench-orthogonal stretching/ shortening ($\epsilon_h$) have measurement uncertainties on the order of 0.15%-0.3%. We also deduce the mean vertical deformation of the upper plate ($\epsilon_v$) by measuring changes in its surface, assuming that the plate is incompressible, with a measurement uncertainty on the order of 0.2%. Time-evolution of velocities are smoothed for graphical representation (Figure S1 in Supporting Information S1) using a simple moving average (window of 6 time units), but mean values and associated standard deviations presented hereafter for each stage are calculated from raw data. We let the experiments run for a duration equivalent to 50 Myr.

# 3. Results

## 3.1. Initial Model (SA06): Free Plates and No Mantle Flow

The initial model (SA06) has been run to determine the intrinsic kinematic and geometric parameters of the subduction system in the absence of mantle flow as plate geometries (width, length, thickness) and physical (density) and rheological (viscosity) properties are identical to those in other models with mantle flow. In this model, both plates are left free. Subduction in this model follows a classic evolution, as largely described in previous studies (e.g., Espurt et al., 2008; Funiciello et al., 2003; Guillaume et al., 2021; Xue et al., 2022).

During the first stage of the experiment, the slab sinks into the upper mantle with a rollback motion until it reaches the 660-km discontinuity after ~9 Myr (Figure 4). Afterward, during stage 2, the slab continues to roll back at velocities of $1 \pm 0.55$ cm/yr at upper mantle middepth until the end of the experiment after 45 Myr (Figure 4). During this second stage, the slab dip progressively decreases from a nearly vertical slab to a more shallowly dipping slab (Figure 4). At the same time, the trench retreats at a rate of $1.14 \pm 0.31$ cm/yr that is slightly faster than the absolute velocity of the overriding plate ($1.11 \pm 0.22$ cm/yr), resulting in the finite trench-orthogonal stretching (by 1.85%) and thinning (by 1.1%) of the upper continental plate after 45 Myr (Figure 5).

During the first stage, the oceanic plate moves toward the overriding plate at a velocity of $2.11 \pm 1.23$ cm/yr, while in the second stage it remains almost fixed ($-0.03 \pm 0.4$ cm/yr) (Figure S1 in Supporting Information S1 and Figure 6c). Subduction velocity also decreases after the arrival of the slab at the 660-km discontinuity from $2.87 \pm 1.27$ cm/yr to $1.11 \pm 0.49$ cm/yr (Figure S1 in Supporting Information S1 and Figure 6d).

## 3.2. Forward Mantle Flow

A horizontal forward mantle flow is imposed, with velocities (considered positive in this direction) varying between +1.66 cm/yr (model SA09) and +10 cm/yr (model SA13). Models with intermediate velocities of +3.33 cm/yr (model SA10), +5 cm/yr (model SA12), and +7.5 cm/yr (model SA14) were additionally performed. In this series of models, the upper plate is fixed to a backwall, that is, it has a zero velocity relative to the bottom of the box.

For all the models, the slab arrives at the 660 km discontinuity with an almost identical steep angle (Figure 4). Instead, during the stage 2 of subduction, the geometry of the slab significantly varies as a function of mantle flow velocity. In particular, for mantle flow velocities $> +3.33$ cm/yr, the shallow part of the slab tends to flatten. For the maximum mantle flow velocity (+10 cm/yr), the slab even exhibits a roll-over geometry (Figure 4).

The absolute velocity of the free plate - here the oceanic plate - increases almost linearly with increasing flow velocity for both subduction stages. During stage 1, it varies from $1.27 \pm 0.43$ cm/yr for the slowest mantle flow to $5.19 \pm 0.68$ cm/yr for the fastest mantle flow ($R^2 = 0.98$). It accounts for ~50% of the imposed background mantle flow velocity (Figure 6c). During stage 2, it varies from $0.41 \pm 0.26$ cm/yr to $2.67 \pm 1.02$ cm/yr ($R^2 = 0.86$), accounting for ~25% of the imposed flow velocity (Figure 6c). At depth, horizontal displacements of the slab and slab tip also linearly correlate with imposed mantle flow during both stages (Figures 6e and 6f). It is only for the model with mantle flow at +10 cm/yr that the slab tip velocity does not increase during stage 2, and it is the only model in which the slab rolls over. The subduction velocity linearly increases with mantle flow during stage 1 (from $0.86 \pm 0.38$ to $2.37 \pm 0.99$ cm/yr, $R^2 = 0.97$) while there is no clear systematic variation during stage 2 (Figure 6d).
of subduction velocity ($v_s$), trench velocity ($v_t$), velocity of the free plate $v_{fp}$ (lower plate for forward mantle flow or upper plate for backward mantle flow), and upper plate length changes by visually tracking the position of points corresponding to the plate's edges along the section, over successive frames. Velocities are reported within the model box reference frame for all models, that is, equivalent to a fixed lower mantle reference frame. Measurement uncertainties vary from 0.075 cm/yr to 0.15 cm/yr, and upper plate trench-orthogonal stretching/ shortening ($\epsilon_h$) have measurement uncertainties on the order of 0.15%-0.3%. We also deduce the mean vertical deformation of the upper plate ($\epsilon_v$) by measuring changes in its surface, assuming that the plate is incompressible, with a measurement uncertainty on the order of 0.2%. Time-evolution of velocities are smoothed for graphical representation (Figure S1 in Supporting Information S1) using a simple moving average (window of 6 time units), but mean values and associated standard deviations presented hereafter for each stage are calculated from raw data. We let the experiments run for a duration equivalent to 50 Myr.

# 3. Results

## 3.1. Initial Model (SA06): Free Plates and No Mantle Flow

The initial model (SA06) has been run to determine the intrinsic kinematic and geometric parameters of the subduction system in the absence of mantle flow as plate geometries (width, length, thickness) and physical (density) and rheological (viscosity) properties are identical to those in other models with mantle flow. In this model, both plates are left free. Subduction in this model follows a classic evolution, as largely described in previous studies (e.g., Espurt et al., 2008; Funiciello et al., 2003; Guillaume et al., 2021; Xue et al., 2022).

During the first stage of the experiment, the slab sinks into the upper mantle with a rollback motion until it reaches the 660-km discontinuity after ~9 Myr (Figure 4). Afterward, during stage 2, the slab continues to roll back at velocities of $1 \pm 0.55$ cm/yr at upper mantle middepth until the end of the experiment after 45 Myr (Figure 4). During this second stage, the slab dip progressively decreases from a nearly vertical slab to a more shallowly dipping slab (Figure 4). At the same time, the trench retreats at a rate of $1.14 \pm 0.31$ cm/yr that is slightly faster than the absolute velocity of the overriding plate ($1.11 \pm 0.22$ cm/yr), resulting in the finite trench-orthogonal stretching (by 1.85%) and thinning (by 1.1%) of the upper continental plate after 45 Myr (Figure 5).

During the first stage, the oceanic plate moves toward the overriding plate at a velocity of $2.11 \pm 1.23$ cm/yr, while in the second stage it remains almost fixed ($-0.03 \pm 0.4$ cm/yr) (Figure S1 in Supporting Information S1 and Figure 6c). Subduction velocity also decreases after the arrival of the slab at the 660-km discontinuity from $2.87 \pm 1.27$ cm/yr to $1.11 \pm 0.49$ cm/yr (Figure S1 in Supporting Information S1 and Figure 6d).

## 3.2. Forward Mantle Flow

A horizontal forward mantle flow is imposed, with velocities (considered positive in this direction) varying between +1.66 cm/yr (model SA09) and +10 cm/yr (model SA13). Models with intermediate velocities of +3.33 cm/yr (model SA10), +5 cm/yr (model SA12), and +7.5 cm/yr (model SA14) were additionally performed. In this series of models, the upper plate is fixed to a backwall, that is, it has a zero velocity relative to the bottom of the box.

For all the models, the slab arrives at the 660 km discontinuity with an almost identical steep angle (Figure 4). Instead, during the stage 2 of subduction, the geometry of the slab significantly varies as a function of mantle flow velocity. In particular, for mantle flow velocities $> +3.33$ cm/yr, the shallow part of the slab tends to flatten. For the maximum mantle flow velocity (+10 cm/yr), the slab even exhibits a roll-over geometry (Figure 4).

The absolute velocity of the free plate - here the oceanic plate - increases almost linearly with increasing flow velocity for both subduction stages. During stage 1, it varies from $1.27 \pm 0.43$ cm/yr for the slowest mantle flow to $5.19 \pm 0.68$ cm/yr for the fastest mantle flow ($R^2 = 0.98$). It accounts for ~50% of the imposed background mantle flow velocity (Figure 6c). During stage 2, it varies from $0.41 \pm 0.26$ cm/yr to $2.67 \pm 1.02$ cm/yr ($R^2 = 0.86$), accounting for ~25% of the imposed flow velocity (Figure 6c). At depth, horizontal displacements of the slab and slab tip also linearly correlate with imposed mantle flow during both stages (Figures 6e and 6f). It is only for the model with mantle flow at +10 cm/yr that the slab tip velocity does not increase during stage 2, and it is the only model in which the slab rolls over. The subduction velocity linearly increases with mantle flow during stage 1 (from $0.86 \pm 0.38$ to $2.37 \pm 0.99$ cm/yr, $R^2 = 0.97$) while there is no clear systematic variation during stage 2 (Figure 6d).
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Slab geometry and horizontal position of the plates when subduction is initiated (blue), when the slab arrives at the 660-km discontinuity (green), and after ~50 Myr of subduction (orange), interpreted from lateral (left) and top (right) views of the 11 models. Stage 1 corresponds to the evolution of subduction from initiation (blue) to when the slab reaches the 660-km discontinuity (green), while stage 2 encompasses the subsequent time until the end of the experiment (orange). Initial model SA06 is reported on top, while other models, with a fixed plate and with backward mantle flow (left half of the figure) or forward mantle flow (right half of the figure), are reported below. The down-pointed triangles on lateral views represent the trench position at different experimental times, and pin-points indicate when a plate is fixed.</figcaption></figure>

The overriding plate always records finite trench-orthogonal shortening, except during the second stage of the model with the lowest mantle flow velocity ($v_f = +1.66 \text{ cm/yr}$) as the deformation regime becomes neutral (Figure 5a, right panel, and Figure 7a). During stage 1, the shortening rate linearly increases with mantle flow velocities from $-0.6 \times 10^{-16} \text{ to } -4.5 \times 10^{-16} \text{ s}^{-1}$ ($R^2 = 0.97$) (Figure 7a). During stage 2, the same linear relation holds, albeit at slower rates, from $-3.2 \times 10^{-18} \text{ to } -1.9 \times 10^{-16} \text{ s}^{-1}$. The final cumulative amount of shortening ranges between -2.8% and -39.7% (Figure 5a, right panel).

The increase in mantle flow velocity also has an impact on the vertical deformation of the continental plate (Figure 5b right panel). During stage 1, trench-orthogonal shortening is accompanied by thickening at rates varying from $1.9 \times 10^{-17} \text{ s}^{-1}$ for $v_f = +1.66 \text{ cm/yr}$ to $3.5 \times 10^{-16} \text{ s}^{-1}$ for $v_f = +10 \text{ cm/yr}$ (Figure 7b). During stage 2, the thickening rate also increases with mantle flow velocity, albeit at slower rates too, from $-1.3 \times 10^{-17} \text{ to } 1.9 \times 10^{-16} \text{ s}^{-1}$. For all models, thickening rates are systematically lower than shortening rates, evidencing that deformation of the upper plate is also accommodated by trench-parallel stretching during subduction. The impact of trench-parallel stretching remains nevertheless limited. We only observe upper plate limited finite thinning
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Slab geometry and horizontal position of the plates when subduction is initiated (blue), when the slab arrives at the 660-km discontinuity (green), and after ~50 Myr of subduction (orange), interpreted from lateral (left) and top (right) views of the 11 models. Stage 1 corresponds to the evolution of subduction from initiation (blue) to when the slab reaches the 660-km discontinuity (green), while stage 2 encompasses the subsequent time until the end of the experiment (orange). Initial model SA06 is reported on top, while other models, with a fixed plate and with backward mantle flow (left half of the figure) or forward mantle flow (right half of the figure), are reported below. The down-pointed triangles on lateral views represent the trench position at different experimental times, and pin-points indicate when a plate is fixed.</figcaption></figure>

The overriding plate always records finite trench-orthogonal shortening, except during the second stage of the model with the lowest mantle flow velocity ($v_f = +1.66$ cm/yr) as the deformation regime becomes neutral (Figure 5a, right panel, and Figure 7a). During stage 1, the shortening rate linearly increases with mantle flow velocities from $-0.6 \times 10^{-16}$ to $-4.5 \times 10^{-16}$ s$^{-1}$ ($R^2 = 0.97$) (Figure 7a). During stage 2, the same linear relation holds, albeit at slower rates, from $-3.2 \times 10^{-18}$ to $-1.9 \times 10^{-16}$ s$^{-1}$. The final cumulative amount of shortening ranges between -2.8% and -39.7% (Figure 5a, right panel).

The increase in mantle flow velocity also has an impact on the vertical deformation of the continental plate (Figure 5b right panel). During stage 1, trench-orthogonal shortening is accompanied by thickening at rates varying from $1.9 \times 10^{-17}$ s$^{-1}$ for $v_f = +1.66$ cm/yr to $3.5 \times 10^{-16}$ s$^{-1}$ for $v_f = +10$ cm/yr (Figure 7b). During stage 2, the thickening rate also increases with mantle flow velocity, albeit at slower rates too, from $-1.3 \times 10^{-17}$ to $1.9 \times 10^{-16}$ s$^{-1}$. For all models, thickening rates are systematically lower than shortening rates, evidencing that deformation of the upper plate is also accommodated by trench-parallel stretching during subduction. The impact of trench-parallel stretching remains nevertheless limited. We only observe upper plate limited finite thinning
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Time-evolution of upper plate deformation. Models with backward mantle flow and forward mantle flow are illustrated to the left and right, respectively. Each colored curve corresponds to a single model, as indicated at the top, with the corresponding equivalent mantle flow velocity given in cm/yr (sign convention as in Figure 2). The down-pointed triangles represent the transition from stage 1 to stage 2, when the slab reaches the 660-km interface. (a) Cumulative amount of upper plate trench-orthogonal deformation measured in the center of the plate (positive values for stretching, negative values for shortening). (b) Cumulative amount of upper plate vertical deformation measured from changes in the measured top surface of the upper plate (positive values for thickening, negative for thinning).</figcaption></figure>

(final amount of -0.6%) in the model with the slowest mantle velocity of +1.66 cm/yr (Figure 5b, right panel), despite upper plate trench-orthogonal shortening (final amount of -2.8%).

### 3.3. Backward Mantle Flow

We apply a mantle flow in the opposite direction, that is, a backward mantle flow, with velocities varying between -1.66 cm/yr (model SA19) and -10 cm/yr (model SA24) (velocities considered negative in this direction, Figure 2b). Three intermediate models with mantle flow velocities of -3.33 cm/yr (model SA20), -5 cm/yr (model SA21) and -7.5 cm/yr (model SA22) were performed. In these models, the oceanic plate is fixed to a backwall, that is, it has a zero velocity relative to the bottom of the box.

For all these models, the slab arrives at the 660 km discontinuity with an almost identical vertical geometry (Figure 4). At that time, the dip of the slab is nevertheless slightly higher than for equivalent models with opposite mantle flow (Figure 4). During stage 2, unlike for models with forward mantle flow, there are no clear changes of slab geometry with varying mantle flow velocity.

For this series of models, the upper plate is the free plate. The absolute velocity of the free plate increases linearly with increasing mantle flow velocity, but at significantly lower rates than equivalent models with mantle flow in the opposite direction (Figure 6c). During stage 1, upper plate velocities linearly correlate with mantle flow velocities, with values ranging between -0.42 ± 0.19 cm/yr for the slowest model and -1.59 ± 0.32 cm/yr for the fastest model (R² = 0.94) (Figure 6c), accounting on average for ~15% of the imposed flow (vs. 50% when the flow is in the opposite direction). For stage 2, upper plate velocity remains correlated with mantle flow velocity
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Time-evolution of upper plate deformation. Models with backward mantle flow and forward mantle flow are illustrated to the left and right, respectively. Each colored curve corresponds to a single model, as indicated at the top, with the corresponding equivalent mantle flow velocity given in cm/yr (sign convention as in Figure 2). The down-pointed triangles represent the transition from stage 1 to stage 2, when the slab reaches the 660-km interface. (a) Cumulative amount of upper plate trench-orthogonal deformation measured in the center of the plate (positive values for stretching, negative values for shortening). (b) Cumulative amount of upper plate vertical deformation measured from changes in the measured top surface of the upper plate (positive values for thickening, negative for thinning).</figcaption></figure>

(final amount of -0.6%) in the model with the slowest mantle velocity of +1.66 cm/yr (Figure 5b, right panel), despite upper plate trench-orthogonal shortening (final amount of -2.8%).

### 3.3. Backward Mantle Flow

We apply a mantle flow in the opposite direction, that is, a backward mantle flow, with velocities varying between -1.66 cm/yr (model SA19) and -10 cm/yr (model SA24) (velocities considered negative in this direction, Figure 2b). Three intermediate models with mantle flow velocities of -3.33 cm/yr (model SA20), -5 cm/yr (model SA21) and -7.5 cm/yr (model SA22) were performed. In these models, the oceanic plate is fixed to a backwall, that is, it has a zero velocity relative to the bottom of the box.

For all these models, the slab arrives at the 660 km discontinuity with an almost identical vertical geometry (Figure 4). At that time, the dip of the slab is nevertheless slightly higher than for equivalent models with opposite mantle flow (Figure 4). During stage 2, unlike for models with forward mantle flow, there are no clear changes of slab geometry with varying mantle flow velocity.

For this series of models, the upper plate is the free plate. The absolute velocity of the free plate increases linearly with increasing mantle flow velocity, but at significantly lower rates than equivalent models with mantle flow in the opposite direction (Figure 6c). During stage 1, upper plate velocities linearly correlate with mantle flow velocities, with values ranging between -0.42 ± 0.19 cm/yr for the slowest model and -1.59 ± 0.32 cm/yr for the fastest model (R² = 0.94) (Figure 6c), accounting on average for ~15% of the imposed flow (vs. 50% when the flow is in the opposite direction). For stage 2, upper plate velocity remains correlated with mantle flow velocity
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Evolution or variation of measured parameters in the different implemented models. (a) Time-evolution of subducting plate velocity for the initial model SA06 explaining how mean values are calculated for each of the two stages (before and after slab interaction with the 660-km interface between the upper and lower mantle) in the following diagrams. Values are averaged for the entire duration of stage 1 (blue) and stage 2 (orange). Reported error bars indicate the standard deviation. (b) Mean trench velocity ($v_t$) for the two stages of subduction as a function of the imposed mantle flow velocity. Positive values indicate trench motion toward the subducting plate, that is, trench retreat. (c) Mean free plate velocity ($v_{fp}$) (velocity of lower plate for positive mantle flow velocities, and of upper plate for negative mantle flow velocities). Positive velocity indicates plate motion toward the upper plate. (d) Mean subduction velocity. (e) Mean slab horizontal velocity measured at 330-km depth ($v_{st}$). Positive velocity indicates slab motion toward the upper plate. (f) Mean slab tip horizontal velocity ($v_{st(x)}$). Positive velocity indicates slab tip motion toward the upper plate.</figcaption></figure>

but exhibits lower rates with values ranging between −0.26 ± 0.13 cm/yr and −0.54 ± 0.2 cm/yr, accounting on average for ~3% of the imposed velocity (vs. 25% when the flow is in the opposite direction).

At depth, horizontal displacements of the slab at upper mantle middepth and slab tip are correlated with imposed mantle flow during the first stage with $R^2 = 0.74$ and $R^2 = 0.66$, respectively (Figures 6e and 6f). After the slab interacts with the 660-km discontinuity, these velocities decrease down to absolute values lower than 0.5 cm/yr. For both series of models with opposite mantle flow direction, upper mantle middepth and slab tip velocities are correlated with the free plate velocity, suggesting that stresses associated with mantle flow favor the horizontal displacement of the slab rather than internal deformation.
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Evolution or variation of measured parameters in the different implemented models. (a) Time-evolution of subducting plate velocity for the initial model SA06 explaining how mean values are calculated for each of the two stages (before and after slab interaction with the 660-km interface between the upper and lower mantle) in the following diagrams. Values are averaged for the entire duration of stage 1 (blue) and stage 2 (orange). Reported error bars indicate the standard deviation. (b) Mean trench velocity ($v_t$) for the two stages of subduction as a function of the imposed mantle flow velocity. Positive values indicate trench motion toward the subducting plate, that is, trench retreat. (c) Mean free plate velocity ($v_{fp}$) (velocity of lower plate for positive mantle flow velocities, and of upper plate for negative mantle flow velocities). Positive velocity indicates plate motion toward the upper plate. (d) Mean subduction velocity. (e) Mean slab horizontal velocity measured at 330-km depth ($v_{st}$). Positive velocity indicates slab motion toward the upper plate. (f) Mean slab tip horizontal velocity ($v_{st(x)}$). Positive velocity indicates slab tip motion toward the upper plate.</figcaption></figure>

but exhibits lower rates with values ranging between −0.26 ± 0.13 cm/yr and −0.54 ± 0.2 cm/yr, accounting on average for ~3% of the imposed velocity (vs. 25% when the flow is in the opposite direction).

At depth, horizontal displacements of the slab at upper mantle middepth and slab tip are correlated with imposed mantle flow during the first stage with $R^2 = 0.74$ and $R^2 = 0.66$, respectively (Figures 6e and 6f). After the slab interacts with the 660-km discontinuity, these velocities decrease down to absolute values lower than 0.5 cm/yr. For both series of models with opposite mantle flow direction, upper mantle middepth and slab tip velocities are correlated with the free plate velocity, suggesting that stresses associated with mantle flow favor the horizontal displacement of the slab rather than internal deformation.
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Deformation of the upper plate in the different models. (a) Mean upper plate trench-orthogonal horizontal deformation rate ($\epsilon_h$) for the two stages of subduction as a function of the imposed mantle flow velocity (positive values indicate stretching). (b) Mean upper plate vertical deformation rate ($\epsilon_v$) as a function of mantle flow velocity (positive values indicate thickening). Values are averaged for the entire duration of stage 1 (blue symbols) and stage 2 (orange symbols), error bars indicate the standard deviation (Figure 6a).</figcaption></figure>

For subduction velocity, we observe the same behavior as for models with opposite mantle flow, that is, the subduction velocity linearly increases with mantle flow during stage 1 (from $0.42 \pm 0.17$ to $1.34 \pm 0.41$ cm/yr, $R^2 = 0.94$) while there is no clear systematic variation during stage 2 (Figure 6d).

The upper plate records trench-orthogonal deformation, but at rates slower than when mantle flow is directed toward the upper plate. The finite deformation after 50 Myr varies between weak stretching (+0.07%) for the model with the slowest mantle flow velocity, to shortening for the models with faster mantle flow velocities $|v_f| \ge 3.33$ cm/yr, with amounts of shortening reaching up to -4.6% for the model with $|v_f| = 10$ cm/yr (Figure 5a). These values are almost one order of magnitude lower than those observed with mantle flow in the opposite direction. Deformation rates correlate with mantle flow velocity ($R^2 = 0.91$ for stage 1, $R^2 = 0.97$ for stage 2), and are systematically lower for stage 2 than those for stage 1 (Figure 7a).

There is no thickening of the continental plate associated with shortening of the upper plate, even for the model with the largest horizontal finite shortening (-4.6% for model with $v_f = -10$ cm/yr) (Figure 7b). Models with fast mantle flow velocities show no vertical deformation, while the finite amount of thinning after 50 Myr is only of -3.5% for the model with the slowest mantle flow (Figure 5b). This is in agreement with previous models with forward mantle flow, where finite thickening of the upper plate was only observed for a finite amount of shortening larger than -5% (Figure 5).

# 4. Discussion

## 4.1. Evaluating Mantle-Plate Coupling

Our models include viscous coupling between the sub-lithospheric mantle and the plates. This coupling occurs at the base of each plate along the horizontal portions of both plates, but also along the subducted part of the lower plate. The addition of a background mantle flow modifies this coupling, which has strong implications on the kinematics and dynamics of the subduction system.

At the base of the lithosphere, the total shear force induced by the mantle is directly proportional to its viscosity and to the velocity gradient. Given that in our models the mantle has a Newtonian rheology and thus a strain-rate independent viscosity, the traction exerted at the base of the plate should increase linearly with increasing mantle velocity, considering that the distance over which the gradient applies does not change. To explore the sensitivity of the free plate motion to imposed mantle flow velocity, we define $v_{fp}/v_f$, the ratio between the free plate and the mantle flow velocities. In the absence of subduction, in our linearly viscous system, $v_{fp}/v_f$ is found to be 1 (complete plate-mantle coupling, experiment not shown in this manuscript). When subduction and background flow are active, $v_{fp}/v_f$ drops below 1 with values varying between 0.05 and 0.77 depending on the subduction stage and on the flow direction (Figure 8). For models with forward flow, one could assume that $v_{fp}/v_f$ would reach values larger than 1, as the subducting plate is not only dragged toward the trench by the imposed mantle flow,
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Deformation of the upper plate in the different models. (a) Mean upper plate trench-orthogonal horizontal deformation rate ($\epsilon_h$) for the two stages of subduction as a function of the imposed mantle flow velocity (positive values indicate stretching). (b) Mean upper plate vertical deformation rate ($\epsilon_v$) as a function of mantle flow velocity (positive values indicate thickening). Values are averaged for the entire duration of stage 1 (blue symbols) and stage 2 (orange symbols), error bars indicate the standard deviation (Figure 6a).</figcaption></figure>

For subduction velocity, we observe the same behavior as for models with opposite mantle flow, that is, the subduction velocity linearly increases with mantle flow during stage 1 (from $0.42 \pm 0.17$ to $1.34 \pm 0.41$ cm/yr, $R^2 = 0.94$) while there is no clear systematic variation during stage 2 (Figure 6d).

The upper plate records trench-orthogonal deformation, but at rates slower than when mantle flow is directed toward the upper plate. The finite deformation after 50 Myr varies between weak stretching (+0.07%) for the model with the slowest mantle flow velocity, to shortening for the models with faster mantle flow velocities $|v_f| \ge 3.33$ cm/yr, with amounts of shortening reaching up to -4.6% for the model with $|v_f| = 10$ cm/yr (Figure 5a). These values are almost one order of magnitude lower than those observed with mantle flow in the opposite direction. Deformation rates correlate with mantle flow velocity ($R^2 = 0.91$ for stage 1, $R^2 = 0.97$ for stage 2), and are systematically lower for stage 2 than those for stage 1 (Figure 7a).

There is no thickening of the continental plate associated with shortening of the upper plate, even for the model with the largest horizontal finite shortening (-4.6% for model with $v_f = -10$ cm/yr) (Figure 7b). Models with fast mantle flow velocities show no vertical deformation, while the finite amount of thinning after 50 Myr is only of -3.5% for the model with the slowest mantle flow (Figure 5b). This is in agreement with previous models with forward mantle flow, where finite thickening of the upper plate was only observed for a finite amount of shortening larger than -5% (Figure 5).

# 4. Discussion

## 4.1. Evaluating Mantle-Plate Coupling

Our models include viscous coupling between the sub-lithospheric mantle and the plates. This coupling occurs at the base of each plate along the horizontal portions of both plates, but also along the subducted part of the lower plate. The addition of a background mantle flow modifies this coupling, which has strong implications on the kinematics and dynamics of the subduction system.

At the base of the lithosphere, the total shear force induced by the mantle is directly proportional to its viscosity and to the velocity gradient. Given that in our models the mantle has a Newtonian rheology and thus a strain-rate independent viscosity, the traction exerted at the base of the plate should increase linearly with increasing mantle velocity, considering that the distance over which the gradient applies does not change. To explore the sensitivity of the free plate motion to imposed mantle flow velocity, we define $v_{fp}/v_f$, the ratio between the free plate and the mantle flow velocities. In the absence of subduction, in our linearly viscous system, $v_{fp}/v_f$ is found to be 1 (complete plate-mantle coupling, experiment not shown in this manuscript). When subduction and background flow are active, $v_{fp}/v_f$ drops below 1 with values varying between 0.05 and 0.77 depending on the subduction stage and on the flow direction (Figure 8). For models with forward flow, one could assume that $v_{fp}/v_f$ would reach values larger than 1, as the subducting plate is not only dragged toward the trench by the imposed mantle flow,
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Ratio $v_{fp}/v_f$ between free plate velocity and imposed mantle flow velocity as a function of mantle flow velocity for the two stages of subduction (stage 1 in blue, stage 2 in orange). For models with a single plate and without subduction, $v_{fp}/v_f$ is found to be 1.</figcaption></figure>

but also by the slab pull. However, $v_{fp}/v_f$ never exceeds 0.8 and is in general of ~0.5 during stage 1 and of ~0.2 during stage 2 (Figure 8). For models with backward flow, $v_{fp}/v_f$ values are even lower, reaching minimum values of 0.05 during stage 2. These relatively low values may be to some extent related to the choice of a Newtonian rheology for the upper mantle. Indeed, models including upper mantle composite rheology showed that dislocation creep can produce weakening in the sub-plate sheared layer relative to the underlying ambient mantle thus facilitating trenchward subducting plate motion (Cerpa et al., 2022; Holt & Becker, 2016), which could result in larger values of $v_{fp}/v_f$.

Another important observation is that $v_{fp}/v_f$ is almost constant within each subduction configuration (stage, flow direction) despite varying mantle flow velocities (Figure 8). This implies that mantle drag increases almost linearly with mantle flow velocity. While having composite rheology may imply larger values of $v_{fp}/v_f$, it is difficult to anticipate if the relation between mantle flow velocity and mantle drag would depart from being linear. This would require further investigation that is beyond the scope of this study. Further-

more, the nearly constant $v_{fp}/v_f$ ratio suggests that the variations in the slab to mantle viscosity ratio resulting from variable laboratory experimental conditions (Table 1) does not significantly influence the plate-mantle coupling.

Our initial model, in agreement with observations made in previous studies (Guillaume et al., 2021), shows that the trench retreats if the mantle beneath the slab moves around the edges of the subduction zone toward the area above the slab (toroidal flow). When backward mantle flow is imposed, this slab-induced toroidal flow is reduced (see Guillaume et al., 2021), preventing fast slab rollback and trench retreat, and in turn preventing fast upper plate displacement. Instead, forward mantle flow is in the same direction as subduction and accompanies the toroidal flow around the slab. The faster plate velocity in this configuration therefore results from the combination of mantle drag at the base of the plate, but also from the facilitated direct traction of the slab that is pulling the subducting plate toward the trench.

## 4.2. Stage 1 (Free Slab Sinking) Versus Stage 2 (Resisted Slab Sinking) Subductions

Our results indicate that most model parameters evolve differently before and after the slab reaches the impermeable 660-km discontinuity of the experimental setup (Figures 5–8). One may keep in mind that, over the duration of the models (50 Myr), most of the models with mantle flow remain in a transient stage where the slab never settles flat at the bottom of the tank.

We observe that the free plate is moving toward the trench at systematically lower rates after the slab reaches the 660-km discontinuity (Figure 6c). However, displacement of the free plate remains correlated with the imposed mantle flow during this second stage (Figure 8). It indicates that horizontal displacement of the free plate is mostly driven by the mantle drag at the base of the lithosphere, which increases with increasing mantle flow velocity.

Deformation of the overriding plate also generally records significantly lower rates of horizontal and vertical deformation after the slab has reached the 660-km discontinuity for models with forward mantle flow (Figure 7). For backward mantle flow, deformation rates are similar within measurement uncertainties between the two stages. In both cases, deformation rates within the overriding plate are well correlated with the free plate velocity ($R^2$ = 0.95 and 0.76 for forward and backward mantle flow, respectively) (Figure 9a), indicating that within the chosen boundary conditions the horizontal displacement of the free plate - and thus mantle drag - controls the deformation rate of the overriding plate during the entire duration of the experiments. For models with backward mantle flow, these results are in agreement with those by Guillaume et al. (2018) and Cerpa et al. (2018) that both showed through analog and numerical models that upper plate deformation is controlled by the velocity of the overriding plate. However, in these previous studies, the motion of the plate was imposed by a piston pushing the upper plate, in contrast with our study in which motion of the overriding plate only results from mantle traction exerted by the subducting slab and the imposed mantle flow.

Our models show that the interaction of the slab with the 660-km discontinuity has profound effects on the kinematics and tectonics of the modeled subduction system. This deep interaction indeed changes the force
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Ratio $v_{fp}/v_f$ between free plate velocity and imposed mantle flow velocity as a function of mantle flow velocity for the two stages of subduction (stage 1 in blue, stage 2 in orange). For models with a single plate and without subduction, $v_{fp}/v_f$ is found to be 1.</figcaption></figure>

but also by the slab pull. However, $v_{fp}/v_f$ never exceeds 0.8 and is in general of ~0.5 during stage 1 and of ~0.2 during stage 2 (Figure 8). For models with backward flow, $v_{fp}/v_f$ values are even lower, reaching minimum values of 0.05 during stage 2. These relatively low values may be to some extent related to the choice of a Newtonian rheology for the upper mantle. Indeed, models including upper mantle composite rheology showed that dislocation creep can produce weakening in the sub-plate sheared layer relative to the underlying ambient mantle thus facilitating trenchward subducting plate motion (Cerpa et al., 2022; Holt & Becker, 2016), which could result in larger values of $v_{fp}/v_f$.

Another important observation is that $v_{fp}/v_f$ is almost constant within each subduction configuration (stage, flow direction) despite varying mantle flow velocities (Figure 8). This implies that mantle drag increases almost linearly with mantle flow velocity. While having composite rheology may imply larger values of $v_{fp}/v_f$, it is difficult to anticipate if the relation between mantle flow velocity and mantle drag would depart from being linear. This would require further investigation that is beyond the scope of this study. Further-

more, the nearly constant $v_{fp}/v_f$ ratio suggests that the variations in the slab to mantle viscosity ratio resulting from variable laboratory experimental conditions (Table 1) does not significantly influence the plate-mantle coupling.

Our initial model, in agreement with observations made in previous studies (Guillaume et al., 2021), shows that the trench retreats if the mantle beneath the slab moves around the edges of the subduction zone toward the area above the slab (toroidal flow). When backward mantle flow is imposed, this slab-induced toroidal flow is reduced (see Guillaume et al., 2021), preventing fast slab rollback and trench retreat, and in turn preventing fast upper plate displacement. Instead, forward mantle flow is in the same direction as subduction and accompanies the toroidal flow around the slab. The faster plate velocity in this configuration therefore results from the combination of mantle drag at the base of the plate, but also from the facilitated direct traction of the slab that is pulling the subducting plate toward the trench.

## 4.2. Stage 1 (Free Slab Sinking) Versus Stage 2 (Resisted Slab Sinking) Subductions

Our results indicate that most model parameters evolve differently before and after the slab reaches the impermeable 660-km discontinuity of the experimental setup (Figures 5–8). One may keep in mind that, over the duration of the models (50 Myr), most of the models with mantle flow remain in a transient stage where the slab never settles flat at the bottom of the tank.

We observe that the free plate is moving toward the trench at systematically lower rates after the slab reaches the 660-km discontinuity (Figure 6c). However, displacement of the free plate remains correlated with the imposed mantle flow during this second stage (Figure 8). It indicates that horizontal displacement of the free plate is mostly driven by the mantle drag at the base of the lithosphere, which increases with increasing mantle flow velocity.

Deformation of the overriding plate also generally records significantly lower rates of horizontal and vertical deformation after the slab has reached the 660-km discontinuity for models with forward mantle flow (Figure 7). For backward mantle flow, deformation rates are similar within measurement uncertainties between the two stages. In both cases, deformation rates within the overriding plate are well correlated with the free plate velocity ($R^2$ = 0.95 and 0.76 for forward and backward mantle flow, respectively) (Figure 9a), indicating that within the chosen boundary conditions the horizontal displacement of the free plate - and thus mantle drag - controls the deformation rate of the overriding plate during the entire duration of the experiments. For models with backward mantle flow, these results are in agreement with those by Guillaume et al. (2018) and Cerpa et al. (2018) that both showed through analog and numerical models that upper plate deformation is controlled by the velocity of the overriding plate. However, in these previous studies, the motion of the plate was imposed by a piston pushing the upper plate, in contrast with our study in which motion of the overriding plate only results from mantle traction exerted by the subducting slab and the imposed mantle flow.

Our models show that the interaction of the slab with the 660-km discontinuity has profound effects on the kinematics and tectonics of the modeled subduction system. This deep interaction indeed changes the force
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Deformation of the upper plate and free plate velocity. (a) Horizontal deformation rate for models with mantle flow as a function of the free plate velocity during stage 1 (blue symbols) and stage 2 (orange symbols) of subduction. S is the slope of the linear regression. (b) Finite amount of upper plate trench-orthogonal deformation ($ε_h$) as a function of the finite amount of upper plate vertical deformation ($ε_v$), after 50 Myr for models with mantle flow. Upper plate thickening ($ε_v$ > 0) only develops when upper plate trench-orthogonal shortening ($ε_h$) overpasses ~4% (gray dashed lines).</figcaption></figure>

equilibrium in the system by increasing resisting forces (basal resistance force and basal frictional force at the 660-km discontinuity, Figure 2), which in turn slows down subduction velocity, plate velocity and upper plate deformation for a given mantle flow. It is worth noticing that in our experimental set-up the 660-km discontinuity is an impermeable boundary, while in nature it rather corresponds to a jump in viscosity that does not necessarily impede the slab to penetrate into the lower mantle (e.g., Forte & Mitrovica, 1996; Hager, 1984; Hager & Richards, 1989; King & Hager, 1994). As such, the effect of the basal resistance force and the modification of the mantle flow and associated viscous resistive forces are possibly overestimated in our models. It means that the decrease observed in our models for most measured parameters between stage 1 and stage 2 should be seen as an upper bound. Previous studies even suggest that the penetration of the slab into the lower mantle could result in a change of the size of the convective cells that ultimately leads to a shift from “slab pull” type mountain belts to “slab suction” type (Faccenna et al., 2013), characterized by a change of deformational regime with larger strain rates in the overriding plate. Such a model has been used in particular to explain the initiation of Andean orogeny following slab penetration into the lower mantle (Chen et al., 2019; Faccenna et al., 2017; Schellart, 2017). While we cannot directly compare these results with our models that do not include a lower mantle, we notice that for a given free plate velocity, the mean horizontal deformation rates are similar or slightly larger for stage 2 (Figure 9a). In addition, it has been proposed that penetration of the slab into the lower mantle results in increased upper mantle flow (Faccenna et al., 2017; Schellart, 2017), which according to our set of models should indeed result in increased strain rates.

### 4.3. Impact of Mantle Flow on the Deformation of the Continental Upper Plate

Habel et al. (2023) recently showed that in models with an unidirectional background mantle flow and a fixed opposite plate exerting a reaction force, upper plate deformation was mainly accommodated by horizontal trench-orthogonal shortening regardless of the direction of mantle flow with respect to subduction. However, they only tested a mantle flow velocity of 10 cm/yr, either in a forward or backward direction, and as such were not able to identify when the deformation regime shifted from trench-orthogonal stretching (initial model with no mantle flow) to trench-orthogonal shortening, and if the deformation rate was proportional to the applied mantle flow velocity. With our systematic approach, we have shown that the stresses that apply on the upper plate in the presence of mantle flow and a fixed opposite plate almost always result in trench-orthogonal shortening. It is only for the models with the slowest mantle flow ($v_f = \pm 1.66$ cm/yr) that no deformation is observed, at least in one of the two subduction stages (Figure 7a). In the absence of mantle flow, the trench is mobile and retreats toward the subducting plate at a rate of $\pm 0.4$ cm/yr, accompanied by trench-orthogonal stretching. For forward mantle flow, the slab and the lower plate are displaced toward the upper plate, resulting in trench advance, with a velocity that increases with increasing mantle flow. Significant trench-orthogonal shortening is only recorded when mantle flow velocity exceeds the intrinsic trench retreat velocity. For backward mantle flow, the trench retreats toward the subducting plate. Again, this is only when the imposed mantle flow velocity
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Deformation of the upper plate and free plate velocity. (a) Horizontal deformation rate for models with mantle flow as a function of the free plate velocity during stage 1 (blue symbols) and stage 2 (orange symbols) of subduction. S is the slope of the linear regression. (b) Finite amount of upper plate trench-orthogonal deformation ($ε_h$) as a function of the finite amount of upper plate vertical deformation ($ε_v$), after 50 Myr for models with mantle flow. Upper plate thickening ($ε_v$ > 0) only develops when upper plate trench-orthogonal shortening ($ε_h$) overpasses ~4% (gray dashed lines).</figcaption></figure>

equilibrium in the system by increasing resisting forces (basal resistance force and basal frictional force at the 660-km discontinuity, Figure 2), which in turn slows down subduction velocity, plate velocity and upper plate deformation for a given mantle flow. It is worth noticing that in our experimental set-up the 660-km discontinuity is an impermeable boundary, while in nature it rather corresponds to a jump in viscosity that does not necessarily impede the slab to penetrate into the lower mantle (e.g., Forte & Mitrovica, 1996; Hager, 1984; Hager & Richards, 1989; King & Hager, 1994). As such, the effect of the basal resistance force and the modification of the mantle flow and associated viscous resistive forces are possibly overestimated in our models. It means that the decrease observed in our models for most measured parameters between stage 1 and stage 2 should be seen as an upper bound. Previous studies even suggest that the penetration of the slab into the lower mantle could result in a change of the size of the convective cells that ultimately leads to a shift from “slab pull” type mountain belts to “slab suction” type (Faccenna et al., 2013), characterized by a change of deformational regime with larger strain rates in the overriding plate. Such a model has been used in particular to explain the initiation of Andean orogeny following slab penetration into the lower mantle (Chen et al., 2019; Faccenna et al., 2017; Schellart, 2017). While we cannot directly compare these results with our models that do not include a lower mantle, we notice that for a given free plate velocity, the mean horizontal deformation rates are similar or slightly larger for stage 2 (Figure 9a). In addition, it has been proposed that penetration of the slab into the lower mantle results in increased upper mantle flow (Faccenna et al., 2017; Schellart, 2017), which according to our set of models should indeed result in increased strain rates.

### 4.3. Impact of Mantle Flow on the Deformation of the Continental Upper Plate

Habel et al. (2023) recently showed that in models with an unidirectional background mantle flow and a fixed opposite plate exerting a reaction force, upper plate deformation was mainly accommodated by horizontal trench-orthogonal shortening regardless of the direction of mantle flow with respect to subduction. However, they only tested a mantle flow velocity of 10 cm/yr, either in a forward or backward direction, and as such were not able to identify when the deformation regime shifted from trench-orthogonal stretching (initial model with no mantle flow) to trench-orthogonal shortening, and if the deformation rate was proportional to the applied mantle flow velocity. With our systematic approach, we have shown that the stresses that apply on the upper plate in the presence of mantle flow and a fixed opposite plate almost always result in trench-orthogonal shortening. It is only for the models with the slowest mantle flow ($v_f = \pm 1.66$ cm/yr) that no deformation is observed, at least in one of the two subduction stages (Figure 7a). In the absence of mantle flow, the trench is mobile and retreats toward the subducting plate at a rate of $\pm 0.4$ cm/yr, accompanied by trench-orthogonal stretching. For forward mantle flow, the slab and the lower plate are displaced toward the upper plate, resulting in trench advance, with a velocity that increases with increasing mantle flow. Significant trench-orthogonal shortening is only recorded when mantle flow velocity exceeds the intrinsic trench retreat velocity. For backward mantle flow, the trench retreats toward the subducting plate. Again, this is only when the imposed mantle flow velocity
exceeds the intrinsic trench retreat velocity that the deformation regime in the upper plate shifts from stretching to shortening, as displacement of the upper plate becomes forced by the traction exerted by the mantle. This is consistent with previous analog models in which the shift was observed when the upper plate velocity exceeded the intrinsic trench retreat velocity (Guillaume et al., 2018). As imposed mantle flow and free plate velocity are linearly correlated in our models, the shift from extensional deformation regime to shortening, which occurs at flow velocities of ~1–2 cm/yr, corresponds to free plate velocities of ~0.4 cm/yr (Figure 9a).

Our series of experiments also allows us to determine how upper plate finite strain and strain rate evolve with increas-
ing absolute velocity of the imposed mantle flow in the box (lower mantle) reference frame. We show that whatever
the direction of the flow, this relation is linear with good correlation coefficients (Figures 7 and 9). This demonstrates
that mantle drag exerts a strong control on plate deformation, as mantle drag scales with mantle flow in our linearly
viscous subduction system. Instead, the slab pull force only marginally varies in between experiments because of
minor variations of slab geometries. Slab pull force control on upper plate deformation is therefore only important
for small mantle flow velocities (<1.66 cm/yr), that is, for subduction systems isolated from large-scale mantle flow.

Deformation of the upper plate is not plane strain, as trench-orthogonal deformation is not entirely compensated by vertical deformation. Indeed, we also observe trench-parallel stretching, which in some cases can lead to upper plate thinning. By correlating finite upper plate trench-orthogonal shortening and thickening after 50 Myr, it appears that continental plate thickening only occurs for amounts of bulk trench-orthogonal shortening above ~4% (Figure 9b). This result can be explained by the presence of unconfined lateral boundaries for the upper plate in our models, which allow lateral viscous flow. If the trench-orthogonal shortening is not sufficient to compensate for the trench-parallel stretching, then no thickening is observed. We do not necessarily expect the natural systems to record thickening only after 4% of shortening because in nature there are generally no free lateral edges for the upper plate. However, deformation associated with plate convergence may also be accommodated to a certain extent by trench-parallel creep of the continental plate. In this case, the combination of trench-orthogonal shortening and trench-parallel stretching may lead to the absence of relief and to the partitioning of deformation.

4.4. Comparison With Previous Models and Implications for Natural Subduction Zones

4.4.1. Slab Geometry

During the first stage of subduction, the slab vertical geometry appears to not vary much between models, despite varying mantle flow velocities (Figure 4). This indicates that the mantle flow is compensated by the trench motion, that is, the slab kinematics and geometry remains almost identical in the upper mantle reference frame. Instead, in the box reference frame, the trench migrates at rates increasing with increasing mantle flow (Figure 6b). Correlation between free plate velocity, slab horizontal velocity at upper mantle middepth and slab tip horizontal velocity (Figure 6) demonstrates that under the applied mantle flow, the slab laterally migrates without significant internal deformation. This is in general agreement with previous results by Guillaume et al. (2021). The main difference is that the most vertical slabs in the models by Guillaume et al. (2021) were obtained for forward mantle flow, while the opposite can be observed here. This can readily be explained by the previous different adopted set-up that did not include an upper plate and considered the velocity of the lower plate as null. Both factors favor trench advance and lower plate thinning for forward mantle flow, which in turn promote steepening of the slab. Chertova et al. (2018) also explored the role of horizontal mantle flow on slab geometry with 3D thermo-mechanical numerical models including a composite nonlinear rheology. They showed that for models with boundary conditions close to ours, that is, trench-orthogonal mantle flow, the slab geometry remained almost similar to the model without flow, highlighting the small impact of mantle flow in this case on slab deformation. Despite different boundary conditions, the subducting plate velocity being prescribed at 1.5 cm/yr in their models, the limited impact of mantle flow on slab geometry suggests that the use of a linear or composite nonlinear rheology does not affect first-order conclusions. Interestingly, in models with trench-oblique mantle flow, larger differences in slab geometry have been reported (Chertova et al., 2018), as a response to more complex interactions between imposed and slab-induced mantle flows. The effect of trench-oblique mantle flow, although possible on Earth especially in areas with complex geometries (e.g., Mediterranean, Rappisi et al., 2022), was not addressed in our study. We can nevertheless anticipate that it may lead to more variability of along-strike slab geometry and in turn to more complex patterns of deformation of the upper plate.

We also show that during the second stage of subduction, there is no major variation in slab geometry with
increasing mantle flow velocity for backward flow (Figure 4, left panels). For forward mantle flow, the slab
exceeds the intrinsic trench retreat velocity that the deformation regime in the upper plate shifts from stretching to shortening, as displacement of the upper plate becomes forced by the traction exerted by the mantle. This is consistent with previous analog models in which the shift was observed when the upper plate velocity exceeded the intrinsic trench retreat velocity (Guillaume et al., 2018). As imposed mantle flow and free plate velocity are linearly correlated in our models, the shift from extensional deformation regime to shortening, which occurs at flow velocities of ~1–2 cm/yr, corresponds to free plate velocities of ~0.4 cm/yr (Figure 9a).

Our series of experiments also allows us to determine how upper plate finite strain and strain rate evolve with increas-
ing absolute velocity of the imposed mantle flow in the box (lower mantle) reference frame. We show that whatever
the direction of the flow, this relation is linear with good correlation coefficients (Figures 7 and 9). This demonstrates
that mantle drag exerts a strong control on plate deformation, as mantle drag scales with mantle flow in our linearly
viscous subduction system. Instead, the slab pull force only marginally varies in between experiments because of
minor variations of slab geometries. Slab pull force control on upper plate deformation is therefore only important
for small mantle flow velocities (<1.66 cm/yr), that is, for subduction systems isolated from large-scale mantle flow.

Deformation of the upper plate is not plane strain, as trench-orthogonal deformation is not entirely compensated by vertical deformation. Indeed, we also observe trench-parallel stretching, which in some cases can lead to upper plate thinning. By correlating finite upper plate trench-orthogonal shortening and thickening after 50 Myr, it appears that continental plate thickening only occurs for amounts of bulk trench-orthogonal shortening above ~4% (Figure 9b). This result can be explained by the presence of unconfined lateral boundaries for the upper plate in our models, which allow lateral viscous flow. If the trench-orthogonal shortening is not sufficient to compensate for the trench-parallel stretching, then no thickening is observed. We do not necessarily expect the natural systems to record thickening only after 4% of shortening because in nature there are generally no free lateral edges for the upper plate. However, deformation associated with plate convergence may also be accommodated to a certain extent by trench-parallel creep of the continental plate. In this case, the combination of trench-orthogonal shortening and trench-parallel stretching may lead to the absence of relief and to the partitioning of deformation.

4.4. Comparison With Previous Models and Implications for Natural Subduction Zones

4.4.1. Slab Geometry

During the first stage of subduction, the slab vertical geometry appears to not vary much between models, despite varying mantle flow velocities (Figure 4). This indicates that the mantle flow is compensated by the trench motion, that is, the slab kinematics and geometry remains almost identical in the upper mantle reference frame. Instead, in the box reference frame, the trench migrates at rates increasing with increasing mantle flow (Figure 6b). Correlation between free plate velocity, slab horizontal velocity at upper mantle middepth and slab tip horizontal velocity (Figure 6) demonstrates that under the applied mantle flow, the slab laterally migrates without significant internal deformation. This is in general agreement with previous results by Guillaume et al. (2021). The main difference is that the most vertical slabs in the models by Guillaume et al. (2021) were obtained for forward mantle flow, while the opposite can be observed here. This can readily be explained by the previous different adopted set-up that did not include an upper plate and considered the velocity of the lower plate as null. Both factors favor trench advance and lower plate thinning for forward mantle flow, which in turn promote steepening of the slab. Chertova et al. (2018) also explored the role of horizontal mantle flow on slab geometry with 3D thermo-mechanical numerical models including a composite nonlinear rheology. They showed that for models with boundary conditions close to ours, that is, trench-orthogonal mantle flow, the slab geometry remained almost similar to the model without flow, highlighting the small impact of mantle flow in this case on slab deformation. Despite different boundary conditions, the subducting plate velocity being prescribed at 1.5 cm/yr in their models, the limited impact of mantle flow on slab geometry suggests that the use of a linear or composite nonlinear rheology does not affect first-order conclusions. Interestingly, in models with trench-oblique mantle flow, larger differences in slab geometry have been reported (Chertova et al., 2018), as a response to more complex interactions between imposed and slab-induced mantle flows. The effect of trench-oblique mantle flow, although possible on Earth especially in areas with complex geometries (e.g., Mediterranean, Rappisi et al., 2022), was not addressed in our study. We can nevertheless anticipate that it may lead to more variability of along-strike slab geometry and in turn to more complex patterns of deformation of the upper plate.

We also show that during the second stage of subduction, there is no major variation in slab geometry with
increasing mantle flow velocity for backward flow (Figure 4, left panels). For forward mantle flow, the slab
tends to flatten for mantle flow velocity ≥5 cm/yr, and even adopts a roll-over shape when mantle flow remains at velocities of 10 cm/yr for several My (model SA13) (Figure 4, right panels). This is made possible because the slab arrives almost vertically at the 660-km discontinuity and the upper plate is shortened, favoring trench advance and slab roll-over. Velocities for horizontal upper mantle flow exceeding 5 cm/yr have been proposed to result for instance for the South Atlantic region from pressure-induced upper mantle flow related to plume activity, explaining changes in South America directional plate motion in the Paleocene and in the Late Eocene/ Oligocene (Colli et al., 2014; Stotz et al., 2023). Slab geometry and upper plate tectonic regime do not evolve together in our models. Indeed, during stage 1, slab dip is almost similar for all models whatever the direction and magnitude of the mantle flow (Figure 4) while deformation rates increase with mantle flow (Figure 7). In addition, after the arrival of the slab at the 660-km discontinuity the slab tends to flatten for models with forward flow and to steepen for models with backward flow while in both cases the shortening rate decreases between stage 1 and stage 2 (Figure 7). These observations suggest that the correlation between slab geometry and upper plate deformation is weak in our models, hence that slab geometry only marginally affects upper plate tectonics within the chosen boundary conditions. This is in agreement with observations of current subduction zones that show no evident correlation between shallow or deep slab dip and the trench-normal component of the upper plate deformation rate (R² values <0.06, after SUBMAP database: https://submap.fr) (Lallemand & Heuret, 2017). Nevertheless, previous studies proposed that changes of slab geometry can locally impact upper plate tectonics regime. In particular, it has been suggested that flat slabs can generate a stronger coupling between the subducting and overriding plates and episodes of intense shortening (e.g., Capitanio et al., 2011; Martinod et al., 2010, 2020; Molnar & Atwater, 1978), indicating that other local mechanisms (subduction of oceanic plateaus, aseismic ridges) could be active in these areas. The apparent absence of correlation between slab geometry and upper plate tectonics in our models and for current natural subduction zones may indicate that subduction systems are three-dimensional and as such that slab geometry and upper plate deformation are not affected in the same way by lateral mantle flow.

### 4.4.2. Upper Plate Strain Versus Mean Forces at the Trench: A Comparison to Natural Subductions

Our models show that the upper plate deformation rate linearly correlates with the velocity of the free plate (Figure 9a), which itself depends on the imposed mantle flow velocity for the chosen linear rheology (Figures 6c, 8, and 10).

Husson (2012) proposed that upper plate strain varies as a function of the mean trench normal integrated mantle drag force from both the upper and lower plates (Σmd/2) for present-day subduction zones. This earlier study suggests that areas undergoing significant shortening (upper plate strain categories 2 and 3, after Heuret and Lallemand (2005)) mostly correspond to values of Σmd/2 larger than 2.5 × 10¹² N/m. We cannot directly compute (Σmd/2) in our models to compare with these findings, but a rough estimate can be established. Indeed, Husson (2012) indicate that the mean force applying on the trench and causing upper plate deformation writes:

$$ \Sigma F/2 = (F_{up} + F_{lp})/2 = (md_{up} + rp_{up} + md_{lp} + sp^* + rp_{lp})/2 = h \eta_{up} \dot{\epsilon} \quad (2) $$

where subscripts up and lp stand for upper plate and lower plate, respectively, md is the drag from the underlying flowing mantle, rp the ridge push, sp* the effective slab pull (slab pull minus viscous dissipation in the slab and mantle), h the upper plate thickness, η_up its viscosity, and ξ its deformation rate.

The initial model (SA06) for which there is no ridge push and no imposed mantle flow constrains the effective slab pull sp* for the other models. The average value for ξ in this model is ~1.27 × 10⁻¹⁷ s⁻¹ (Figure 5), which gives ΣF/2 for this model of ~0.8 × 10¹¹ N/m and therefore implies that sp* is low in the other models and thus negligible. For the other models, the imposed unidirectional flow implies that mantle drag at the base of the free plate, on one side of the subduction, largely dominates the force balance, together with the reaction force exerted by the backwall on the fixed plate, which can be assimilated to a force pushing this plate toward the trench (e.g., ridge push) (Habel et al., 2023) (Figure 3). Equation 2 can then be simplified to:

$$ \Sigma F/2 = (rp_{up} + md_{lp})/2 \quad (3) $$

for forward mantle flow, and:

$$ \Sigma F/2 = (rp_{lp} + md_{up})/2 \quad (4) $$
tends to flatten for mantle flow velocity ≥5 cm/yr, and even adopts a roll-over shape when mantle flow remains at velocities of 10 cm/yr for several My (model SA13) (Figure 4, right panels). This is made possible because the slab arrives almost vertically at the 660-km discontinuity and the upper plate is shortened, favoring trench advance and slab roll-over. Velocities for horizontal upper mantle flow exceeding 5 cm/yr have been proposed to result for instance for the South Atlantic region from pressure-induced upper mantle flow related to plume activity, explaining changes in South America directional plate motion in the Paleocene and in the Late Eocene/ Oligocene (Colli et al., 2014; Stotz et al., 2023). Slab geometry and upper plate tectonic regime do not evolve together in our models. Indeed, during stage 1, slab dip is almost similar for all models whatever the direction and magnitude of the mantle flow (Figure 4) while deformation rates increase with mantle flow (Figure 7). In addition, after the arrival of the slab at the 660-km discontinuity the slab tends to flatten for models with forward flow and to steepen for models with backward flow while in both cases the shortening rate decreases between stage 1 and stage 2 (Figure 7). These observations suggest that the correlation between slab geometry and upper plate deformation is weak in our models, hence that slab geometry only marginally affects upper plate tectonics within the chosen boundary conditions. This is in agreement with observations of current subduction zones that show no evident correlation between shallow or deep slab dip and the trench-normal component of the upper plate deformation rate (R² values <0.06, after SUBMAP database: https://submap.fr) (Lallemand & Heuret, 2017). Nevertheless, previous studies proposed that changes of slab geometry can locally impact upper plate tectonics regime. In particular, it has been suggested that flat slabs can generate a stronger coupling between the subducting and overriding plates and episodes of intense shortening (e.g., Capitanio et al., 2011; Martinod et al., 2010, 2020; Molnar & Atwater, 1978), indicating that other local mechanisms (subduction of oceanic plateaus, aseismic ridges) could be active in these areas. The apparent absence of correlation between slab geometry and upper plate tectonics in our models and for current natural subduction zones may indicate that subduction systems are three-dimensional and as such that slab geometry and upper plate deformation are not affected in the same way by lateral mantle flow.

### 4.4.2. Upper Plate Strain Versus Mean Forces at the Trench: A Comparison to Natural Subductions

Our models show that the upper plate deformation rate linearly correlates with the velocity of the free plate (Figure 9a), which itself depends on the imposed mantle flow velocity for the chosen linear rheology (Figures 6c, 8, and 10).

Husson (2012) proposed that upper plate strain varies as a function of the mean trench normal integrated mantle drag force from both the upper and lower plates (Σmd/2) for present-day subduction zones. This earlier study suggests that areas undergoing significant shortening (upper plate strain categories 2 and 3, after Heuret and Lallemand (2005)) mostly correspond to values of Σmd/2 larger than 2.5 × 10¹² N/m. We cannot directly compute (Σmd/2) in our models to compare with these findings, but a rough estimate can be established. Indeed, Husson (2012) indicate that the mean force applying on the trench and causing upper plate deformation writes:

$$ \Sigma F/2 = (F_{up} + F_{lp})/2 = (md_{up} + rp_{up} + md_{lp} + sp^* + rp_{lp})/2 = h \eta_{up} \dot{\epsilon} \quad (2) $$

where subscripts up and lp stand for upper plate and lower plate, respectively, md is the drag from the underlying flowing mantle, rp the ridge push, sp* the effective slab pull (slab pull minus viscous dissipation in the slab and mantle), h the upper plate thickness, η_up its viscosity, and ξ its deformation rate.

The initial model (SA06) for which there is no ridge push and no imposed mantle flow constrains the effective slab pull sp* for the other models. The average value for ξ in this model is ~1.27 × 10⁻¹⁷ s⁻¹ (Figure 5), which gives ΣF/2 for this model of ~0.8 × 10¹¹ N/m and therefore implies that sp* is low in the other models and thus negligible. For the other models, the imposed unidirectional flow implies that mantle drag at the base of the free plate, on one side of the subduction, largely dominates the force balance, together with the reaction force exerted by the backwall on the fixed plate, which can be assimilated to a force pushing this plate toward the trench (e.g., ridge push) (Habel et al., 2023) (Figure 3). Equation 2 can then be simplified to:

$$ \Sigma F/2 = (rp_{up} + md_{lp})/2 \quad (3) $$

for forward mantle flow, and:

$$ \Sigma F/2 = (rp_{lp} + md_{up})/2 \quad (4) $$
for backward mantle flow.

$\Sigma F/2$ can be computed from the deformation of the plate (Equation 2), which, considering the reported mean strain rate values ($\epsilon_{\text{dot}}$) (Figure 7), gives values ranging between 0 and ~0.3 $\times 10^{12}$ N/m for backward flow (Figure 10). These values are around one order of magnitude below the mean trench normal integrated mantle drag force proposed by Husson (2012) at $2.5 \times 10^{12}$ N/m to shift from extensional/neutral regime to compressional regime. For forward mantle flow, the largest $\Sigma F/2$ values are found during stage 1 for the model with the fastest flow with a value of $3.9 \times 10^{12}$ N/m. However, $\Sigma F/2$ decreases rapidly with decreasing velocity of mantle flow and drops below $1.5 \times 10^{12}$ N/m for stage 2 of all models and for stage 1 of models with mantle flow velocities $\le 5$ cm/yr (Figure 10). The direct comparison of our models with the results by Husson (2012) on present-day subduction zones is not straightforward because he only considers mantle drag in his calculations, while ours also include the reaction force exerted by the back wall. However, the fact that we only achieved forces higher than $2.5 \times 10^{12}$ N/m for the fastest model ($v_f = 10$ cm/yr) may indicate that our linearly viscous model does not encompass all the rheological complexity found in nature.

Our modeling results can nevertheless give us some insights on the effect of mantle drag when compared to present-day subduction zones (Figures 1 and 10). For the Tonga-Kermadec subduction zone, present-day models of mantle flow suggest that the mean trench normal integrated mantle drag force is small ($\Sigma \text{md}/2$ on the order of $0 \text{ to } -1 \times 10^{12}$ N/m) and directed toward the west (i.e., toward the upper plate) beneath both the lower and upper

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Synthetic diagrams summarizing the evolution of various model parameters, during stage 2, as a function of the direction of the imposed mantle flow with respect to the subduction, as well as of its velocity. Left: backward mantle flow; right: forward flow. (a) Geometry of the slab after 50 Myr in our models, with mantle flow direction. (b) Horizontal trench-orthogonal deformation rate. (c) Vertical deformation rate. (d) Mean force applied on the trench ($\Sigma F/2$). (e) The absolute value of free plate velocity. (f) Trench velocity. (g) Subduction velocity.</figcaption></figure>
for backward mantle flow.

$\Sigma F/2$ can be computed from the deformation of the plate (Equation 2), which, considering the reported mean strain rate values ($\text{edot}$) (Figure 7), gives values ranging between 0 and ~0.3 $\times 10^{12}$ N/m for backward flow (Figure 10). These values are around one order of magnitude below the mean trench normal integrated mantle drag force proposed by Husson (2012) at $2.5 \times 10^{12}$ N/m to shift from extensional/neutral regime to compressional regime. For forward mantle flow, the largest $\Sigma F/2$ values are found during stage 1 for the model with the fastest flow with a value of $3.9 \times 10^{12}$ N/m. However, $\Sigma F/2$ decreases rapidly with decreasing velocity of mantle flow and drops below $1.5 \times 10^{12}$ N/m for stage 2 of all models and for stage 1 of models with mantle flow velocities $\le 5$ cm/yr (Figure 10). The direct comparison of our models with the results by Husson (2012) on present-day subduction zones is not straightforward because he only considers mantle drag in his calculations, while ours also include the reaction force exerted by the back wall. However, the fact that we only achieved forces higher than $2.5 \times 10^{12}$ N/m for the fastest model ($v_f = 10$ cm/yr) may indicate that our linearly viscous model does not encompass all the rheological complexity found in nature.

Our modeling results can nevertheless give us some insights on the effect of mantle drag when compared to present-day subduction zones (Figures 1 and 10). For the Tonga-Kermadec subduction zone, present-day models of mantle flow suggest that the mean trench normal integrated mantle drag force is small ($\Sigma \text{md}/2$ on the order of $0 \text{ to } -1 \times 10^{12}$ N/m) and directed toward the west (i.e., toward the upper plate) beneath both the lower and upper

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Synthetic diagrams summarizing the evolution of various model parameters, during stage 2, as a function of the direction of the imposed mantle flow with respect to the subduction, as well as of its velocity. Left: backward mantle flow; right: forward flow. (a) Geometry of the slab after 50 Myr in our models, with mantle flow direction. (b) Horizontal trench-orthogonal deformation rate. (c) Vertical deformation rate. (d) Mean force applied on the trench ($\Sigma F/2$). (e) The absolute value of free plate velocity. (f) Trench velocity. (g) Subduction velocity.</figcaption></figure>
plates (Husson, 2012) (Figure 1b). Upper plate deformation is characterized by extensional deformation (upper plate strain classes -3 and -2; Heuret & Lallemand, 2005). This configuration resembles that of our models with the slowest forward mantle flow (models SA09) where ΣF/2 is close to zero (Figure 10). We show that weak mantle drag is the most favorable situation to obtain a neutral to extensional regime within the upper plate, thus confirming that the effect of mantle drag is probably limited in this region.

The Kuril subduction zone is a peculiar zone along the Western Pacific margin. Indeed, this is one of the few places where the regime of deformation varies between slight and moderate shortening (classes 1 to 2; Heuret & Lallemand, 2005), while it is generally extensional elsewhere (Figure 1). Interestingly, this is also a place where mantle drag beneath both the lower and upper plates is directed toward the W-NW (i.e., toward the upper plate), with a strong mantle drag beneath the lower plate and a weaker mantle drag beneath the upper plate, giving a mean mantle drag force Σmd/2 on the order of 1 to 3 × 10¹² N/m (Husson, 2012; Figure 1c). This situation resembles that of our models with a unidirectional forward mantle flow with velocities larger than 5 cm/yr, where we show that ΣF/2 ≥ 1 × 10¹² N/m and the upper plate accommodates deformation through trench-orthogonal shortening at rates larger than 10⁻¹⁶ s⁻¹ (Figure 7a). This suggests that the deformational regime for the Kuril subduction zone could be controlled by mantle drag at the base of the lithosphere. The Kuril subduction zone is also characterized by a slab with a flatter geometry than for other areas along the Western Pacific (e.g., Boutelier & Cruden, 2008). This is in agreement with our findings for models with moderate to large mantle flow (≥5 cm/yr) that also show a flatter slab in these conditions (Figure 4). Our model with the largest mantle drag shows a roll-over geometry that is not described by tomographic and seismic slab geometry reconstructions (e.g., Hayes et al., 2018), thus suggesting that mantle drag in this area may not be high enough to produce such geometry or would not have maintained these high values for a long enough period of time.

The South America subduction zone is the active prototype of a subduction zone exhibiting a subduction-related orogen. Upper plate strain is indeed highly compressive at present-day (upper plate strain class 2 to 3; Heuret & Lallemand, 2005), and it has been the case for several Myr as building of the Andes may have started as early as middle to late Cretaceous (e.g., Arriagada et al., 2006; Boschman, 2021; Jaillard, 2022). Models of mantle flow suggest that present-day mean mantle drag is high, largely exceeding 3 × 10¹² N/m in the Central Andes (Husson, 2012; Figure 1d), as a result of the large mantle drag toward the trench beneath both the lower and the upper plates. Our models show that irrespective of the mantle flow direction, this is the model with the largest mantle flow and thus the largest mantle drag that produces the largest upper plate shortening rates, with mean values up to 4.5 × 10⁻¹⁶ s⁻¹ for a forward mantle flow velocity of 10 cm/yr, corresponding to ΣF/2 of 3.9 × 10¹² N/m (Figure 10). We cannot directly compare our models with the South America subduction zone, as we can only impose mantle drag in one direction in our models. However, we impose a no-velocity condition to the upper plate that would otherwise move away from the trench. As a result, we introduce a reaction force exerted by the backwall onto the fixed plate. Such force can be assimilated to a ridge push force and/or to mantle drag toward the trench, which in turn promotes upper plate shortening (Habel et al., 2023). While imposing mantle drag in opposite directions is currently technically infeasible, we can nevertheless hypothesize that the combination of large mantle flow directed toward the trench from opposite directions would result in large values for the mean trench normal integrated mantle drag force, and in turn in large deformation rates, allowing the building of subduction-related orogens by shortening.

**5. Conclusions**

We performed a series of three-dimensional subduction models, in which we systematically study the effects of the direction and magnitude of an imposed unidirectional horizontal mantle flow on subduction dynamics and upper plate deformation, and in which the plate opposite to the flow direction (either upper or lower plate) is fixed to a backwall. In our experiments, mantle flow exerts a mantle drag at the base of the lithosphere that tends to displace the slab rather than to produce internal deformation. In the upper mantle reference frame, the slab kinematics and geometry remain almost identical. This results in almost similar vertical slab geometries for mantle flow velocities lower than 5 cm/yr. For larger mantle flow velocities and when considering forward flow, the slab tends to flatten and can even adopt a roll-over shape if mantle flow is maintained for several My. In our experiments, the plate velocity and the upper plate deformation rates linearly correlate with the imposed mantle flow, indicating that mantle drag exerts a key control on plate motion and associated deformation.

For the boundary conditions tested here where plates are systematically dragged toward the trench, the upper plate almost always accommodates deformation through trench-orthogonal shortening. However, shortening is
plates (Husson, 2012) (Figure 1b). Upper plate deformation is characterized by extensional deformation (upper plate strain classes −3 and −2; Heuret & Lallemand, 2005). This configuration resembles that of our models with the slowest forward mantle flow (models SA09) where ΣF/2 is close to zero (Figure 10). We show that weak mantle drag is the most favorable situation to obtain a neutral to extensional regime within the upper plate, thus confirming that the effect of mantle drag is probably limited in this region.

The Kuril subduction zone is a peculiar zone along the Western Pacific margin. Indeed, this is one of the few places where the regime of deformation varies between slight and moderate shortening (classes 1 to 2; Heuret & Lallemand, 2005), while it is generally extensional elsewhere (Figure 1). Interestingly, this is also a place where mantle drag beneath both the lower and upper plates is directed toward the W-NW (i.e., toward the upper plate), with a strong mantle drag beneath the lower plate and a weaker mantle drag beneath the upper plate, giving a mean mantle drag force Σmd/2 on the order of 1 to 3 × 10¹² N/m (Husson, 2012; Figure 1c). This situation resembles that of our models with a unidirectional forward mantle flow with velocities larger than 5 cm/yr, where we show that ΣF/2 ≥ 1 × 10¹² N/m and the upper plate accommodates deformation through trench-orthogonal shortening at rates larger than 10⁻¹⁶ s⁻¹ (Figure 7a). This suggests that the deformational regime for the Kuril subduction zone could be controlled by mantle drag at the base of the lithosphere. The Kuril subduction zone is also characterized by a slab with a flatter geometry than for other areas along the Western Pacific (e.g., Boutelier & Cruden, 2008). This is in agreement with our findings for models with moderate to large mantle flow (≥5 cm/yr) that also show a flatter slab in these conditions (Figure 4). Our model with the largest mantle drag shows a roll-over geometry that is not described by tomographic and seismic slab geometry reconstructions (e.g., Hayes et al., 2018), thus suggesting that mantle drag in this area may not be high enough to produce such geometry or would not have maintained these high values for a long enough period of time.

The South America subduction zone is the active prototype of a subduction zone exhibiting a subduction-related orogen. Upper plate strain is indeed highly compressive at present-day (upper plate strain class 2 to 3; Heuret & Lallemand, 2005), and it has been the case for several Myr as building of the Andes may have started as early as middle to late Cretaceous (e.g., Arriagada et al., 2006; Boschman, 2021; Jaillard, 2022). Models of mantle flow suggest that present-day mean mantle drag is high, largely exceeding 3 × 10¹² N/m in the Central Andes (Husson, 2012; Figure 1d), as a result of the large mantle drag toward the trench beneath both the lower and the upper plates. Our models show that irrespective of the mantle flow direction, this is the model with the largest mantle flow and thus the largest mantle drag that produces the largest upper plate shortening rates, with mean values up to 4.5 × 10⁻¹⁶ s⁻¹ for a forward mantle flow velocity of 10 cm/yr, corresponding to ΣF/2 of 3.9 × 10¹² N/m (Figure 10). We cannot directly compare our models with the South America subduction zone, as we can only impose mantle drag in one direction in our models. However, we impose a no-velocity condition to the upper plate that would otherwise move away from the trench. As a result, we introduce a reaction force exerted by the backwall onto the fixed plate. Such force can be assimilated to a ridge push force and/or to mantle drag toward the trench, which in turn promotes upper plate shortening (Habel et al., 2023). While imposing mantle drag in opposite directions is currently technically infeasible, we can nevertheless hypothesize that the combination of large mantle flow directed toward the trench from opposite directions would result in large values for the mean trench normal integrated mantle drag force, and in turn in large deformation rates, allowing the building of subduction-related orogens by shortening.

**5. Conclusions**

We performed a series of three-dimensional subduction models, in which we systematically study the effects of the direction and magnitude of an imposed unidirectional horizontal mantle flow on subduction dynamics and upper plate deformation, and in which the plate opposite to the flow direction (either upper or lower plate) is fixed to a backwall. In our experiments, mantle flow exerts a mantle drag at the base of the lithosphere that tends to displace the slab rather than to produce internal deformation. In the upper mantle reference frame, the slab kinematics and geometry remain almost identical. This results in almost similar vertical slab geometries for mantle flow velocities lower than 5 cm/yr. For larger mantle flow velocities and when considering forward flow, the slab tends to flatten and can even adopt a roll-over shape if mantle flow is maintained for several My. In our experiments, the plate velocity and the upper plate deformation rates linearly correlate with the imposed mantle flow, indicating that mantle drag exerts a key control on plate motion and associated deformation.

For the boundary conditions tested here where plates are systematically dragged toward the trench, the upper plate almost always accommodates deformation through trench-orthogonal shortening. However, shortening is
associated with thickening only for trench-orthogonal strain rates larger than $2.5 \times 10^{-17} s^{-1}$ as the upper plate also accommodates deformation through mild trench-parallel stretching. We also show that upper plate strain rates are three times larger when the lower plate is the one that is dragged toward the trench, as the mantle flow drives subduction instead of opposing it. In addition, the reaction force exerted by the backwall on the fixed upper plate, which may be taken as an analog to ridge push or opposite mantle drag, also favors trench-orthogonal shortening.

We finally show that the interaction of the slab with an impermeable boundary at 660-km depth has profound consequences on subduction dynamics, and in turn on upper plate deformation. For similar free plate velocities, slightly larger strain rates are recorded when the slab has already reached the 660-km discontinuity. Deformation regimes observed along present-day subduction zones are compatible with the findings from our models and from computations of present-day integrated drag forces at the trench. This suggests that the drag exerted by the convective mantle may be a key ingredient in controlling the upper plate deformation regime.

## Data Availability Statement

Raw data associated with this paper (movies of the top and lateral views of each experiment) are available following this link: https://doi.org/10.57932/e21b6394-6557-4ba7-ac6a-6ecce61d5c69.

## Acknowledgments

We would like to express our gratitude to the Associate Editor, Ernst Willingshofer, and the reviewers, Wim Spakman and the anonymous reviewer, for their constructive comments on this manuscript. Their feedback has greatly contributed to enhancing the quality of this paper.

## References

Arriagada, C., Cobbold, P. R., & Roperch, P. (2006). Salar de Atacama Basin: A record of compressional tectonics in the Central Andes since the mid-Cretaceous. *Tectonics*, 25(1), TC1008. https://doi.org/10.1029/2004TC001770

Becker, T. W., Kustowski, B., & Ekström, G. (2008). Radial seismic anisotropy as a constraint for upper mantle rheology. *Earth and Planetary Science Letters*, 267(1–2), 213–227. https://doi.org/10.1016/j.epsl.2007.11.038

Behr, W. M., Holt, A. F., Becker, T. W., & Faccenna, C. (2022). The effects of plate interface rheology on subduction kinematics and dynamics. *Geophysical Journal International*, 230(2), 796–812. https://doi.org/10.1093/gji/ggac075

Boschman, L. M. (2021). Andean mountain building since the late cretaceous: A paleoelevation reconstruction. *Earth-Science Reviews*, 220, 103640. https://doi.org/10.1016/j.earscirev.2021.103640

Boutelier, D. A., & Cruden, A. R. (2008). Impact of regional mantle flow on subducting plate geometry and interplate stress: Insights from physical modelling. *Geophysical Journal International*, 174(2), 719–732. https://doi.org/10.1111/j.1365-246x.2008.03826.x

Bull, A. L., McNamara, A. K., Becker, T. W., & Ritsema, J. (2010). Global scale models of the mantle flow field predicted by synthetic tomography models. *Physics of the Earth and Planetary Interiors*, 182(3–4), 129–138. https://doi.org/10.1016/j.pepi.2010.03.004

Cande, S. C., & Stegman, D. R. (2011). Indian and African plate motions driven by the push force of the réunion plume head. *Nature*, 475(7354), 47–52. https://doi.org/10.1038/nature10174

Capitanio, F. A., Faccenna, C., Zlotnik, S., & Stegman, D. R. (2011). Subduction dynamics and the origin of Andean orogeny and the Bolivian orocline. *Nature*, 480(7375), 83–86. https://doi.org/10.1038/nature10596

Cerpa, N. G., Guillaume, B., & Martinod, J. (2018). The interplay between overriding plate kinematics, slab dip and tectonics. *Geophysical Journal International*, 215(3), 1789–1802. https://doi.org/10.1093/gji/ggy365

Cerpa, N. G., Sigloch, K., Garel, F., Heuret, A., Davies, D. R., & Mihalynuk, M. G. (2022). The effect of a weak asthenospheric layer on surface kinematics, subduction dynamics and slab morphology in the lower mantle. *Journal of Geophysical Research: Solid Earth*, 127(8), e2022JB024494. https://doi.org/10.1029/2022JB024494

Chang, S.-J., Ferreira, A. M., & Faccenda, M. (2016). Upper- and mid-mantle interaction between the Samoan plume and the Tonga–Kermadec slabs. *Nature Communications*, 7(1), 10799. https://doi.org/10.1038/ncomms10799

Chen, Y.-W., Wu, J., & Suppe, J. (2019). Southward propagation of Nazca subduction along the Andes. *Nature*, 565(7740), 441–447. https://doi.org/10.1038/s41586-018-0860-1

Chertova, M. V., Spakman, W., & Steinberger, B. (2018). Mantle flow influence on subduction evolution. *Earth and Planetary Science Letters*, 489, 258–266. https://doi.org/10.1016/j.epsl.2018.02.038

Colli, L., Stotz, I., Bunge, H. P., Smethurst, M., Clark, S., Iaffaldano, G., et al. (2014). Rapid South Atlantic spreading changes and coeval vertical motion in surrounding continents: Evidence for temporal changes of pressure-driven upper mantle flow. *Tectonics*, 33(7), 1304–1321. https://doi.org/10.1002/2014TC003612

Coltice, N., Husson, L., Faccenna, C., & Arnould, M. (2019). What drives tectonic plates? *Science Advances*, 5(10), eaax4295. https://doi.org/10.1126/sciadv.aax4295

Conrad, C. P., & Behn, M. D. (2010). Constraints on lithosphere net rotation and asthenospheric viscosity from global mantle flow models and seismic anisotropy. *Geochemistry, Geophysics, Geosystems*, 11(5), Q05W05. https://doi.org/10.1029/2009gc002970

Conrad, C. P., & Lithgow-Bertelloni, C. (2002). How mantle slabs drive plate tectonics. *Science*, 298(5591), 207–209. https://doi.org/10.1126/science.1074161

Crespi, M., Cuffaro, M., Doglioni, C., Giannone, F., & Riguzzi, F. (2007). Space geodesy validation of the global lithospheric flow. *Geophysical Journal International*, 168(2), 491–506. https://doi.org/10.1111/j.1365-246x.2006.03226.x

Cuffaro, M., & Doglioni, C. (2007). Global kinematics in deep versus shallow hotspot reference frames. *Special Paper 430: Plates, Plumes and Planetary Processes*, 359–374. https://doi.org/10.1130/2007.2430(18)

Di Leo, J. F., Wookey, J., Hammond, J. O., Kendall, J.-M., Kaneshima, S., Inoue, H., et al. (2012). Mantle flow in regions of complex tectonics: Insights from Indonesia. *Geochemistry, Geophysics, Geosystems*, 13(12), Q12008. https://doi.org/10.1029/2012gc004417

Duarte, J. C., Schellart, W. P., & Cruden, A. R. (2013). Three-dimensional dynamic laboratory models of subduction with an overriding plate and variable interplate rheology. *Geophysical Journal International*, 195(1), 47–66. https://doi.org/10.1093/gji/ggt257
associated with thickening only for trench-orthogonal strain rates larger than $2.5 \times 10^{-17} s^{-1}$ as the upper plate also accommodates deformation through mild trench-parallel stretching. We also show that upper plate strain rates are three times larger when the lower plate is the one that is dragged toward the trench, as the mantle flow drives subduction instead of opposing it. In addition, the reaction force exerted by the backwall on the fixed upper plate, which may be taken as an analog to ridge push or opposite mantle drag, also favors trench-orthogonal shortening.

We finally show that the interaction of the slab with an impermeable boundary at 660-km depth has profound consequences on subduction dynamics, and in turn on upper plate deformation. For similar free plate velocities, slightly larger strain rates are recorded when the slab has already reached the 660-km discontinuity. Deformation regimes observed along present-day subduction zones are compatible with the findings from our models and from computations of present-day integrated drag forces at the trench. This suggests that the drag exerted by the convective mantle may be a key ingredient in controlling the upper plate deformation regime.

## Data Availability Statement

Raw data associated with this paper (movies of the top and lateral views of each experiment) are available following this link: https://doi.org/10.57932/e21b6394-6557-4ba7-ac6a-6ecce61d5c69.

## Acknowledgments

We would like to express our gratitude to the Associate Editor, Ernst Willingshofer, and the reviewers, Wim Spakman and the anonymous reviewer, for their constructive comments on this manuscript. Their feedback has greatly contributed to enhancing the quality of this paper.

## References

Arriagada, C., Cobbold, P. R., & Roperch, P. (2006). Salar de Atacama Basin: A record of compressional tectonics in the Central Andes since the mid-Cretaceous. *Tectonics*, 25(1), TC1008. https://doi.org/10.1029/2004TC001770

Becker, T. W., Kustowski, B., & Ekström, G. (2008). Radial seismic anisotropy as a constraint for upper mantle rheology. *Earth and Planetary Science Letters*, 267(1–2), 213–227. https://doi.org/10.1016/j.epsl.2007.11.038

Behr, W. M., Holt, A. F., Becker, T. W., & Faccenna, C. (2022). The effects of plate interface rheology on subduction kinematics and dynamics. *Geophysical Journal International*, 230(2), 796–812. https://doi.org/10.1093/gji/ggac075

Boschman, L. M. (2021). Andean mountain building since the late cretaceous: A paleoelevation reconstruction. *Earth-Science Reviews*, 220, 103640. https://doi.org/10.1016/j.earscirev.2021.103640

Boutelier, D. A., & Cruden, A. R. (2008). Impact of regional mantle flow on subducting plate geometry and interplate stress: Insights from physical modelling. *Geophysical Journal International*, 174(2), 719–732. https://doi.org/10.1111/j.1365-246x.2008.03826.x

Bull, A. L., McNamara, A. K., Becker, T. W., & Ritsema, J. (2010). Global scale models of the mantle flow field predicted by synthetic tomography models. *Physics of the Earth and Planetary Interiors*, 182(3–4), 129–138. https://doi.org/10.1016/j.pepi.2010.03.004

Cande, S. C., & Stegman, D. R. (2011). Indian and African plate motions driven by the push force of the réunion plume head. *Nature*, 475(7354), 47–52. https://doi.org/10.1038/nature10174

Capitanio, F. A., Faccenna, C., Zlotnik, S., & Stegman, D. R. (2011). Subduction dynamics and the origin of Andean orogeny and the Bolivian orocline. *Nature*, 480(7375), 83–86. https://doi.org/10.1038/nature10596

Cerpa, N. G., Guillaume, B., & Martinod, J. (2018). The interplay between overriding plate kinematics, slab dip and tectonics. *Geophysical Journal International*, 215(3), 1789–1802. https://doi.org/10.1093/gji/ggy365

Cerpa, N. G., Sigloch, K., Garel, F., Heuret, A., Davies, D. R., & Mihalynuk, M. G. (2022). The effect of a weak asthenospheric layer on surface kinematics, subduction dynamics and slab morphology in the lower mantle. *Journal of Geophysical Research: Solid Earth*, 127(8), e2022JB024494. https://doi.org/10.1029/2022JB024494

Chang, S.-J., Ferreira, A. M., & Faccenda, M. (2016). Upper- and mid-mantle interaction between the Samoan plume and the Tonga–Kermadec slabs. *Nature Communications*, 7(1), 10799. https://doi.org/10.1038/ncomms10799

Chen, Y.-W., Wu, J., & Suppe, J. (2019). Southward propagation of Nazca subduction along the Andes. *Nature*, 565(7740), 441–447. https://doi.org/10.1038/s41586-018-0860-1

Chertova, M. V., Spakman, W., & Steinberger, B. (2018). Mantle flow influence on subduction evolution. *Earth and Planetary Science Letters*, 489, 258–266. https://doi.org/10.1016/j.epsl.2018.02.038

Colli, L., Stotz, I., Bunge, H. P., Smethurst, M., Clark, S., Iaffaldano, G., et al. (2014). Rapid South Atlantic spreading changes and coeval vertical motion in surrounding continents: Evidence for temporal changes of pressure-driven upper mantle flow. *Tectonics*, 33(7), 1304–1321. https://doi.org/10.1002/2014TC003612

Coltice, N., Husson, L., Faccenna, C., & Arnould, M. (2019). What drives tectonic plates? *Science Advances*, 5(10), eaax4295. https://doi.org/10.1126/sciadv.aax4295

Conrad, C. P., & Behn, M. D. (2010). Constraints on lithosphere net rotation and asthenospheric viscosity from global mantle flow models and seismic anisotropy. *Geochemistry, Geophysics, Geosystems*, 11(5), Q05W05. https://doi.org/10.1029/2009gc002970

Conrad, C. P., & Lithgow-Bertelloni, C. (2002). How mantle slabs drive plate tectonics. *Science*, 298(5591), 207–209. https://doi.org/10.1126/science.1074161

Crespi, M., Cuffaro, M., Doglioni, C., Giannone, F., & Riguzzi, F. (2007). Space geodesy validation of the global lithospheric flow. *Geophysical Journal International*, 168(2), 491–506. https://doi.org/10.1111/j.1365-246x.2006.03226.x

Cuffaro, M., & Doglioni, C. (2007). Global kinematics in deep versus shallow hotspot reference frames. *Special Paper 430: Plates, Plumes and Planetary Processes*, 359–374. https://doi.org/10.1130/2007.2430(18)

Di Leo, J. F., Wookey, J., Hammond, J. O., Kendall, J.-M., Kaneshima, S., Inoue, H., et al. (2012). Mantle flow in regions of complex tectonics: Insights from Indonesia. *Geochemistry, Geophysics, Geosystems*, 13(12), Q12008. https://doi.org/10.1029/2012gc004417

Duarte, J. C., Schellart, W. P., & Cruden, A. R. (2013). Three-dimensional dynamic laboratory models of subduction with an overriding plate and variable interplate rheology. *Geophysical Journal International*, 195(1), 47–66. https://doi.org/10.1093/gji/ggt257
Duarte, J. C., Schellart, W. P., & Cruden, A. R. (2014). Rheology of petrolatum–paraffin oil mixtures: Applications to analogue modelling of geological processes. *Journal of Structural Geology*, **63**, 1–11. https://doi.org/10.1016/j.jsg.2014.02.004

Duarte, J. C., Schellart, W. P., & Cruden, A. R. (2015). How weak is the subduction zone interface? *Geophysical Research Letters*, **42**(8), 2664–2673. https://doi.org/10.1002/2014gl062876

Espurt, N., Funiciello, F., Martinod, J., Guillaume, B., Regard, V., Faccenna, C., & Brusset, S. (2008). Flat subduction dynamics and deformation of the South American Plate: Insights from analog modeling. *Tectonics*, **27**(3). https://doi.org/10.1029/2007tc002175

Faccenna, C., Becker, T. W., Conrad, C. P., & Husson, L. (2013). Mountain building and mantle dynamics. *Tectonics*, **32**(1), 80–93. https://doi.org/10.1029/2012tc003176

Faccenna, C., Oncken, O., Holt, A. F., & Becker, T. W. (2017). Initiation of the Andean orography by lower mantle subduction. *Earth and Planetary Science Letters*, **463**, 189–201. https://doi.org/10.1016/j.epsl.2017.01.041

Fernández-García, C. (2019). *Processus de déformation lors de la rupture latérale d'un slab: Apports de la modélisation analogique et perspectives pour la tectonique méditerranéenne* (PhD thesis) (p. 200). Université de Rennes 1.

Ficini, E., Dal Zilio, L., Doglioni, C., & Gerya, T. V. (2017). Horizontal mantle flow controls subduction dynamics. *Scientific Reports*, **7**(1), 7550. https://doi.org/10.1038/s41598-017-06551-y

Forsyth, D., & Uyeda, S. (1975). On the relative importance of the driving forces of plate motion. *Geophysical Journal International*, **43**(1), 163–200. https://doi.org/10.1111/j.1365-246x.1975.tb00631.x

Forte, A. M., & Mitrovica, J. X. (1996). New inferences of mantle viscosity from joint inversion of long-wavelength mantle convection and post-glacial rebound data. *Geophysical Research Letters*, **23**(10), 1147–1150. https://doi.org/10.1029/96g100964

Funiciello, F., Faccenna, C., Giardini, D., & Regenauer-Lieb, K. (2003). Dynamics of retreating slabs: 2. Insights from three-dimensional laboratory experiments. *Journal of Geophysical Research*, **108**(B4), 2207. https://doi.org/10.1029/2001jb000896

Funiciello, F., Faccenna, C., Heuret, A., Lallemand, S., Di Giuseppe, E., & Becker, T. W. (2008). Trench migration, net rotation and slab–mantle coupling. *Earth and Planetary Science Letters*, **271**(1–4), 233–240. https://doi.org/10.1016/j.epsl.2008.04.006

Guillaume, B., Funiciello, F., & Faccenna, C. (2012). Interplays between mantle flow and slab pull at subduction zones in 3D. *Journal of Geophysical Research: Solid Earth*, **126**(5), e2020JB021574. https://doi.org/10.1029/2020JB021574

Guillaume, B., Hertgen, S., Martinod, J., & Cerpa, N. G. (2018). Slab dip, surface tectonics: How and when do they change following an acceleration/slow down of the overriding plate? *Tectonophysics*, **726**, 110–120. https://doi.org/10.1016/j.tecto.2018.01.030

Guillaume, B., Martinod, J., & Espurt, N. (2009). Variations of slab dip and overriding plate tectonics during subduction: Insights from analogue modelling. *Tectonophysics*, **463**(1–4), 167–174. https://doi.org/10.1016/j.tecto.2008.09.043

Guillaume, B., Moroni, M., Funiciello, F., Martinod, J., & Faccenna, C. (2010). Mantle flow and dynamic topography associated with slab window opening: Insights from laboratory models. *Tectonophysics*, **496**(1–4), 83–98. https://doi.org/10.1016/j.tecto.2010.10.014

Habel, T., Replumaz, A., Guillaume, B., Simoes, M., Geffroy, T., Kermarrec, J. J., & Lacassin, R. (2023). Upper-plate shortening and mountain-building in the context of mantle-driven oceanic subduction. *Tektonika*, **1**(2), 158–176. https://doi.org/10.55575/tektonika2023.1.2.39

Hager, B. H. (1984). Subducted slabs and the geoid: Constraints on mantle rheology and flow. *Journal of Geophysical Research*, **89**(B7), 6003–6015. https://doi.org/10.1029/jb089ib07p06003

Hager, B. H., & O'Connell, R. J. (1979). Kinematic models of large-scale flow in the Earth's mantle. *Journal of Geophysical Research*, **84**(B3), 1031–1048. https://doi.org/10.1029/jb084ib03p01031

Hager, B. H., & Richards, M. A. (1989). Long-wavelength variations in Earth's geoid: Physical models and dynamical implications. *Philosophical Transactions of the Royal Society of London - Series A: Mathematical and Physical Sciences*, **328**(1599), 309–327.

Hayes, G. P., Moore, G. L., Portner, D. E., Hearne, M., Flamme, H., Furtney, M., & Smoczyk, G. M. (2018). Slab2, a comprehensive subduction zone geometry model. *Science*, **362**(6410), 58–61. https://doi.org/10.1126/science.aat4723

Hertgen, S., Yamato, P., Guillaume, B., Magni, V., Schliffke, N., & Hunen, J. (2020). Influence of the thickness of the overriding plate on Convergence zone dynamics. *Geochemistry, Geophysics, Geosystems*, **21**(2). https://doi.org/10.1029/2019gc008678

Heuret, A., & Lallemand, S. (2005). Plate motions, slab dynamics and back-arc deformation. *Physics of the Earth and Planetary Interiors*, **149**(1–2), 31–51. https://doi.org/10.1016/j.pepi.2004.08.022

Höink, T., Jellinek, A. M., & Lenardic, A. (2011). Viscous coupling at the lithosphere-asthenosphere boundary. *Geochemistry, Geophysics, Geosystems*, **12**(10), Q0AK02. https://doi.org/10.1029/2011GC003698

Holt, A. F., & Becker, T. W. (2016). The effect of a power-law mantle viscosity on trench retreat rate. *Geophysical Journal International*, **208**(1), 491–507. https://doi.org/10.1093/gji/ggw392

Holt, A. F., Becker, T. W., & Buffett, B. A. (2015). Trench migration and overriding plate stress in dynamic subduction models. *Geophysical Journal International*, **201**(1), 172–192. https://doi.org/10.1093/gji/ggw011

Husson, L. (2012). Trench migration and upper plate strain over a convecting mantle. *Physics of the Earth and Planetary Interiors*, **212**–213, 32–43. https://doi.org/10.1016/j.pepi.2012.09.006

Jaillard, E. (2022). Late Cretaceous-Paleogene orogenic build-up of the Ecuadorian Andes: Review and discussion. *Earth-Science Reviews*, **230**, 104033. https://doi.org/10.1016/j.earscirev.2022.104033

Karato, S. I. (2010). Rheology of the deep upper mantle and its implications for the preservation of the continental roots: A review. *Tectonophysics*, **481**(1–4), 82–98. https://doi.org/10.1016/j.tecto.2009.04.011

Kincaid, C., & Griffiths, R. W. (2003). Laboratory models of the thermal evolution of the mantle during rollback subduction. *Nature*, **425**(6953), 58–62. https://doi.org/10.1038/nature01923

King, S. D. (2016). An evolving view of transition zone and midmantle viscosity. *Geochemistry, Geophysics, Geosystems*, **17**(3), 1234–1237. https://doi.org/10.1002/2016gc006279

King, S. D., & Hager, B. H. (1994). Subducted slabs and the geoid: 1. Numerical experiments with temperature-dependent viscosity. *Journal of Geophysical Research*, **99**(B10), 19843–19852. https://doi.org/10.1029/94JB01552

Király, Á., Holt, A. F., Funiciello, F., Faccenna, C., & Capitanio, F. A. (2018). Modeling slab-slab interactions: Dynamics of outward dipping double-sided subduction systems. *Geochemistry, Geophysics, Geosystems*, **19**(3), 693–714. https://doi.org/10.1002/2017gc007199

Kreemer, C. (2009). Absolute plate motions constrained by shear wave splitting orientations with implications for hot spot motions and Mantle flow. *Journal of Geophysical Research*, **114**(B10), B10405. https://doi.org/10.1029/2009JB006416

Lallemand, S., & Heuret, A. (2017). Subduction zones parameters. *Reference Module in Earth Systems and Environmental Sciences*. https://doi.org/10.1016/b978-0-12-409548-9.09495-1

Lallemand, S., Heuret, A., & Boutelier, D. (2005). On the relationships between slab dip, back-arc stress, Upper Plate absolute motion, and crustal nature in subduction zones. *Geochemistry, Geophysics, Geosystems*, **6**(9). https://doi.org/10.1029/2005gc000917

Lithgow-Bertelloni, C., & Richards, M. A. (1998). The dynamics of Cenozoic and Mesozoic plate motions. *Reviews of Geophysics*, **36**(1), 27–78. https://doi.org/10.1029/97rg02282
Duarte, J. C., Schellart, W. P., & Cruden, A. R. (2014). Rheology of petrolatum–paraffin oil mixtures: Applications to analogue modelling of geological processes. *Journal of Structural Geology*, **63**, 1–11. https://doi.org/10.1016/j.jsg.2014.02.004

Duarte, J. C., Schellart, W. P., & Cruden, A. R. (2015). How weak is the subduction zone interface? *Geophysical Research Letters*, **42**(8), 2664–2673. https://doi.org/10.1002/2014gl062876

Espurt, N., Funiciello, F., Martinod, J., Guillaume, B., Regard, V., Faccenna, C., & Brusset, S. (2008). Flat subduction dynamics and deformation of the South American Plate: Insights from analog modeling. *Tectonics*, **27**(3). https://doi.org/10.1029/2007tc002175

Faccenna, C., Becker, T. W., Conrad, C. P., & Husson, L. (2013). Mountain building and mantle dynamics. *Tectonics*, **32**(1), 80–93. https://doi.org/10.1029/2012tc003176

Faccenna, C., Oncken, O., Holt, A. F., & Becker, T. W. (2017). Initiation of the Andean orography by lower mantle subduction. *Earth and Planetary Science Letters*, **463**, 189–201. https://doi.org/10.1016/j.epsl.2017.01.041

Fernández-García, C. (2019). *Processus de déformation lors de la rupture latérale d'un slab: Apports de la modélisation analogique et perspectives pour la tectonique méditerranéenne* (PhD thesis) (p. 200). Université de Rennes 1.

Ficini, E., Dal Zilio, L., Doglioni, C., & Gerya, T. V. (2017). Horizontal mantle flow controls subduction dynamics. *Scientific Reports*, **7**(1), 7550. https://doi.org/10.1038/s41598-017-06551-y

Forsyth, D., & Uyeda, S. (1975). On the relative importance of the driving forces of plate motion. *Geophysical Journal International*, **43**(1), 163–200. https://doi.org/10.1111/j.1365-246x.1975.tb00631.x

Forte, A. M., & Mitrovica, J. X. (1996). New inferences of mantle viscosity from joint inversion of long-wavelength mantle convection and post-glacial rebound data. *Geophysical Research Letters*, **23**(10), 1147–1150. https://doi.org/10.1029/96g100964

Funiciello, F., Faccenna, C., Giardini, D., & Regenauer-Lieb, K. (2003). Dynamics of retreating slabs: 2. Insights from three-dimensional laboratory experiments. *Journal of Geophysical Research*, **108**(B4), 2207. https://doi.org/10.1029/2001jb000896

Funiciello, F., Faccenna, C., Heuret, A., Lallemand, S., Di Giuseppe, E., & Becker, T. W. (2008). Trench migration, net rotation and slab–mantle coupling. *Earth and Planetary Science Letters*, **271**(1–4), 233–240. https://doi.org/10.1016/j.epsl.2008.04.006

Guillaume, B., Funiciello, F., & Faccenna, C. (2012). Interplays between mantle flow and slab pull at subduction zones in 3D. *Journal of Geophysical Research: Solid Earth*, **126**(5), e2020JB021574. https://doi.org/10.1029/2020JB021574

Guillaume, B., Hertgen, S., Martinod, J., & Cerpa, N. G. (2018). Slab dip, surface tectonics: How and when do they change following an acceleration/slow down of the overriding plate? *Tectonophysics*, **726**, 110–120. https://doi.org/10.1016/j.tecto.2018.01.030

Guillaume, B., Martinod, J., & Espurt, N. (2009). Variations of slab dip and overriding plate tectonics during subduction: Insights from analogue modelling. *Tectonophysics*, **463**(1–4), 167–174. https://doi.org/10.1016/j.tecto.2008.09.043

Guillaume, B., Moroni, M., Funiciello, F., Martinod, J., & Faccenna, C. (2010). Mantle flow and dynamic topography associated with slab window opening: Insights from laboratory models. *Tectonophysics*, **496**(1–4), 83–98. https://doi.org/10.1016/j.tecto.2010.10.014

Habel, T., Replumaz, A., Guillaume, B., Simoes, M., Geffroy, T., Kermarrec, J. J., & Lacassin, R. (2023). Upper-plate shortening and mountain-building in the context of mantle-driven oceanic subduction. *Tektonika*, **1**(2), 158–176. https://doi.org/10.55575/tektonika2023.1.2.39

Hager, B. H. (1984). Subducted slabs and the geoid: Constraints on mantle rheology and flow. *Journal of Geophysical Research*, **89**(B7), 6003–6015. https://doi.org/10.1029/jb089ib07p06003

Hager, B. H., & O'Connell, R. J. (1979). Kinematic models of large-scale flow in the Earth's mantle. *Journal of Geophysical Research*, **84**(B3), 1031–1048. https://doi.org/10.1029/jb084ib03p01031

Hager, B. H., & Richards, M. A. (1989). Long-wavelength variations in Earth's geoid: Physical models and dynamical implications. *Philosophical Transactions of the Royal Society of London - Series A: Mathematical and Physical Sciences*, **328**(1599), 309–327.

Hayes, G. P., Moore, G. L., Portner, D. E., Hearne, M., Flamme, H., Furtney, M., & Smoczyk, G. M. (2018). Slab2, a comprehensive subduction zone geometry model. *Science*, **362**(6410), 58–61. https://doi.org/10.1126/science.aat4723

Hertgen, S., Yamato, P., Guillaume, B., Magni, V., Schliffke, N., & Hunen, J. (2020). Influence of the thickness of the overriding plate on Convergence zone dynamics. *Geochemistry, Geophysics, Geosystems*, **21**(2). https://doi.org/10.1029/2019gc008678

Heuret, A., & Lallemand, S. (2005). Plate motions, slab dynamics and back-arc deformation. *Physics of the Earth and Planetary Interiors*, **149**(1–2), 31–51. https://doi.org/10.1016/j.pepi.2004.08.022

Höink, T., Jellinek, A. M., & Lenardic, A. (2011). Viscous coupling at the lithosphere-asthenosphere boundary. *Geochemistry, Geophysics, Geosystems*, **12**(10), Q0AK02. https://doi.org/10.1029/2011GC003698

Holt, A. F., & Becker, T. W. (2016). The effect of a power-law mantle viscosity on trench retreat rate. *Geophysical Journal International*, **208**(1), 491–507. https://doi.org/10.1093/gji/ggw392

Holt, A. F., Becker, T. W., & Buffett, B. A. (2015). Trench migration and overriding plate stress in dynamic subduction models. *Geophysical Journal International*, **201**(1), 172–192. https://doi.org/10.1093/gji/ggw011

Husson, L. (2012). Trench migration and upper plate strain over a convecting mantle. *Physics of the Earth and Planetary Interiors*, **212–213**, 32–43. https://doi.org/10.1016/j.pepi.2012.09.006

Jaillard, E. (2022). Late Cretaceous-Paleogene orogenic build-up of the Ecuadorian Andes: Review and discussion. *Earth-Science Reviews*, **230**, 104033. https://doi.org/10.1016/j.earscirev.2022.104033

Karato, S. I. (2010). Rheology of the deep upper mantle and its implications for the preservation of the continental roots: A review. *Tectonophysics*, **481**(1–4), 82–98. https://doi.org/10.1016/j.tecto.2009.04.011

Kincaid, C., & Griffiths, R. W. (2003). Laboratory models of the thermal evolution of the mantle during rollback subduction. *Nature*, **425**(6953), 58–62. https://doi.org/10.1038/nature01923

King, S. D. (2016). An evolving view of transition zone and midmantle viscosity. *Geochemistry, Geophysics, Geosystems*, **17**(3), 1234–1237. https://doi.org/10.1002/2016gc006279

King, S. D., & Hager, B. H. (1994). Subducted slabs and the geoid: 1. Numerical experiments with temperature-dependent viscosity. *Journal of Geophysical Research*, **99**(B10), 19843–19852. https://doi.org/10.1029/94JB01552

Király, Á., Holt, A. F., Funiciello, F., Faccenna, C., & Capitanio, F. A. (2018). Modeling slab-slab interactions: Dynamics of outward dipping double-sided subduction systems. *Geochemistry, Geophysics, Geosystems*, **19**(3), 693–714. https://doi.org/10.1002/2017gc007199

Kreemer, C. (2009). Absolute plate motions constrained by shear wave splitting orientations with implications for hot spot motions and Mantle flow. *Journal of Geophysical Research*, **114**(B10), B10405. https://doi.org/10.1029/2009JB006416

Lallemand, S., & Heuret, A. (2017). Subduction zones parameters. *Reference Module in Earth Systems and Environmental Sciences*. https://doi.org/10.1016/b978-0-12-409548-9.09495-1

Lallemand, S., Heuret, A., & Boutelier, D. (2005). On the relationships between slab dip, back-arc stress, Upper Plate absolute motion, and crustal nature in subduction zones. *Geochemistry, Geophysics, Geosystems*, **6**(9). https://doi.org/10.1029/2005gc000917

Lithgow-Bertelloni, C., & Richards, M. A. (1998). The dynamics of Cenozoic and Mesozoic plate motions. *Reviews of Geophysics*, **36**(1), 27–78. https://doi.org/10.1029/97rg02282
Loiselet, C., Husson, L., & Braun, J. (2009). From longitudinal slab curvature to slab rheology. *Geology*, *37*(8), 747–750. https://doi.org/10.1130/g30052a.1

Martinod, J., Gérault, M., Husson, L., & Regard, V. (2020). Widening of the Andes: An interplay between subduction dynamics and crustal wedge tectonics. *Earth-Science Reviews*, *204*, 103170. https://doi.org/10.1016/j.earscirev.2020.103170

Martinod, J., Husson, L., Roperch, P., Guillaume, B., & Espurt, N. (2010). Horizontal subduction zones, convergence velocity and the building of the Andes. *Earth and Planetary Science Letters*, *299*(3–4), 299–309. https://doi.org/10.1016/j.epsl.2010.09.010

Menant, A., Sternai, P., Jolivet, L., Guillou-Frottier, L., & Gerya, T. (2016). 3D numerical modeling of mantle flow, crustal dynamics and magma genesis associated with slab roll-back and tearing: The eastern Mediterranean case. *Earth and Planetary Science Letters*, *442*, 93–107. https://doi.org/10.1016/j.epsl.2016.03.002

Mitrovica, J. X., & Forte, A. M. (2004). A new inference of mantle viscosity based upon joint inversion of convection and glacial isostatic adjustment data. *Earth and Planetary Science Letters*, *225*(1–2), 177–189. https://doi.org/10.1016/j.epsl.2004.06.005

Molnar, P., & Atwater, T. (1978). Interarc spreading and Cordilleran tectonics as alternates related to the age of subducted oceanic lithosphere. *Earth and Planetary Science Letters*, *41*(3), 330–340. https://doi.org/10.1016/0012-821x(78)90187-5

Molnar, P., & Gray, D. (1979). Subduction of continental lithosphere: Some constraints and uncertainties. *Geology*, *7*(1), 58–62. https://doi.org/10.1130/0091-7613(1979)7<58:soclsc>2.0.co;2

Müller, R. D., Seton, M., Zahirovic, S., Williams, S. E., Matthews, K. J., Wright, N. M., et al. (2016). Ocean basin evolution and global-scale plate reorganization events since Pangea Breakup. *Annual Review of Earth and Planetary Sciences*, *44*(1), 107–138. https://doi.org/10.1146/annurev-earth-060115-012211

Obayashi, M., Yoshimitsu, J., Nolet, G., Fukao, Y., Shiobara, H., Sugioka, H., et al. (2013). Finite frequency whole mantle *P* wave tomography: Improvement of subducted slab images. *Geophysical Research Letters*, *40*(21), 5652–5657. https://doi.org/10.1002/2013gl057401

Obrebski, M., Allen, R. M., Xue, M., & Hung, S.-H. (2010). Slab-plume interaction beneath the Pacific Northwest. *Geophysical Research Letters*, *37*(14), L14305. https://doi.org/10.1029/2010gl043489

Palano, M., Piromallo, C., & Chiarabba, C. (2017). Surface imprint of toroidal flow at retreating slab edges: The first geodetic evidence in the Calabrian subduction system. *Geophysical Research Letters*, *44*(2), 845–853. https://doi.org/10.1002/2016gl071452

Pasyanos, M. E., Masters, T. G., Laske, G., & Ma, Z. (2014). LITHO1.0: An updated crust and lithospheric model of the Earth. *Journal of Geophysical Research: Solid Earth*, *119*(3), 2153–2173. https://doi.org/10.1002/2013jb010626

Peral, M., Király, Á., Zlotnik, S., Funiciello, F., Fernández, M., Faccenna, C., & Vergés, J. (2018). Opposite subduction polarity in adjacent plate segments. *Tectonics*, *37*(9), 3285–3302. https://doi.org/10.1029/2017tc004896

Piromallo, C., Becker, T. W., Funiciello, F., & Faccenna, C. (2006). Three-dimensional instantaneous mantle flow induced by subduction. *Geophysical Research Letters*, *33*(8), L08304. https://doi.org/10.1029/2005gl025390

Rappisi, F., VanderBeek, B. P., Faccenda, M., Morelli, A., & Molinari, I. (2022). Slab geometry and upper mantle flow patterns in the Central Mediterranean from 3D anisotropic *P*-wave tomography. *Journal of Geophysical Research: Solid Earth*, *127*(5), e2021JB023488. https://doi.org/10.1029/2021jb023488

Ricard, Y., Doglioni, C., & Sabadini, R. (1991). Differential rotation between lithosphere and mantle: A consequence of lateral mantle viscosity variations. *Journal of Geophysical Research*, *96*(B5), 8407–8415. https://doi.org/10.1029/91jb00204

Rodríguez-González, J., Negredo, A. M., & Carminati, E. (2014). Slab-mantle flow interaction: Influence on subduction dynamics and duration. *Terra Nova*, *26*(4), 265–272. https://doi.org/10.1111/ter.12095

Rudolf, M., Boutelier, D., Rosenau, M., Schreurs, G., & Oncken, O. (2016). Rheological benchmark of silicone oils used for analog modeling of short- and long-term lithospheric deformation. *Tectonophysics*, *684*, 12–22. https://doi.org/10.1016/j.tecto.2015.11.028

Schellart, W. P. (2005). Influence of the subducting plate velocity on the geometry of the slab and migration of the subduction hinge. *Earth and Planetary Science Letters*, *231*(3–4), 197–219. https://doi.org/10.1016/j.epsl.2004.12.019

Schellart, W. P. (2010). Evolution of subduction zone curvature and its dependence on the trench velocity and the slab to upper mantle viscosity ratio. *Journal of Geophysical Research*, *115*(B11), B11406. https://doi.org/10.1029/2009jb006643

Schellart, W. P. (2011). Rheology and density of glucose syrup and honey: Determining their suitability for usage in analogue and fluid dynamic models of geological processes. *Journal of Structural Geology*, *33*(6), 1079–1088. https://doi.org/10.1016/j.jsg.2011.03.013

Schellart, W. P. (2017). Andean mountain building and magmatic arc migration driven by subduction-induced whole mantle flow. *Nature Communications*, *8*(1), 2010. https://doi.org/10.1038/s41467-017-01847-z

Schellart, W. P., & Moresi, L. (2013). A new driving mechanism for backarc extension and backarc shortening through slab sinking induced toroidal and poloidal mantle flow: Results from dynamic subduction models with an overriding plate. *Journal of Geophysical Research: Solid Earth*, *118*(6), 3221–3248. https://doi.org/10.1002/jgrb.50173

Schepers, G., Van Hinsbergen, D. J., Spakman, W., Kosters, M. E., Boschman, L. M., & McQuarrie, N. (2017). South-American plate advance and forced Andean trench retreat as drivers for transient flat subduction episodes. *Nature Communications*, *8*(1), 15249. https://doi.org/10.1038/ncomms15249

Semple, A., & Lenardic, A. (2020). The robustness of pressure-driven asthenospheric flow in mantle convection models with plate-like behavior. *Geophysical Research Letters*, *47*(17), e2020GL089556. https://doi.org/10.1029/2020GL089556

Stotz, I. L., Vilacís, B., Hayek, J. N., Carena, S., & Bunge, H. P. (2023). Plume driven plate motion changes: New insights from the South Atlantic realm. *Journal of South American Earth Sciences*, *124*, 104257. https://doi.org/10.1016/j.jsames.2023.104257

Turcotte, D. L., & Schubert, G. (2002). *Geodynamics*. Cambridge University Press. https://doi.org/10.1017/cbo9780511807442

Weismüller, J., Gmeiner, B., Ghelichkhan, S., Huber, M., John, L., Wohlhmuth, B., et al. (2015). Fast asthenosphere motion in high-resolution global mantle flow models. *Geophysical Research Letters*, *42*(18), 7429–7435. https://doi.org/10.1002/2015GL063727

Wu, B., Conrad, C. P., Heuret, A., Lithgow-Bertelloni, C., & Lallemand, S. (2008). Reconciling strong slab pull and weak plate bending: The plate motion constraint on the strength of mantle slabs. *Earth and Planetary Science Letters*, *272*(1–2), 412–421. https://doi.org/10.1016/j.epsl.2008.05.009

Xue, K., Schellart, W. P., & Strak, V. (2022). Overriding plate deformation and topography during slab rollback and slab rollover: Insights from subduction experiments. *Tectonics*, *41*(2), e2021TC007089. https://doi.org/10.1029/2021tc007089

Zhong, S., & Gurnis, M. (1994). Controls on trench topography from dynamic models of subducted slabs. *Journal of Geophysical Research*, *99*(B8), 15683–15695. https://doi.org/10.1029/94jb00809
Loiselet, C., Husson, L., & Braun, J. (2009). From longitudinal slab curvature to slab rheology. *Geology*, *37*(8), 747–750. https://doi.org/10.1130/g30052a.1

Martinod, J., Gérault, M., Husson, L., & Regard, V. (2020). Widening of the Andes: An interplay between subduction dynamics and crustal wedge tectonics. *Earth-Science Reviews*, *204*, 103170. https://doi.org/10.1016/j.earscirev.2020.103170

Martinod, J., Husson, L., Roperch, P., Guillaume, B., & Espurt, N. (2010). Horizontal subduction zones, convergence velocity and the building of the Andes. *Earth and Planetary Science Letters*, *299*(3–4), 299–309. https://doi.org/10.1016/j.epsl.2010.09.010

Menant, A., Sternai, P., Jolivet, L., Guillou-Frottier, L., & Gerya, T. (2016). 3D numerical modeling of mantle flow, crustal dynamics and magma genesis associated with slab roll-back and tearing: The eastern Mediterranean case. *Earth and Planetary Science Letters*, *442*, 93–107. https://doi.org/10.1016/j.epsl.2016.03.002

Mitrovica, J. X., & Forte, A. M. (2004). A new inference of mantle viscosity based upon joint inversion of convection and glacial isostatic adjustment data. *Earth and Planetary Science Letters*, *225*(1–2), 177–189. https://doi.org/10.1016/j.epsl.2004.06.005

Molnar, P., & Atwater, T. (1978). Interarc spreading and Cordilleran tectonics as alternates related to the age of subducted oceanic lithosphere. *Earth and Planetary Science Letters*, *41*(3), 330–340. https://doi.org/10.1016/0012-821x(78)90187-5

Molnar, P., & Gray, D. (1979). Subduction of continental lithosphere: Some constraints and uncertainties. *Geology*, *7*(1), 58–62. https://doi.org/10.1130/0091-7613(1979)7<58:soclsc>2.0.co;2

Müller, R. D., Seton, M., Zahirovic, S., Williams, S. E., Matthews, K. J., Wright, N. M., et al. (2016). Ocean basin evolution and global-scale plate reorganization events since Pangea Breakup. *Annual Review of Earth and Planetary Sciences*, *44*(1), 107–138. https://doi.org/10.1146/annurev-earth-060115-012211

Obayashi, M., Yoshimitsu, J., Nolet, G., Fukao, Y., Shiobara, H., Sugioka, H., et al. (2013). Finite frequency whole mantle P wave tomography: Improvement of subducted slab images. *Geophysical Research Letters*, *40*(21), 5652–5657. https://doi.org/10.1002/2013gl057401

Obrebski, M., Allen, R. M., Xue, M., & Hung, S.-H. (2010). Slab-plume interaction beneath the Pacific Northwest. *Geophysical Research Letters*, *37*(14), L14305. https://doi.org/10.1029/2010gl043489

Palano, M., Piromallo, C., & Chiarabba, C. (2017). Surface imprint of toroidal flow at retreating slab edges: The first geodetic evidence in the Calabrian subduction system. *Geophysical Research Letters*, *44*(2), 845–853. https://doi.org/10.1002/2016gl071452

Pasyanos, M. E., Masters, T. G., Laske, G., & Ma, Z. (2014). LITHO1.0: An updated crust and lithospheric model of the Earth. *Journal of Geophysical Research: Solid Earth*, *119*(3), 2153–2173. https://doi.org/10.1002/2013jb010626

Peral, M., Király, Á., Zlotnik, S., Funiciello, F., Fernández, M., Faccenna, C., & Vergés, J. (2018). Opposite subduction polarity in adjacent plate segments. *Tectonics*, *37*(9), 3285–3302. https://doi.org/10.1029/2017tc004896

Piromallo, C., Becker, T. W., Funiciello, F., & Faccenna, C. (2006). Three-dimensional instantaneous mantle flow induced by subduction. *Geophysical Research Letters*, *33*(8), L08304. https://doi.org/10.1029/2005gl025390

Rappisi, F., VanderBeek, B. P., Faccenda, M., Morelli, A., & Molinari, I. (2022). Slab geometry and upper mantle flow patterns in the Central Mediterranean from 3D anisotropic P-wave tomography. *Journal of Geophysical Research: Solid Earth*, *127*(5), e2021JB023488. https://doi.org/10.1029/2021jb023488

Ricard, Y., Doglioni, C., & Sabadini, R. (1991). Differential rotation between lithosphere and mantle: A consequence of lateral mantle viscosity variations. *Journal of Geophysical Research*, *96*(B5), 8407–8415. https://doi.org/10.1029/91jb00204

Rodríguez-González, J., Negredo, A. M., & Carminati, E. (2014). Slab-mantle flow interaction: Influence on subduction dynamics and duration. *Terra Nova*, *26*(4), 265–272. https://doi.org/10.1111/ter.12095

Rudolf, M., Boutelier, D., Rosenau, M., Schreurs, G., & Oncken, O. (2016). Rheological benchmark of silicone oils used for analog modeling of short- and long-term lithospheric deformation. *Tectonophysics*, *684*, 12–22. https://doi.org/10.1016/j.tecto.2015.11.028

Schellart, W. P. (2005). Influence of the subducting plate velocity on the geometry of the slab and migration of the subduction hinge. *Earth and Planetary Science Letters*, *231*(3–4), 197–219. https://doi.org/10.1016/j.epsl.2004.12.019

Schellart, W. P. (2010). Evolution of subduction zone curvature and its dependence on the trench velocity and the slab to upper mantle viscosity ratio. *Journal of Geophysical Research*, *115*(B11), B11406. https://doi.org/10.1029/2009jb006643

Schellart, W. P. (2011). Rheology and density of glucose syrup and honey: Determining their suitability for usage in analogue and fluid dynamic models of geological processes. *Journal of Structural Geology*, *33*(6), 1079–1088. https://doi.org/10.1016/j.jsg.2011.03.013

Schellart, W. P. (2017). Andean mountain building and magmatic arc migration driven by subduction-induced whole mantle flow. *Nature Communications*, *8*(1), 2010. https://doi.org/10.1038/s41467-017-01847-z

Schellart, W. P., & Moresi, L. (2013). A new driving mechanism for backarc extension and backarc shortening through slab sinking induced toroidal and poloidal mantle flow: Results from dynamic subduction models with an overriding plate. *Journal of Geophysical Research: Solid Earth*, *118*(6), 3221–3248. https://doi.org/10.1002/jgrb.50173

Schepers, G., Van Hinsbergen, D. J., Spakman, W., Kosters, M. E., Boschman, L. M., & McQuarrie, N. (2017). South-American plate advance and forced Andean trench retreat as drivers for transient flat subduction episodes. *Nature Communications*, *8*(1), 15249. https://doi.org/10.1038/ncomms15249

Semple, A., & Lenardic, A. (2020). The robustness of pressure-driven asthenospheric flow in mantle convection models with plate-like behavior. *Geophysical Research Letters*, *47*(17), e2020GL089556. https://doi.org/10.1029/2020GL089556

Stotz, I. L., Vilacís, B., Hayek, J. N., Carena, S., & Bunge, H. P. (2023). Plume driven plate motion changes: New insights from the South Atlantic realm. *Journal of South American Earth Sciences*, *124*, 104257. https://doi.org/10.1016/j.jsames.2023.104257

Turcotte, D. L., & Schubert, G. (2002). *Geodynamics*. Cambridge University Press. https://doi.org/10.1017/cbo9780511807442

Weismüller, J., Gmeiner, B., Ghelichkhan, S., Huber, M., John, L., Wohlmuth, B., et al. (2015). Fast asthenosphere motion in high-resolution global mantle flow models. *Geophysical Research Letters*, *42*(18), 7429–7435. https://doi.org/10.1002/2015GL063727

Wu, B., Conrad, C. P., Heuret, A., Lithgow-Bertelloni, C., & Lallemand, S. (2008). Reconciling strong slab pull and weak plate bending: The plate motion constraint on the strength of mantle slabs. *Earth and Planetary Science Letters*, *272*(1–2), 412–421. https://doi.org/10.1016/j.epsl.2008.05.009

Xue, K., Schellart, W. P., & Strak, V. (2022). Overriding plate deformation and topography during slab rollback and slab rollover: Insights from subduction experiments. *Tectonics*, *41*(2), e2021TC007089. https://doi.org/10.1029/2021tc007089

Zhong, S., & Gurnis, M. (1994). Controls on trench topography from dynamic models of subducted slabs. *Journal of Geophysical Research*, *99*(B8), 15683–15695. https://doi.org/10.1029/94jb00809