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

HAL
open science

# Overriding Plate Deformation and Topography During Slab Rollback and Slab Rollover: Insights From Subduction Experiments

K. Xue, W. P. Schellart, V Strak

► To cite this version:

K. Xue, W. P. Schellart, V Strak. Overriding Plate Deformation and Topography During Slab Rollback and Slab Rollover: Insights From Subduction Experiments. Tectonics, 2022, 41 (2), 10.1029/2021TC007089. [https://doi.org/10.1029/2021TC007089](https://doi.org/10.1029/2021TC007089)

HAL Id: hal-03662532

https://hal.science/hal-03662532

Submitted on 11 May 2022

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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Tectonics®

RESEARCH ARTICLE

10.1029/2021TC007089

**Key Points:**

* 4D dynamic analog subduction models investigate how end-member subduction modes affect the overriding plate deformation and topography

* Continuous trench retreat produces overall overriding plate extension, while long term trench advance promotes overall shortening

* Continuous trench retreat facilitates trenchward overriding plate tilting, while long term trench advance promotes landward tilting

**Correspondence to:**

K. Xue,
k.xue@vu.nl

**Citation:**

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*, e2021TC007089. https://doi.org/10.1029/2021TC007089

Received 27 SEP 2021
Accepted 17 JAN 2022

**Author Contributions:**

**Conceptualization:** K. Xue, W. P.

Schellart, V. Strak

**Data curation:** K. Xue

**Formal analysis:** K. Xue

**Funding acquisition:** W. P. Schellart

**Investigation:** K. Xue

**Methodology:** K. Xue, W. P. Schellart, V. Strak

**Project Administration:** W. P. Schellart

**Resources:** W. P. Schellart

**Supervision:** W. P. Schellart, V. Strak

**Validation:** K. Xue, W. P. Schellart, V. Strak

**Visualization:** K. Xue, V. Strak

**Writing – original draft:** K. Xue

**Writing – review & editing:** W. P. Schellart, V. Strak

© Wiley Periodicals LLC. The Authors.

This is an open access article under
the terms of the Creative Commons
Attribution License, which permits use,
distribution and reproduction in any
medium, provided the original work is
properly cited.

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

# Overriding Plate Deformation and Topography During Slab Rollback and Slab Rollover: Insights From Subduction Experiments

K. Xue¹ ID, W. P. Schellart¹ ID, and V. Strak¹,² ID

¹Department of Earth Sciences, Vrije Universiteit Amsterdam, Amsterdam, The Netherlands, ²Aix-Marseille University, CEREGE, CNRS, IRD, INRAE, Aix-en-Provence, France

**Abstract** Some subduction zones in nature show mainly overriding plate (OP) extension and low topography, and others show mainly shortening and elevated topography. Here we investigate how end-member subduction modes (trench retreat with slab rollback and trench advance with slab rollover) affect overriding plate deformation (OPD), topography, and mantle flow with time-evolving three-dimensional fully-dynamic analog models using particle image velocimetry. We conduct two sets of experiments, one of which is characterized by trench retreat, and the other characterized by trench advance. Experiments showing continuous trench retreat experience overall OP extension, while experiments dominated by trench advance experience overall shortening. Both subduction modes present fore-arc shortening and intra-arc extension. Our experiments indicate that the overall OPD is mainly driven by the horizontal mantle flow at the base of the OP inducing a viscous drag force ($F_D$), and is determined by the horizontal gradient of the horizontal mantle shear rate $d\dot{\gamma}/dx$, which controls the horizontal trench-normal gradient in $F_D$. Furthermore, a large-scale trenchward OP tilting and overall subsidence are observed in the experiments showing continuous trench retreat, while a landward OP tilting and an overall uplift are observed during long-term trench advance. The two types of topography during the two different subduction modes can be ascribed to the downward component of the large-scale trenchward mantle flow and the upward component of the landward mantle flow, respectively, and thus represent forms of dynamic topography. Our models showing trench advance provide a possible mechanism for OPD and topography at the Makran subduction zone.

## 1. Introduction

Subduction zones are generally considered to be the main driver of plate tectonics and mantle flow, and they also have a major impact on the Earth's surface topography. One way that subduction zones affect topography is by deforming the overriding plate, which can involve extension, causing subsidence, or shortening, causing uplift. Previous studies have investigated overriding plate deformation (OPD) using statistical parametric investigations (Heuret & Lallemand, 2005; Jarrard, 1986; Lallemand et al., 2008; Schellart, 2008b) and geodynamic models (Arcay et al., 2008; Boutelier & Cruden, 2013; Capitanio et al., 2010; Clark et al., 2008; Garel et al., 2014; Hertgen et al., 2020; Holt et al., 2015; Schellart & Moresi, 2013). Some of these works (Hertgen et al., 2020; Schellart & Moresi, 2013) concluded that the rate of overriding plate deformation is particularly dependent on the trench velocity.

Previous geodynamic models (Bellahsen et al., 2005; Di Giuseppe et al., 2008; Funiciello et al., 2008; Heuret et al., 2007; Ribe, 2010; Schellart, 2008a; Xue et al., 2020) and tomographic models (Van Der Voo et al., 1999; Widiyantoro et al., 1999; Wortel & Spakman, 2000) of subducting slabs have distinguished three main subduction styles as determined by the trench motion and the slab geometry: (1) continuous trench retreat with a slab rollback geometry (lazy “S” or “Z” slab geometry, Xue et al. [2020]; e.g., Calabria subduction zone, Scotia subduction zone, Tonga subduction zone), (2) long-term trench advance and slab rollover forming a “U” shaped slab geometry rotated 90° (e.g., India-Eurasia collision zone, Makran subduction zone), and (3) intermittent trench retreat and trench advance forming a steep folded slab pile (e.g., Mariana subduction zone). The first and the third subduction styles are relatively common in nature, and the patterns and mechanisms of OPD in these subduction styles have been investigated in buoyancy-driven subduction modeling works (Alsaif et al., 2020; Chen et al., 2015, 2016; Duarte et al., 2013; Hertgen et al., 2020; Holt et al., 2015; Yang et al., 2018). The second subduction style with slab rollover is less common in nature and the corresponding OPD has not been investigated
using buoyancy-driven models with a single subduction zone. However, there are at least three examples of subduction zones with a roll-over slab geometry in nature, that is, the Makran subduction zone (Amaru, 2007), the India-Eurasia continental subduction zone (Van Der Voo et al., 1999), and the subduction zone at the Caroline microplate (Fuji et al., 2021). With this contribution, we strive to fill the gap in knowledge that currently exists on slab rollover-overriding plate interaction. The aim of our study is to build buoyancy-driven models of progressive subduction in 3D space to investigate the two end-member subduction styles, slab rollover with trench advance and slab rollback with trench retreat, and to test the influence of subduction style on the variability of overriding plate deformation and topography. Previous experimental subduction studies investigated overriding plate deformation, but not topography, during slab rollover (Heuret et al., 2007) and overriding plate topography and deformation during trench retreat but not during trench advance (Martinod et al., 2013). Furthermore, these models used kinematic boundary conditions. Our models are buoyancy-driven and exclude external force or velocity boundary conditions, which can result in critically different model outcomes (Schellart & Strak, 2016).

In this study, we will focus on two research questions: (1) how do the deformation and topography of the overriding plate evolve during subduction with the rollover and rollback subduction styles? (2) what are the driving mechanisms that control OPD and topography during slab rollover and slab rollback? We present four dynamic upper-mantle subduction analog models with an overriding plate that evolve in three-dimensional space and which are driven only by the negative buoyancy of the slab, following an approach developed in earlier works (Chen et al., 2016, 2017; Duarte et al., 2013). Different subduction styles can be achieved by varying different geophysical parameters, for example, viscosity ratio of the subducting plate and ambient mantle material (Di Giuseppe et al., 2008; Funiciello et al., 2008; Garel et al., 2014; Ribe, 2010; Schellart, 2008a), plate thickness (Bellahsen et al., 2005; Ribe, 2010; Schellart, 2008a), plate length (Xue et al., 2020) and inclusion of lateral continental margins (Magni et al., 2014). In this study, we include two experimental sets that differ by their subducting plate thickness, to produce the two different subduction styles. We investigate the full subduction evolution, during which the slab geometry and dip of the slab change. We calculate the mantle flow velocity field in a vertical cross-section at the center of the subduction zone, as well as the OPD and topography using a particle image velocimetry (PIV) technique. Mantle flow, OPD, and the topographic evolution will be presented and described, and the forces responsible for the deformation and topography will be discussed.

## 2. Methodology

In this study, four experiments are conducted in an 80 by 60 cm transparent plexiglas tank filled with 8.25 cm of glucose syrup to model an upper mantle reservoir (Figure 1). Highly viscous layers of silicone putty mixed with iron powder are placed on top of the syrup simulating the overriding plate and subducting plate. Their dimensions are given in Figure 1. A weak coupling at the subduction zone interface is obtained using a mixture of paraffin oil (90wt%) and petrolatum (10wt%) as a lubricant (Duarte et al., 2013). The glucose syrup has a density of 1,408 kg/m³ and the overriding plate has a neutral buoyancy ($\rho_{OP} = 1,408$ kg/m³) relative to the sublithospheric mantle material. Because our experiments cannot model metamorphic reactions in nature, which increase the slab density during subduction, we use a relatively high density contrast (100 kg/m³) with high subducting plate density ($\rho_{SP} = 1,508$ kg/m³) following earlier works (Chen et al., 2015; Duarte et al., 2013). This is slightly higher than a density contrast of 80 kg/m³ between a mature subducting plate and sub-lithospheric mantle as determined by Cloos (1993) to account for the surface tension forces in the experiments that are negligible in nature (Schellart, 2008a). The viscosity of the glucose syrup is ~50 Pa s at 20°C and the viscosity ratio between the subducting plate and the glucose syrup is ~1300 in the experiments. We use a time scale of 1 s in the experiments representing 8300 years in nature, and a length scale of 1 cm representing 80 km in nature. To achieve dynamic similarity between model and nature, we assume that the slab sinking velocity follows the Stokes' settling law (Jacoby, 1973). With the dynamic scaling using the Stokes' settling law described in earlier works (Duarte et al., 2013; Strak & Schellart, 2016; Xue et al., 2020), we simulate the upper mantle with a viscosity that represents ~8.4 × 10¹⁹ Pa s in nature, which is comparable to estimates of 10¹⁹–10²¹ Pa s (Harig et al., 2010; Peltier, 2004), and a subducting oceanic plate with a viscosity representing 1.1 × 10²³ Pa s in nature.

In this study, two physical parameters are varied between experiments (Table 1): (1) the viscosity ratio of the overriding plate to the upper mantle ($\eta_{OP}/\eta_{UM}$) to consistently investigate OPD in a relatively weaker and stronger plate and (2) the thickness of the subducting plate ($T_{SP}$) to generate the two different subduction modes.
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Three-dimensional schematic model setup. The subducting plate thickness is 1.2 cm for Exps. 1 and 2, and 1.4 cm for Exps. 3 and 4. The laser produces a light sheet with a width of ~1 cm that is shot perpendicular to the trench axis along the mid line of the model. The experiments are recorded at 10 s intervals with one particle image velocimetry camera recording from the side to compute mantle flow, one normal camera recording from the side, and two stereoscopic cameras recording from the top to compute overriding plate deformation and topography.</figcaption></figure>

The experiments all have a free top surface, no-slip side walls, and a no-slip bottom boundary. The bottom boundary represents an impenetrable upper-lower mantle boundary. Both the overriding and subducting plates have two free lateral sides and a free trailing edge, simulating plates that are relatively mobile with strike-slip faults at the lateral sides and a mid-oceanic ridge at the trailing edge. To avoid the boundary effects from the side walls, the plates are placed in the middle of the tank and the subduction zone is far from the lateral side walls.

We scale our models to nature using density contrasts between the plates and the ambient mantle rather than using densities, and therefore we need to apply a topographic correction factor ($C_{Topo}$) when scaling the model topography to nature (Schellart & Strak, 2016). Since the density contrasts are the same in the experiments and in nature, the topographic correction factor is:

$$C_{Topo} = \frac{\rho_{m\_UM}}{\rho_{n\_UM}} \qquad (1)$$

Where $\rho_{m\_UM}$ and $\rho_{n\_UM}$ represent the density of the upper mantle in the models and nature, respectively. With this correction factor, we scale the topography in our model to nature as follows:

$$\frac{h_n}{h_m} = \frac{\rho_{m\_UM} \cdot l_n}{\rho_{n\_UM} \cdot l_m} \qquad (2)$$

Where $h_n$ and $h_m$ represent the topography in nature and the model, respectively, and $l_n/l_m = 8 \times 10^6$ as mentioned earlier. With an upper mantle density of 3250 kg/m³ in nature, we obtain a scaling ratio for topography between nature and experiments of $h_n/h_m = 3.47 \times 10^6$.

The experiments are initiated by pouring ~8 ml of syrup on top of the first ~3 cm of the subducting plate leading edge, forming an initial slab perturbation with a dip angle of ~30°. PIV was introduced to analog modeling in the 2000s (e.g., Adam et al., 2005; Hampel et al., 2004). Here we use

<table><caption>**Table 1**
The Difference in Model Setup Between Experiments</caption><thead><tr><th>Experiment number</th><th>η<sub>OP</sub>/η<sub>UM</sub></th><th>SP thickness (cm)</th><th>Slab geometry</th></tr></thead><tbody><tr><td>Exp. 1</td><td>~520</td><td>1.2</td><td>Slab rollover</td></tr><tr><td>Exp. 2</td><td>~1160</td><td>1.2</td><td>Slab rollover</td></tr><tr><td>Exp. 3</td><td>~520</td><td>1.4</td><td>Slab rollback</td></tr><tr><td>Exp. 4</td><td>~1160</td><td>1.4</td><td>Slab rollback</td></tr></tbody></table>

Note. $\eta_{OP}/\eta_{UM}$ = overriding plate to sublithospheric upper mantle viscosity ratio. SP = subducting plate.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Sequential side-view photographs of four experiments with a different subducting plate thickness ($T_{SP}$) and overriding plate to sublithospheric upper mantle viscosity ratio ($\eta_{OP}/\eta_{UM}$), showing the progressive subduction evolution. $\theta_B$: slab tip dip angle.</figcaption></figure>

two stereoscopic PIV cameras to record the top view of the experiments and to compute the strain and topography of the overriding plate (Figure 1). The stereoscopic PIV technique employed here is similar to Chen et al. (2016, 2017). White passive tracers are sprinkled on top of the overriding plate and are used in a cross-correlation technique to compute strain in 2D and topography using stereo photogrammetry with PIV cameras 1 and 2. We record the subduction process from the side with normal camera 4 while we compute the mantle flow using PIV camera 3. This camera visualizes phosphorescent polymer particles with a diameter of 20–50 µm that are mixed with the syrup homogeneously and that are illuminated by a laser sheet at 10 s intervals. The mantle flow velocity field is computed with a seeding density of ~40 particles/cm², a multi-pass interrogation window decreasing from 256 × 256 pixels to 128 × 128 pixels with an overlap of 75%. The overriding plate strain and topography are calculated with a subset size of 25 pixels and a step size of 6 pixels. Both the time interval (here 10 s) and the size of the interrogation window are chosen to reach an optimal signal-to-noise ratio. Pictures of a 3D textured calibration board aligned to the laser sheet position and to the overriding plate surface allow us to scale and correct the images used for the calculations.

# 3. Results

## 3.1. Subduction Kinematics

Our four experiments present two main subduction styles, where Exps. 1 and 2 with a thinner subducting plate ($T_{SP}$ = 1.2 cm) show trench retreat followed by trench advance with a rollover slab geometry (a U-shape on its side) and Exps. 3 and 4 with a thicker subducting plate ($T_{SP}$ = 1.4 cm) show continuous trench retreat during the whole experiment with a backward slab draping geometry (lazy S-shape; Figure 2).

Three main subduction phases are distinguished: the free sinking phase, the transitional phase, and the steady-state phase. We calculate the subducting plate ($v_{SP}$) and trench ($v_T$) velocities during subduction by tracking passive markers in successive photographs. All experiments experience a similar free sinking phase during which the slab subducts rapidly into the ambient upper mantle due to its negative buoyancy, with the trench retreating in the direction away from the overriding plate. During this phase, $v_{SP}$ (trenchward motion is positive), $v_T$ (retreat, i.e., oceanward motion, is positive) and subduction velocity ($v_S = v_{SP} + v_T$) increase rapidly to a maximum before the slab tip reaches the bottom boundary (Figures 3a–3c). Slab bending evolves differently between experiments with more significant slab bending in Exps. 1 and 2 (slab tip dip angle $\theta_B$ is ~83° for Exp. 1 and ~95° for Exp. 2 for a subduction depth of ~6.25 cm, Figures 2a and 2b) compared to Exps. 3 and 4 ($\theta_B$ is ~67° for Exp. 3 and ~64° for Exp. 4, Figures 2c and 2d). When the slab tip reaches the bottom boundary, the difference in slab tip dip angle ($\theta_B$) between the two sets of experiments is even more pronounced and, importantly, $\theta_B$ > 90° in Exps. 1-2
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Diagrams showing the evolution of (a) $v_{SP}$, (b) $v_T$, (c) $v_S$, (d) $v_{SP}/v_S$, (e) change in length of the overriding plate compared to the length at the time of subduction initiation, and (f) the average velocities and $v_{SP}/v_S$ during the steady-state phase for four experiments with a different $T_{SP}$ and $\eta_{OF}/\eta_{UM}$. Experiments 1 and 2 present trench retreat followed by trench advance, and Exps. 3 and 4 show continuous trench retreat. RO: rollover. RB: rollback. The vertical arrows with different colors in (a) point to the beginning (solid) and the end (dotted) of the transitional phase of the corresponding models. For the lines in (e), a positive and negative slope indicates (on average) overriding plate extension and shortening, respectively.</figcaption></figure>

and $\theta_B < 90^\circ$ in Exps. 3-4. The slab tip interacts with the bottom boundary during a short period and all velocities drop quickly to a minimum, which we define as the transitional phase (the beginning and the end are indicated with solid and dotted arrows, respectively, in Figure 3a). Following this, the experiments reach a steady-state phase with $v_{SP}$, $v_T$, and $v_S$ reaching relatively stable values (Figures 3a-3c), while $v_{SP}/v_S$ remains almost constant (Figure 3d). During the steady-state phase, Exps. 1 and 2 show trench advance and a rollover slab geometry (last panels in Figures 2a and 2b), whereas Exps. 3 and 4 continue with trench retreat and slab rollback (last panels in Figures 2c and 2d).

The velocities $v_{SP}$, $v_T$ and $v_S$ evolve similarly during the free sinking phase for all experiments but are very different during the steady-state subduction phase for the two different subduction modes. Indeed, the trench velocity ($v_T$, positive toward the subducting plate) is negative for Exps. 1 and 2 showing trench advance and positive for Exps. 3 and 4 show trench retreat, and $v_{SP}$ is much higher in Exps. 1 and 2 compared to that in Exps. 3 and 4 (Figures 3a and 3f). Furthermore, the experiments with a lower $\eta_{OF}/\eta_{UM}$ have slightly lower $v_{SP}$, compared to the experiments with higher $\eta_{OF}/\eta_{UM}$ but with the same subduction mode (Figure 3f). $v_T$ is similar for Exps. 1 and 2,
with an average of -0.028 mm/s for both experiments (Figures 3b and 3f). On the other hand, Exp. 4 with higher $\eta_{OP}/\eta_{UM}$ shows a higher $v_T$ (0.089 mm/s) compared to Exp. 3 (0.066 mm/s).

During both free sinking and steady-state phases, $v_S$ is higher for Exps. 3–4 showing trench retreat compared to that in Exps. 1-2 with trench advance (Figures 3c and 3f). In addition, $v_S$ in Exp. 4 with a higher $\eta_{OP}/\eta_{UM}$ is higher than that in Exp. 3, while it is only slightly higher in Exp. 2 than that in Exp. 1.

A subduction partitioning ratio ($v_{SP}/v_S$) that exceeds 1 indicates trench advance, while $v_{SP}/v_S < 1$ indicates trench retreat. The subduction partitioning remains below 1 during the free sinking phase for all experiments, indicating that all experiments experience trench retreat during this phase. During the steady-state phase, $v_{SP}/v_S$ remains stable at ~0.2 for Exps. 3 and 4, while it is above 1 for Exps. 1 (~1.9) and 2 (~1.7).

## 3.2. Deformation of the Overriding Plate

The pattern and degree of OPD vary between experiments with different subduction styles and $\eta_{OP}/\eta_{UM}$ (Figures 3e and 4). The overriding plate length in Exps. 1 and 2 increases during the free sinking and transition phases and decreases during the steady-state phase, while in Exps. 3 and 4 it increases continuously (Figure 3e). All experiments experience a similar pattern of OPD during the free sinking phase (Figures 4a, 4b, 4g and 4h), with overall overriding plate extension in the trench normal direction and only shortening in a narrow region located within ~20 mm from the trench, which we define as the fore-arc region (Figure 4a). In addition, all experiments show a localized area with higher extensional strain, which we define as the intra-arc domain for this study as shown in Figure 4a, at a distance between ~20 and ~45 mm from the trench. Note that we define the area between the intra-arc and the trailing edge of the overriding plate as the back-arc region (Figure 4a). Furthermore, strains are higher for experiments with a lower $\eta_{OP}/\eta_{UM}$ (Exps. 1 and 3) than with a higher $\eta_{OP}/\eta_{UM}$ (Exps. 2 and 4).

After the free sinking phase, fore-arc shortening continues and is more pronounced for experiments showing trench advance ($\epsilon_{XX} \sim -46\%$ (Exp. 1), $\sim -40\%$ (Exp. 2)) compared to experiments showing trench retreat ($\epsilon_{XX} \sim -28\%$ (Exp. 3), $\sim -10\%$ (Exp. 4)). The localized extensional intra-arc area at ~20–45 mm from the trench experiences stronger extension and becomes more clearly defined during the steady-state phase for all experiments (Figures 4c, 4d, 4i and 4j). In addition, there is also a narrow area at the overriding plate trailing edge showing extensional deformation for all experiments. Apart from these local areas that show a common strain style (shortening or extension) in all experiments, the overall style of deformation in the back-arc during the steady-state phase differs for the two different subduction styles. Experiments 1 and 2 with the slab rollover structure experience extensive trench-normal shortening (Figures 4c and 4d), whereas Exps. 3 and 4 with the backward slab draping geometry experience continuous trench-normal extension (Figures 4i and 4j).

The total overriding plate strain is presented in Figures 4e, 4f, 4k and 4l, which is equivalent to the sum of the strain during the free sinking phase, the transitional phase, and the steady-state phase. Considering that the deformation style for Exps. 3 and 4 do not change during progressive subduction evolution, the total finite strain for these experiments shows pronounced fore-arc shortening, intra-arc extension, and back-arc extension. For Exps. 1 and 2, the deformation style in the fore-arc and intra-arc domain remain unchanged, but in the back-arc domain the extensional deformation changes to shortening, and the finite shortening deformation summed over the entire subduction evolution is smaller than that generated during the steady-state phase.

## 3.3. Subduction Induced Mantle Flow and Its Correlation With OPD

The cross section of the mantle flow velocity field and the gradient of the horizontal trench-normal mantle flow velocity ($dv_x/dx$) are presented at two moments for all experiments (lower panels in Figure 5): (1) when the slab tip reaches a depth of 6.25 cm (~2 cm above the bottom boundary) during the free sinking phase (Figures 5a, 5b, 5e and 5f) and (2) half way during the steady-state phase (Figures 5c, 5d, 5g and 5h). The velocity fields and $dv_x/dx$ are plotted with the corresponding overriding plate strain (upper panels in Figures 5a–5h) calculated between that particular amount of subduction and the beginning of the corresponding subduction phase to allow for a direct comparison. We note that it is actually the horizontal gradient of the trench-normal basal shear rate ($d\dot{\gamma}/dx$) multiplied with the sublithospheric mantle viscosity ($\eta_{UM}$, which is constant in our experiments) that determines the actual trench-normal horizontal normal stress ($\sigma_{XX}$), normal strain rate ($\epsilon_{XX}$) and strain ($\epsilon_{XX}$) in the overriding plate. However, considering that $dv_x/dx$ in the overriding plate is one to three orders of magnitude
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Finite strain maps of the overriding plate in the trench normal direction ($\epsilon_{xx}$) at the end of the free sinking phase (a, b, g, and h), for the period that includes the transitional and entire steady-state phases (c, d, i, and j), and at the end of the subduction experiments (entire duration of the experiments) (e, f, k, and l) for four experiments with a different $T_{sp}$ and $\eta_{OF}/\eta_{UM}$. A division between fore-arc, intra-arc, and back-arc regions is illustrated in (a). For each phase, overriding plate deformation is computed for an equivalent amount of subduction.</figcaption></figure>

smaller than $dv_x/dx$ in the underlying mantle, the sign and magnitude of $dv_x/dx$ directly correspond to those of $d\dot{\gamma}/dx$. Nevertheless, we have plotted $d\dot{\gamma}/dx$ for the uppermost ~10 mm of sublithospheric mantle below the base of the overriding plate in the upper panels of Figures 5a–5h.

During the free sinking phase (Figures 5a, 5b, 5e and 5f), the upper mantle material below the overriding plate mainly flows toward the trench, and only in a relatively small area near the slab tip, it directs away from the trench. During this phase, $dv_x/dx$ and $d\dot{\gamma}/dx$ are mostly positive, which indicates that $\dot{\gamma}$ increases toward the trench in most parts of the upper mantle. Relatively high $dv_x/dx$ and $d\dot{\gamma}/dx$ occur at ~50–90 mm from the
<figure><img src="image_8.png" /></figure>

**Figure 5.** Cross-sections of mantle flow velocity field shown with black vectors and $dv_x/dx$ shown with blue-red colors (lower panel in each subfigure), as well as the trench-normal horizontal overriding plate normal strain ($\epsilon_{XX}$, black line in upper panel of each subfigure) along the mid-line of the overriding plate. The trench-normal horizontal component of the mantle flow velocity ($v_x$) is also shown in the upper panels with orange lines (~10 mm below the plate base as indicated with an orange arrow in the lower panel of (a)) and purple lines (~40 mm below the plate base as indicated with a purple arrow in the lower panel of (a)), and $v_T$ at the corresponding moment is shown as a blue dot in the upper panel. The trench-normal horizontal gradient of the shear rate ($d\dot{\gamma}/dx$) for a ~10 mm thick zone below the base of the overriding plate is indicated with a green line in the upper panels of (a)–(h). The overriding plate strain is computed between the moment shown in each corresponding lower panel and the start of the free sinking phase for a, b, e and f, and the start of the steady-state phase for c, d, g and h. In the upper panels showing OP strain, the shaded areas in blue represent shortening and red represent extension. RO: rollover. RB: rollback. Note that an enlargement of the velocity field as shown in panels (c) and (g) is shown in Figure 6.
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Cross-sections through the center of the subduction zone, showing the mantle velocity field (vectors) and vertical mantle velocity component (gray-white-green color scheme) for (a) the trench advance subduction mode (Exp. 1) and (b) the trench retreat subduction mode (Exp. 3). Note that positive $v_z$ (green) indicates upwelling mantle flow and negative $v_z$ (gray) indicates downwelling mantle flow.</figcaption></figure>

trench, generally correlating with higher magnitudes of the overriding plate extension at ~25–100 mm from the trench. The extensional overriding plate strain generally coincides with the distribution of positive $dv_x/dx$ and $d\dot{\gamma}/dx$ right below the overriding plate base.

During the steady-state phase, the mantle flow below the overriding plate in the experiments showing trench advance directs away from the trench ($v_x < 0$), and $|v_x|$ right below the overriding plate gradually increases from the overriding plate trailing edge toward a distance of 50–60 mm from the trench (so $d\dot{\gamma}/dx < 0$) but $|v_x|$ gradually decreases then up to the trench (so $d\dot{\gamma}/dx > 0$, Figures 5c and 5d). In spatial accordance with the distribution of $dv_x/dx$ and $d\dot{\gamma}/dx$, the overriding plate presents negative strain (shortening) from near the overriding plate trailing edge (~160–180 mm from the trench) to ~50 mm from the trench (strain reaches ~-4% For Exp. 1 and ~-2% for Exp. 2) and shows positive strain (extension) between ~50 and ~20 mm from the trench with maximum extensional strains of ~12% for Exp. 1 and ~4% for Exp. 2. In the experiments showing trench retreat during the steady-state phase (Exps. 3 and 4), the mantle flow below the overriding plate is directed toward the trench, and $\epsilon_{XX}$, $dv_x/dx$ and $d\dot{\gamma}/dx$ are generally positive (Figures 5g, 5h and 6b). The magnitude of $dv_x/dx$ and $d\dot{\gamma}/dx$ are especially high at ~30–120 mm from the trench.

### 3.4. Evolution of the Overriding Plate Topography

The overriding plate shows large-scale tilting for all experiments. During the free sinking phase all experiments show an overriding plate with a downward slope toward the trench with an average tilting angle of ~0.04°. The angle is measured by plotting for each experiment the best-fit line of the scaled topography for the back-arc region that is confined by the vertical lines shown in Figure 7. During the steady-state phase, however, the overriding plate tilts in different directions depending on the subduction mode. For the slab rollover mode (Exps. 1 and 2), a large-scale uplift is observed in the back-arc except close to the trailing edge, resulting in an upward slope toward the trench with an average tilting angle of ~0.06° (orange lines in Figures 7a and 7b). In contrast, for the slab rollback mode (Exps. 3 and 4), the entire overriding plate is dragged down and tilts toward the trench side (orange lines in Figures 7c and 7d), with an average tilting angle of ~0.1° at the end of the experiment. We also quantify the vertical force ($F_{Tilt}$) that generates the large-scale tilting (uplift for Exps. 1 and 2 and subsidence for Exps. 3 and 4) of the overriding plate, which equals the product of the total volume of the depression or uplifted region (with respect to the reference level (topography = 0)) due to tilting ($V_{Tilt}$), the density contrast between the upper mantle ($\rho_{UM} = 1,408$ kg/m$^3$) and air ($\rho_{Air} = 1.225$ kg/m$^3$; $\Delta\rho_{UM-Air}$), and the gravitational acceleration (g), giving $F_{Tilt} = V_{Tilt}\Delta\rho_{UM-Air}g$. Note that the subsidence is defined in this study as the long-wavelength subsidence signal (which is in the back-arc region), not the short-wavelength subsidence signal (which we refer to as depression and is in the fore-arc region). This gives a vertical force used to produce the overriding plate tilting of approximately –0.029, –0.047, 0.03, and 0.044 N for Exps. 1-4, respectively (note that negative numbers indicate a vertical compressive force, causing uplift, and positive numbers indicate a vertical tensile force, causing subsidence). We can compare this tilting force with the only driving force present in our experiments, the slab negative buoyancy force $F_{BU}$, which is simply $V_{slab}\Delta\rho_{SLAB-UM}g$, where $V_{slab}$ is the volume of the slab that is sinking in the upper mantle (i.e., the inclined part) and $\Delta\rho_{SLAB-UM}$ is the density contrast between the slab and ambient mantle. Then we get $|F_{Tilt}|/|F_{BU}|$ is 0.28, 0.40, 0.21, 0.29 (+/-30%) for Exps. 1-4, respectively.

From the beginning to the end of the free sinking phase (Figure 7), a depression forms at ~15–20 mm from the trench for all experiments, with a depth ranging between ~0.5 (Exp. 4) and ~1.1 mm (Exp. 1) relative to the
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Diagram showing the overriding plate topography for four experiments at three moments: early free sinking phase (black line), the end of the free sinking phase (aqua line), and the end of the steady-state phase (orange line). The vertical arrows indicate the location of the fore-arc depression for the curves with the corresponding color. The shaded area in magenta determines the area that is above the subduction zone plate boundary interface at the end of the steady-state phase. Averaged tilting is measured by fitting a best-fit line to the observed topography between two vertical dashed lines in each phase with the corresponding color. The x-position in the diagrams has been centered on the trench location for better comparison between the individual curves.</figcaption></figure>

edge on the back-arc side of the local depression. During the steady-state phase, this fore-arc depression becomes deeper for Exps. 3 and 4, and shallower for Exps. 1 and 2. At the end of the entire subduction process, the fore-arc depression for all experiments reaches a depth ranging between ~0.8 (Exps. 1 and 2) and ~1 mm (Exps. 3 and 4), with respect to the edge on the back-arc side of the local depression. The fore-arc depression is often flanked on the trench side by a fore-arc bulge, up to ~0.7 mm high with respect to the base of the fore-arc depression, which is generally most pronounced during the free sinking phase.

# 4. Discussion

## 4.1. Subduction Kinematics and Slab Geometries

The different subduction styles that are observed due to the different $T_{SP}$ can be explained by the increase in slab negative buoyancy force with increasing $T_{SP}$, which is the only driving force in the subduction system. The increased driving force promotes faster subduction and particularly promotes faster $v_T$ during the free sinking phase, while $v_{SP}$ is less affected (Figure 3b). Earlier works (Griffiths et al., 1995; Xue et al., 2020) have shown that an increase in $v_T$ particularly causes a decrease in slab dip angle. In addition, the bending resistance at the subduction zone hinge scales with the cube of $T_{SP}$ (Conrad & Hager, 1999). Therefore, the slab develops a larger bending radius and a smaller bending angle at the slab tip for a larger $T_{SP}$ during the free sinking phase (Figure 2), in agreement with earlier subduction modeling works (Capitanio et al., 2007; Irvine & Schellart, 2012). The larger bending radius, together with a lower slab dip angle due to a faster trench retreat, promotes a lower slab tip dip angle (~75° and 70° for Exps. 3 and 4, respectively) when the slab tip reaches the bottom boundary, thereby producing trench retreat with a rollback slab geometry (Figures 2c and 2d). On the other hand, a thinner subducting plate will develop a smaller bending radius at the subduction hinge and, together with a lower trench velocity, this will cause a higher slab tip dip angle (~109° and 108° for Exps. 1 and 2, respectively) when the slab tip reaches the bottom boundary (Figures 2a and 2b), facilitating slab rollover and trench advance. Furthermore, Lallemand et al. (2008) showed that trench advance with a rollover slab geometry is associated with a higher $v_{SP}$, whereas trench retreat with a rollback slab geometry corresponds with a relatively lower $v_{SP}$, which is consistent with the experimental results reported in this study.

The effect of plate thickness on the subduction mode in our study is opposite to trends reported in some earlier modeling studies (e.g., Di Giuseppe et al., 2008; Faccenna et al., 2009), which found that a thicker subducting plate promotes trench advance and a thinner plate promotes trench retreat. The apparent conflict can be explained
<figure><img src="image_11.png" /><figcaption>**Figure 8.** (a, b) Illustration of the forces and the velocity directions involved in two different subduction modes. (c, d) Decomposition of forces at the subduction zone interface in two different subduction models. Subduction interface shear drag ($F_{SH}$), trench suction ($F_{SU}$), and basal shear drag ($F_D$) are shown by purple, green, and blue arrows, respectively. Note that $d\dot{y}/dx$ is the trench-normal horizontal gradient of the horizontal shear rate in the sublithospheric mantle, while $v_T$ is the trench velocity. Abbreviations: IAE: Intra-arc extension; FAS: Fore-arc shortening.</figcaption></figure>

when we consider previously published regime diagrams showing the dependence of the subduction style on the mantle depth/slab thickness ratio and slab/mantle viscosity ratio (Li & Ribe, 2012; Ribe, 2010; Schellart, 2008a). For one particular viscosity ratio, one might observe, when increasing the slab thickness (and thus decreasing the mantle depth/slab thickness ratio), that the subduction style changes from rollover to rollback, while for another viscosity ratio, one might observe, when increasing the slab thickness, that the subduction style changes from rollback to rollover (see e.g., Figure 13 in Schellart [2008a] or Figure 11 in Ribe [2010]). We add to this the complexity that the regime diagrams and the studies from Di Giuseppe et al. (2008) and Faccenna et al. (2009) were based on subduction models without an overriding plate, while the current models include an overriding plate, and the complexity that plate length also influences the subduction style (Xue et al., 2020). Indeed, earlier studies have shown that a variety of physical parameters can affect the subduction mode (e.g., Bellahsen et al., 2005; Di Giuseppe et al., 2008; Funiciello et al., 2008; Garel et al., 2014; Magni et al., 2014; Schellart, 2008a; Xue et al., 2020). A detailed discussion on how the different parameters affect the subduction mode, however, is beyond the scope of this study, as we focus our investigation on how the two end-member subduction modes (rollback and rollover) affect OPD and topography.

Another parameter that is varied in this study is $\eta_{OP}/\eta_{UM}$. The experiments with a lower $\eta_{OP}/\eta_{UM}$ (~520) experience almost double maximum overriding plate strain (Figure 4) compared to that with a higher $\eta_{OP}/\eta_{UM}$ (~1160), which can be explained by the linear viscous rheology of the overriding plate and the relationship where stress scales with the product of the effective viscosity and strain rate. Hertgen et al. (2020) found that with a stronger overriding plate, the back-arc basin forms later. In addition, we find that $\eta_{OP}/\eta_{UM}$ has no impact on the slab geometry and the subduction mode, in agreement with Holt et al. (2015), who showed that overriding plate viscosity affects overriding plate stress state, but has no obvious impact on the slab morphology.

## 4.2. Overriding Plate Deformation and Topography Evolution

Some of the main forces in plate tectonics and mantle dynamics include the slab negative buoyancy force, the trench suction force, basal drag forces, and plate boundary forces (Elsasser, 1971; Forsyth & Uyeda, 1975; Shemenda, 1993). The forces responsible for deforming the overriding plate are illustrated in Figure 8 and include: (1) the viscous drag force ($F_D$) at the base of the overriding plate, which is related to the basal shear stress, which is $\eta_{mantle}$ ($dv/dz$). The horizontal gradient in basal shear stress, $\eta_{mantle}$ ($d\dot{y}/dx$), is what causes extension or
shortening in the overriding plate, with a positive trench-directed gradient (i.e., increasing trenchward) causing an extension, and a negative trench-directed gradient (i.e., decreasing trenchward) causing shortening. (2) The shear force ($F_{SH}$) at the subduction zone interface, which is the product of the subduction interface shear rate, the subduction zone plate boundary interface surface area, and the subduction channel effective viscosity, where the shear rate is equal to the subduction rate divided by the channel thickness. And (3) the suction force ($F_{SU}$) between the plates which is oriented perpendicular to the subduction zone interface and points downward during trench retreat and upward during trench advance.

During trench retreat, the mantle flow right below the overriding plate is directed toward the subducting plate, and the magnitude of $lv_x$ and the trenchward basal viscous drag $F_D$ increase toward the mantle wedge corner with $d\dot{\gamma}/dx > 0$ (Figures 5g and 5h). The increased $F_D$ toward the trench, together with the trench-retreat-induced trench suction ($F_{SU}$), exerts trench-normal deviatoric tension in the overriding plate, causing extension (Figures 5a, 5b, 5e, 5f, 5g, and 5h), except for the fore-arc region where subduction interface shear drag ($F_{SH}$) and basal drag ($F_D$) have an opposite shear sense, causing overriding plate shortening (Figure 8b). On the other hand, in the trench advance subduction mode, the mantle flow below the overriding plate is oriented away from the subducting plate. Furthermore, $v_x$ and $F_D$ generally become more negative from the overriding plate trailing edge toward the trench with $d\dot{\gamma}/dx < 0$ before reaching the intra-arc region and then become less negative from the intra-arc region toward the trench with $d\dot{\gamma}/dx > 0$ (Figure 8a). The negative $d\dot{\gamma}/dx$ and the associated landward-decreasing $F_D$, combined with the deviatoric compressive stresses above the subduction zone interface that are caused by the trench advance and slab rollover (inducing negative trench suction), result in shortening in the area between the overriding plate trailing edge and the intra-arc edge, as well as above the subduction interface, while the local positive $d\dot{\gamma}/dx$ in the intra-arc region and the landward-increasing $F_D$ facilitate overriding plate extension in the intra-arc region.

Earlier modeling studies have investigated OPD in a trench retreat mode (Arcay et al., 2008; Chen et al., 2015, 2016, 2017; Duarte et al., 2013; Guillaume et al., 2013; Heuret et al., 2007; Holt et al., 2015; Martinod et al., 2013; Schellart & Moresi, 2013). From these studies, those that presented models using a buoyancy-driven approach showed similar OPD patterns as here during trench retreat and slab rollback, with overall overriding plate extension and mostly local fore-arc shortening (Chen et al., 2015, 2016; Duarte et al., 2013; Schellart & Moresi, 2013). Our results agree with earlier conceptual models and observations, which suggest that the subduction mode generally determines the style of OPD, with overriding plate extension correlating with trench retreat (Elsasser, 1971; Le Pichon et al., 1982; Lonergan & White, 1997; Malinverno & Ryan, 1986; Molnar & Atwater, 1978; Rosenbaum & Lister, 2004; Schellart, 2008b) and shortening correlating with trench advance or slow trench motion (Schellart, 2008b). We note that the boundary condition at the trailing edge of the overriding plate (free or fixed) also has an effect on the motion and deformation style of the overriding plate. For example, models from Chen et al. (2015) showed that an overriding plate will experience significantly larger extensional strains when it is fixed at the trailing edge compared to when it is free. With our experiments, we show that the different subduction modes (slab rollback and slab rollover) cause different mantle flow patterns, and thereby different mantle velocity gradients and gradients in basal drag force, as well as different trench suction forces. The combination of these forces determines if the overriding plate experiences shortening or extension.

Our experimental models are simplified in that they are isothermal. Including thermal gradients would likely increase mantle flow rates below the overriding plate, but at the same time decrease the effective viscosity in the mantle below the overriding plate. The combined effect of increased velocity and decreased viscosity on basal drag is likely small, as they likely cancel out each other, such that the basal drag force is comparable in magnitude as for an isothermal model set-up. Another limitation of our models is that the viscous overriding plate material slowly spreads laterally because of the unwanted side effect of surface tension. The amount of spreading depends on $\eta_{OP}/\eta_{UM}$ and the duration of the experiments, such that lower $\eta_{OP}/\eta_{UM}$ and longer experimental duration promote larger spreading. The lateral spreading is most evident in the trench-parallel direction (y-direction). Toward the end of each experimental run, the increase in surface area is ~10.4%, 2.4%, 3.5%, and 1.0% for Exps. 1 to 4, respectively, which is consistent with the lower $\eta_{OP}/\eta_{UM}$ in Exps. 1 and 3 and the relatively longer experimental duration of Exps. 1 and 2.

### 4.2.1. Fore-arc Shortening

We observe fore-arc shortening (Figures 4 and 5) and a fore-arc depression (Figure 7) during the whole subduction process for both subduction styles. During trench retreat, both the horizontal component of $F_{SH}$ ($F_{SH\_HOR}$) and
$F_D$ ($F_{D\_HOR}$) work together to promote fore-arc shortening due to their opposite shear drag sense that enhances convergence, while the horizontal component of $F_{SU}$ ($F_{SU\_HOR}$) enhances extension (Figure 8b). To produce fore-arc shortening in this subduction mode, then $|F_{SU\_HOR}| + |F_{D\_HOR}| > |F_{D\_HOR}|$. In the trench advance subduction mode, $F_{SH}$ and $F_{SU}$ promote fore-arc shortening and $F_{D\_HOR}$ promotes fore-arc extension as its drag is directed away from the fore-arc domain. Thus, to produce shortening during trench advance, then $|F_{D\_HOR}| < |F_{SU\_HOR}| + |F_{SU\_HOR}|$. $F_{SU\_HOR}$ depends on the resistance to translate the overriding plate at the trench and correlates with the magnitude of $v_T$ (Chen et al., 2015; Shemenda, 1993), while $F_D$ depends on the shear rate at the base of the overriding plate and effective viscosity of the sublithospheric mantle below the base of the overriding plate. In addition, $F_{SU\_HOR}$ always provides compression to the fore-arc region in both subduction modes, and its magnitude depends on the effective viscosity, shear rate, and the contact surface area between the plates, which is larger (~21%) during long-term trench advance, compared to that during long-term trench retreat.

During the free sinking phase, the fore-arc region experiences shortening and the trench retreats in all experiments, which means $|F_{SU\_HOR}| + |F_{D\_HOR}| > |F_{D\_HOR}|$. During this phase, the fore-arc shortening strain is lower in Exps. 3 and 4 with higher $T_{SP}$ than that in Exps. 1 and 2. A higher $T_{SP}$ produces a higher $v_S$ (Figure 3c) resulting in a higher $F_{SH}$ and a higher trenchward $v_x$ resulting in a higher $F_D$, but it also promotes a higher tensional $F_{SU}$ due to the associated higher $v_T$ (Figure 3b). The experiments indicate that $F_{SU}$ increases more than the other two forces, such that the overall horizontal deviatoric compressive stress in the fore-arc is lower with higher $T_{SP}$, producing a lower trench-normal shortening strain.

During the steady-state phase, the fore-arc shortening is also less pronounced in the trench retreat mode (Exps. 3 and 4) than that in the trench advance mode (Exps. 1 and 2). In the trench advance mode, $v_S$ is lower producing a lower $F_{SH}$, but at the same time the subduction interface is longer (interface length = ~1.9 and 1.8 cm for Exps. 1 and 2, respectively, compared to ~1.6 and 1.3 cm for Exps. 3 and 4, respectively) and the subduction channel is thinner because the lubricant is partly squeezed out from the interface (as observed in the experiments), promoting a higher $F_{SH}$, which together with a compressional $F_{SU}$ due to the trench advance, promotes shortening. In the trench retreat mode, the combination of the three forces remains the same as that during the free sinking phase. However, the combination of a relatively higher $F_{SH}$ and $F_D$ promoting compression and a tensional $F_{SU}$ due to trench retreat, results in less fore-arc shortening than during trench advance.

### 4.2.2. Intra-arc and Back-arc Extension

Both subduction modes produce a confined intra-arc and near-back-arc region of maximum extensional deformation throughout the experiments, and at almost the same distance to the trench of approximately half of the subduction zone width. This region generally coincides spatially with the distribution of $d\dot{\gamma}/dx > 0$ for both subduction modes (Figure 5), and so we propose that $F_D$ is the main factor driving such extension. Some earlier studies suggest that $F_{SU}$ is the main factor driving such extension (Elsasser, 1971; Shemenda, 1993). However, our experiments demonstrate that $F_{SU\_HOR}$ is not even enough to compete with $F_{D\_HOR}$ and $F_{SH\_HOR}$ near the trench to drive continuous extension in the fore-arc region during trench retreat, while it induces compression in the fore-arc region during trench advance. Some recent models also demonstrate that subduction-induced mantle return flow is the dominant driver of intra-arc and back-arc extension (Chen et al., 2016; Holt et al., 2015; Schellart & Moresi, 2013). In addition, models from Chen et al. (2015) indicate that back-arc extension reaches a maximum at a distance to the trench of approximately half of the subduction zone width, which is comparable to the toroidal mantle return flow radius. Our results agree with these earlier studies and imply that $F_D$ resulting from mantle return flow with $d\dot{\gamma}/dx > 0$ (Figure 8) is the dominant controlling factor for intra-arc and back-arc extension.

### 4.2.3. Long-Wavelength Dynamic Topography and Overriding Plate Fore-arc Topography

Our experiments demonstrate an exciting new observation that concerns the long-wavelength (1,000–2,000 km) topographic evolution of the overriding plate, involving whole plate tilting, the direction of which depends on the style of subduction (slab rollback vs. slab rollover). The entire overriding plate is dragged down and tilts toward the trench in the trench retreat mode (~0.04–0.10° averaged over ~1,500 km), while it is pushed up with a maximum upward deflection at ~2 cm (scaling to ~160 km) from the trench and tilts toward the trailing edge in the trench advance mode (~0.04–0.06° averaged over ~1,500 km). Among the forces, the vertical component of $F_{SH}$ ($F_{SH\_VER}$) is in the same direction (downward) for both subduction modes, and only the vertical component of $F_{SU}$ ($F_{SU\_VER}$) is upward for the trench advance mode and downward for the trench retreat mode. If it is $F_{SU}$ that drives the overriding plate tilt, one would expect the maximum uplift during trench advance to be above the subduction
zone interface. However, the model results show that we have a fore-arc depression in Exps. 1 and 2 (Figure 7). In addition, we expect that $F_{SU\_VER}$ only affects a small region (the fore-arc region of a few cm), as it is applied to a small area (the subduction zone interface), but not the entire plate of ~20 cm. On the other hand, the mantle flow during rollover shows an upward velocity component below the overriding plate (Figures 5c, 5d and 6a), while the mantle flow during rollback shows a downward velocity component (Figures 5g, 5h and 6b). Therefore, we suggest that this large-scale plate tilting in both subduction modes is related to the upper mantle flow (the vertical velocity component thereof), and thus can be considered as dynamic topography. The force generating the large-scale overriding plate tilting ($F_{Till}$) accounts for a considerable portion of the slab negative buoyancy force ($F_{BU}$) with $|F_{Till}|/|F_{BU}| = 0.2–0.4$, with an error margin of ~30%. The observation of large-scale trench-ward overriding plate tilting in the trench retreat subduction mode has also been observed in numerical models from Crameri and Lithgow-Bertelloni (2018), who proposed it is caused by the large scale mantle flow when the slab tip reached the more viscous lower mantle and excited a larger mantle flow cell. Their models showed overriding plate tilting ranging between 0.01° and 0.07° over a distance exceeding 1,000 km, which generally agrees with our models (tilting of ~0.04°–0.1° over ~1,500 km).

In the work by Husson (2006) the dynamic topography from theoretical models and for several narrow retreating subduction zones is presented. It shows that subsidence is greatest at a distance from the trench of approximately half of the subduction zone width, reaching depths of the order ~1,500 m. Furthermore, the extent to which overriding plate subsidence is observed extends in a trench-normal direction over a distance of approximately twice the subduction zone width. The dynamic topography observed in our models with a trench retreat mode generally agrees with Husson (2006), with a broad local maximum in subsidence at a distance of ~0.3–0.5 times the subduction zone width of 2–2.5 km, and a comparable lateral extent of subsidence of ~2 times the subduction zone width. Models from Husson et al. (2012) characterized by slab rollback show subsidence in the region where one would have the overriding plate, which is not present in their models, but there is no systematic tilting toward the trench and the subsidence is smaller (maximum of ~0.1 mm compared to our models with 0.6–1.0 mm) at a distance that is ~half the width of the subduction zone. We attribute the lower subsidence to the absence of an overriding plate and the presence of low-viscosity (asthenospheric) material at the surface in the models of Husson et al. (2012).

A fore-arc depression is observed in all experiments with different subduction styles, with a maximum depression located above the subduction interface. This has only been observed in subduction models characterized by trench retreat from earlier studies (Chen et al., 2017; Hassani et al., 1997). Hassani et al. (1997) have proposed that the suction force between the plates is the main driver of the fore-arc depression. Models from Chen et al. (2017) with a similar setup and physical parameters as in our study have shown that $F_{SU}$ dominantly drives the formation of the fore-arc depression, while $F_{SH}$ plays a minor role because of the low $v_s$, and mantle flow does not produce the fore-arc depression because the maximum vertical mantle flow velocity does not coincide with the deepest point of the fore-arc depression. Chen et al. (2017) also found that as the slab dip angle gradually increases during the free sinking phase, $F_{SU\_VER}$ increases, and thus the maximum fore-arc depression increases. Our experimental results show that, during the free sinking phase, Exps. 1 and 2 have a higher maximum fore-arc depression compared to Exps. 3 and 4. This is because of the slab dip angle in Exps. 1 and 2 increases more during the free sinking phase, resulting in a higher $F_{SU\_VER}$. This further confirms earlier findings that $F_{SU}$ is the main driver of the fore-arc depression during the free-sinking phase. The depth of the fore-arc basin is ~0.6–1 mm (~2.1–3.5 km in nature) in the experiments showing continuous trench retreat, which is comparable with that from Chen et al. (2017) (~1.4–4 km in nature). Subduction models from Cerpa and Arcay (2020) show a similar topographic evolution in the fore-arc region, compared to our models. Their models show a bulge next to the trench, which is bordered by a local depression, in agreement with our model results, and the fore-arc bulge and depression are the most pronounced during the free sinking phase, gradually decreasing over time after the free sinking phase, which also agrees with our model results. The free subduction model (no imposed velocity boundary conditions) with a subduction interface friction coefficient of 0.04 in Cerpa and Arcay (2020), which is geodynamically most comparable to our buoyancy-driven subduction models, also shows the most comparable fore-arc topography. Their models show a bulge topography that lies lower than the trailing part of the overriding plate, as in our Exps. 3 and 4 (Figures 7c and 7d), and a fore-arc depression with a depth of ~2.5 km which is comparable to that of our Exp. 3 with a depth of 2.1 km (blue line in Figure 7c).
The fore-arc depression remains above the subduction interface during the steady-state phase in all experiments and is most pronounced for the rollback subduction models, as one would expect from the downward-directed $F_{SU}$ that operates during slab rollback (Figure 8). It thus might appear surprising that the fore-arc depression remains during trench advance and slab rollover. In the trench advance mode, among the two relevant forces, the vertical component of $F_{SU}$ directs upward and thus cannot contribute to the depression. This leaves only $F_{SH\_VER}$ to explain the fore-arc depression. One might expect $F_{SH}$, which depends on $v_s$ and subduction channel thickness at the trench, to be lower in the trench advance mode compared to the trench retreat mode due to a lower $v_s$ during the steady-state phase (Figure 3c). However, $F_{SH}$ is likely comparable or higher because the shear rates in the lubricant at the interface are comparable or greater due to a thinner layer of lubricant in a rollover slab setting (negative $F_{SU}$, i.e., toward the overriding plate, squeezes lubricant out) compared to a rollback slab setting (positive $F_{SU}$, i.e., toward the subducting plate, sucks lubricant in). Furthermore, a larger $F_{SH}$ in the slab rollover setting is consistent with the production of a longer subduction interface. This might also partly explain why $v_s$ is significantly lower in a rollover setting compared to a rollback setting (although part of it is explained by the lower $T_{SP}$ in the rollover models).

### 4.3. Implications for Subduction Zones in Nature

We present a comparison of the overriding plate topography of our rollover experiments with a natural subduction setting that likely presents a rollover slab geometry, namely the Makran subduction zone, indicated by a mantle seismic tomographic model (Amaru, 2007) showing a very steep upper mantle slab and an overall rollover geometry (Figures 9a and 9d). The Makran subduction zone is a relatively narrow subduction zone (~900 km wide), comparable to our experimental subduction zones, that formed at the boundary between the Arabian Plate and the overriding Eurasian plate. Trench migration at the Makran subduction zone has been considered as either close to stationary or advancing with a velocity of ~0–2 cm/yr (Schellart et al., 2008) and the Eurasian overriding plate has experienced overall trench-normal shortening (Burg et al., 2013; Haghipour et al., 2012). The overall slab geometry, trench motion, and overriding plate deformation in the Makran setting are in general agreement with our slab rollover experiments that are dominated by trench advance and overriding plate shortening. In nature, there is also a domain, the coastal Makran region at ~120 and 180 km north of the trench, that is dominated by normal faults and extension. This extensional zone, however, likely does not correspond to the extensional domain found in our rollover experiments, as it is located above the subduction zone plate boundary interface in nature (Normand et al., 2019; Pajang et al., 2021) but above the sublithospheric mantle in our experiments. A possible explanation for the lack of an extensional zone further northward in the Makran domain is the much larger overriding plate (Eurasia), and thus larger resistance to motion, in nature, which would enhance horizontal trench-normal compressive stresses during slab rollover, thereby suppressing the local zone of extension. Regional mantle flow is also one of the factors that can affect subduction style. In the case of the nearby India-Eurasia collision zone, the quasi-toroidal mantle flow that is likely induced at the western syntaxaxis of the collision zone is anticlockwise (see e.g., Figure 1b in Schellart et al. [2019]) which would promote rollback of the Makran slab rather than rollover.

Regarding the large-scale topography, a large part of the Makran accretionary wedge between ~200 and ~400–600 km from the trench (Figure 9b) has a relatively high elevation, and the area further to the north shows a trend with a gradual decrease in elevation. One may consider that this large-scale variation in topography could be caused by the lateral variation in crustal thickness, with higher elevations corresponding to higher crustal thicknesses and vice versa. However, we do not find a good correlation between long-wavelength elevation and crustal thickness (Figures 9b and 9c). Notably, the region of maximum crustal thickness of ~45 km ($x$ = 1,500 km in Figure 9c) is also the region of the lowest overriding plate elevation (Figure 9b). We argue that the discrepancy between observed topography and expected isostatic topographic elevation can be explained by rollover-induced dynamic topography, which has pushed regions closer to the trench with thinner crust to higher elevations. Indeed, the long-wavelength topography of the Makran region, showing a topography that decreases away from the trench (northward) over a distance of ~1,500 km, is consistent with the overriding plate topography of our rollover experiments, showing a comparable topographic decrease over a similar distance, in particular Exp. 1. We thus propose that the long-wavelength topography of the Makran is mostly due to mantle flow and therefore has a dynamic origin.
<figure><img src="image_12.png" /><figcaption>**Figure 9.** (a) Slab geometry as implied by the global P-wave seismic tomography model of Amaru (2007), (b) topographic profile (black line) at the Makran subduction zone (elevation data is from GMTED2010 at USGS website), and the topographic profiles at the end of Exps. 1 (green) and 2 (blue), (c) crustal thickness profile from Crust 1.0 model at http://igppweb.ucsd.edu/~gabi/rem.html (Laske & Masters, 2013). The positions of the tomographic cross-section, topographic profile, and crustal thickness profile at the Makran subduction zone are indicated in (d).</figcaption></figure>

Apart from the long-wavelength topographic signal (~1,000–2,000 km) that we are interested in, the experiments also show a short-wavelength topographic signal (~100 km) that is likely noise, while the Makran topographic profile also shows short-wavelength topography at ~10–100 km wavelength. The latter is due to local crustal deformation processes (e.g., thrust faulting) and surface processes, processes that are not modeled in our analog experiments. From the two experiments showing trench advance and slab rollover, Exp. 1 is likely more applicable to the Makran than Exp. 2, because overriding plate shortening is more extensive in Exp. 1 and thereby is more comparable to the Makran setting, both in distribution and magnitude, and the overriding plate-mantle viscosity ratio in Exp. 1 is more comparable to plate-mantle viscosity estimates in nature (Funiciello et al., 2008; Ribe, 2010; Schellart, 2008a).
## 5. Conclusions

Our buoyancy-driven four-dimensional subduction experiments with an overriding plate produce two subduction styles, namely trench retreat with a backward slab draping (rollback) geometry and trench advance with a forward slab (rollover) geometry. The main findings from the different subduction styles are as follows:

1. Trench advance and slab rollover produce overall shortening in the overriding plate, while trench retreat and slab rollback produce overall extension (Figure 4). Superposed on this overall deformation field, all experiments show fore-arc shortening and intra-arc to near back-arc extension located some 150–300 km from the trench. The main force for OPD is the basal viscous drag force ($F_D$) due to mantle flow, with a trench-directed mantle flow and an increase in trench-directed basal shear rate ($d\dot{\gamma}/dx > 0$) causing overall extension and an OP-directed mantle flow and a decrease in trench-directed basal shear rate ($d\dot{\gamma}/dx < 0$) causing overall shortening (Figure 8). Fore-arc shortening is mostly a result of the subduction interface shear drag force ($F_{SH}$) but is also modulated by the trench suction force ($F_{SU}$) and $F_D$ in the mantle wedge corner

2. The topography of the overriding plate shows large-scale (1,000–2,000 km) overriding plate tilting that depends on the subduction style (Figure 7), with trenchward tilting and subsidence during trench retreat (maximum subsidence of 2–2.5 km), but landward tilting and uplift during trench advance (maximum uplift of 1.5–2 km). These large-scale topographic signatures are a form of dynamic topography as they are caused by the different styles of upper mantle flow during slab rollback and slab rollover (Figures 5 and 6)

3. Our models showing trench advance and slab rollover provide a potential mechanism for the long-wavelength uplift as observed in the overriding plate at the Makran subduction zone (Figure 9)

## Data Availability Statement

The data produced for this research are available in this in-text data citation reference: Xue et al. (2021) under the CC BY 4.0 license.

## References

Adam, J., Urai, J. L., Wieneke, B., Oncken, O., Pfeiffer, K., Kukowski, N., et al. (2005). Shear localisation and strain distribution during tectonic faulting - New insights from granular-flow experiments and high-resolution optical image correlation techniques. *Journal of Structural Geology*, 27(2), 299–301. https://doi.org/10.1016/j.jsg.2004.08.008

Alsaif, M., Garel, F., Gueydan, F., & Davies, D. R. (2020). Upper plate deformation and trench retreat modulated by subduction-driven shallow asthenospheric flows. *Earth and Planetary Science Letters*, 532, 116013. https://doi.org/10.1016/j.epsl.2019.116013

Amaru, M. (2007). Global travel time tomography with 3-D reference models. *Geologica Ultrieaictina*, 274.

Arcay, D., Lallemand, S., & Doin, M. P. (2008). Back-arc strain in subduction zones: Statistical observations versus numerical modeling. *Geochemistry, Geophysics, Geosystems*, 9(5), 1–30. https://doi.org/10.1029/2007GC001875

Bellahsen, N., Faccenna, C., & Funiciello, F. (2005). Dynamics of subduction and plate motion in laboratory experiments: Insights into the “plate tectonics” behavior of the Earth. *Journal of Geophysical Research: Solid Earth*, 110(1), 1–15. https://doi.org/10.1029/2004JB002999

Boutelier, D., & Cruden, A. (2013). Slab rollback rate and trench curvature controlled by arc deformation. *Geology*, 41(8), 911–914. https://doi.org/10.1130/G34338.1

Burg, J.-P., Dolati, A., Bernoulli, D., & Smit, J. (2013). Structural style of the Makran Tertiary accretionary complex in SE-Iran. In K. Al Hosani, F. Roure, R. Ellison, & S. Lokier (Eds.), *Lithosphere dynamics and sedimentary basins: The Arabian Plate and analogues* (pp. 239–259). Springer. https://doi.org/10.1007/978-3-642-30609-9_12

Capitanio, F. A., Morra, G., & Goes, S. (2007). Dynamic models of downgoing plate-buoyancy driven subduction: Subduction motions and energy dissipation. *Earth and Planetary Science Letters*, 262(1–2), 284–297. https://doi.org/10.1016/j.epsl.2007.07.039

Capitanio, F. A., Stegman, D. R., Moresi, L. N., & Sharples, W. (2010). Upper plate controls on deep subduction, trench migrations and deformations at convergent margins. *Tectonophysics*, 483(1–2), 80–92. https://doi.org/10.1016/j.tecto.2009.08.020

Cerpa, N. G., & Arcay, D. (2020). Overriding plate velocity control on surface topography in 2-D models of subduction zones. *Geochemistry, Geophysics, Geosystems*, 21(4), 1–22. https://doi.org/10.1029/2019GC008900

Chen, Z., Schellart, W. P., & Duarte, J. C. (2015). Overriding plate deformation and variability of fore-arc deformation during subduction: Insight from geodynamic models and application to the Calabria subduction zone. *Geochemistry, Geophysics, Geosystems*, 18(1–2), 3697–3715. https://doi.org/10.1002/2015gc005958

Chen, Z., Schellart, W. P., Duarte, J. C., & Strak, V. (2017). Topography of the overriding plate during progressive subduction: A dynamic model to explain fore-arc subsidence. *Geophysical Research Letters*, 44(19), 9632–9643. https://doi.org/10.1002/2017GL074672

Chen, Z., Schellart, W. P., Strak, V., & Duarte, J. C. (2016). Does subduction-induced mantle flow drive backarc extension? *Earth and Planetary Science Letters*, 441, 200–210. https://doi.org/10.1016/j.epsl.2016.02.027

Clark, S. R., Stegman, D., & Müller, R. D. (2008). Episodicity in back-arc tectonic regimes. *Physics of the Earth and Planetary Interiors*, 171(1–4), 265–279. https://doi.org/10.1016/j.pepi.2008.04.012

Cloos, M. (1993). Lithospheric buoyancy and collisional orogenesis - subduction of oceanic plateaus, continental margins, island arcs, spreading ridges, and seamounts. *Geological Society of America Bulletin*, 105, 715–737. https://doi.org/10.1130/0016-7606(1993)105<0715:ibacos>2.3.co;2

Conrad, C. P., & Hager, B. H. (1999). Effects of plate bending and fault strength at subduction zones on plate dynamics. *Journal of Geophysical Research: Solid Earth*, 104(B8), 17551–17571. https://doi.org/10.1029/1999jb900149

## Acknowledgments

We would like to thank reviewer Adam Beall and two anonymous reviewers for their constructive comments that have improved the contents of the paper. This research has been funded by a Vici grant (016. VICI.170.110) from the Dutch National Science Foundation (NWO) awarded to W. P. Schellart. We acknowledge the use of a Lavision stereoscopic PIV equipment and software for visualization.
Crameri, F., & Lithgow-Bertelloni, C. (2018). Abrupt upper-plate tilting during slab-transition-zone collision. *Tectonophysics*, **746**, 199–211.
https://doi.org/10.1016/j.tecto.2017.09.013

Di Giuseppe, E., Van Hunen, J., Funiciello, F., Facenna, C., & Giardini, D. (2008). Slab stiffness control of trench motion: Insights from numerical models. *Geochemistry, Geophysics, Geosystems*, **9**(2), 1. https://doi.org/10.1029/2007GC001776

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

Elsasser, W. M. (1971). Sea-floor spreading as thermal convection. *Journal of Geophysical Research*, **76**(5), 1101–1112. https://doi.org/10.1029/JB076i005p01101

Faccenna, C., Di Giuseppe, E., Funiciello, F., Lallemand, S., & van Hunen, J. (2009). Control of seafloor aging on the migration of the Izu-Bonin-Mariana trench. *Earth and Planetary Science Letters*, **288**(3–4), 386–398. https://doi.org/10.1016/j.epsl.2009.09.042

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

Fuji, N., Jang, H., Nakao, A., Kim, Y. H., Fernández-Blanco, D., Lee, S. M., et al. (2021). A possible roll-over slab geometry under the Caroline plate imaged by Monte Carlo finite-frequency traveltime inversion of teleseismic SS phases. *Frontiers of Earth Science*, **9**, 1–9. https://doi.org/10.3389/feart.2021.593947

Funiciello, F., Facenna, 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

Garel, F., Goes, S., Davies, D. R., Davies, J. H., Kramer, S. C., & Wilson, C. R. (2014). Interaction of subducted slabs with the mantle transition-zone: A regime diagram from 2-D thermo-mechanical models with a mobile trench and an overriding plate. *Geochemistry, Geophysics, Geosystems*, **15**(5), 1739–1765. https://doi.org/10.1002/2014GC005257

Griffiths, R. W., Hackney, R. I., & Van Der Hilst, R. D. (1995). A laboratory investigation of effects of trench migration on the descent of subducted slabs. *Earth and Planetary Science Letters*, **133**, 1–17. https://doi.org/10.1016/0012-821x(95)00027-a

Guillaume, B., Husson, L., Funiciello, F., & Facenna, C. (2013). The dynamics of laterally variable subductions: Laboratory models applied to the Hellenides. *Solid Earth*, **4**(2), 179–200. https://doi.org/10.5194/se-4-179-2013

Haghipour, N., Burg, J. P., Kober, F., Zeilinger, G., Ivy-Ochs, S., Kubik, P. W., & Faridi, M. (2012). Rate of crustal shortening and non-Coulomb behaviour of an active accretionary wedge: The folded fluvial terraces in Makran (SE, Iran). *Earth and Planetary Science Letters*, **355–356**, 187–198. https://doi.org/10.1016/j.epsl.2012.09.001

Hampel, A., Adam, J., & Kukowski, N. (2004). Response of the tectonically erosive south Peruvian fore-arc to subduction of the Nazca Ridge: Analysis of three-dimensional analogue experiments. *Tectonics*, **23**, TC5003. https://doi.org/10.1029/2003TC001585

Harig, C., Zhong, S., & Simons, F. J. (2010). Constraints on upper mantle viscosity from the flow-induced pressure gradient across the Australian continental keel. *Geochemistry, Geophysics, Geosystems*, **11**(6), Q06004. https://doi.org/10.1029/2010GC003038

Hassani, R., Jongmans, D., & Chery, J. (1997). Study of plate deformation and stress in subduction processes using two-dimensional numerical models. *Journal of Geophysical Research*, **102**, 17951–17965. https://doi.org/10.1029/97jb01354

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

Heuret, A., Funiciello, F., Facenna, C., & Lallemand, S. (2007). Plate kinematics, slab shape and back-arc stress: A comparison between laboratory models and current subduction zones. *Earth and Planetary Science Letters*, **256**(3–4), 473–483. https://doi.org/10.1016/j.epsl.2007.02.004

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

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/ggv011

Husson, L. (2006). Dynamic topography above retreating subduction zones. *Geology*, **34**(9), 741–744. https://doi.org/10.1130/G22436.1

Husson, L., Guillaume, B., Funiciello, F., Facenna, C., & Royden, L. H. (2012). Unraveling topography around subduction zones from laboratory models. *Tectonophysics*, **526–529**, 5–15. https://doi.org/10.1016/j.tecto.2011.09.001

Irvine, D. N., & Schellart, W. P. (2012). Effect of plate thickness on bending radius and energy dissipation at the subduction zone hinge. *Journal of Geophysical Research: Solid Earth*, **117**(6), 1. https://doi.org/10.1029/2011JB009113

Jacoby, W. R. (1973). Model experiment of plate movements. *Nature; Physical Science*, **242**, 103–134. https://doi.org/10.1038/physci242130a0

Jarrard, R. D. (1986). Relations among subduction parameters. *Reviews of Geophysics*, **24**(2), 217–284. https://doi.org/10.1029/RG024i002p00217

Lallemand, S., Heuret, A., Facenna, C., & Funiciello, F. (2008). Subduction dynamics as revealed by trench migration. *Tectonics*, **27**(3), 1. https://doi.org/10.1029/2007TC002212

Laske, G., & Masters, G. (2013). Update on CRUST1.0 - a 1-degree global model of Earth's crust. *EGU General Assembly*, **15**, 2658.

Le Pichon, X., Lybéris, N., Angelier, J., & Renard, V. (1982). Strain distribution over the east mediterranean ridge: A synthesis incorporating new Sea-Beam data. *Tectonophysics*, **86**(1–3), 243–274. https://doi.org/10.1016/0040-1951(82)90069-5

Li, Z., & Ribe, N. M. (2012). Dynamics of free subduction from 3-D boundary element modeling. *Journal of Geophysical Research*, **117**, B06408. https://doi.org/10.1029/2012JB009165

Lonergan, L., & White, N. (1997). Origin of the Betic-Rif mountain belt. *Tectonics*, **16**(3), 504–522. https://doi.org/10.1029/96TC03937

Magni, V., Facenna, C., van Hunen, J., & Funiciello, F. (2014). How collision triggers backarc extension: Insight into mediterranean style of extension from 3-d numerical models. *Geology*, **42**(6), 511–514. https://doi.org/10.1130/G35446.1

Malinverno, A., & Ryan, W. B. F. (1986). Extension in the Tyrrenian Sea and shortening in the Apennines as result of arc migration driven by sinking of the lithosphere. *Tectonics*, **5**(2), 227–245. https://doi.org/10.1029/tc005i002p00227

Martinod, J., Guillaume, B., Espurt, N., Facenna, C., Funiciello, F., & Regard, V. (2013). Effect of aseismic ridge subduction on slab geometry and Overriding plate deformation: Insights from analogue modeling. *Tectonophysics*, **588**, 39–55. https://doi.org/10.1016/j.tecto.2012.12.010

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

Normand, R., Simpson, G., & Bahroudi, A. (2019). Extension at the coast of the Makran subduction zone (Iran). *Terra Nova*, **31**(6), 503–510. https://doi.org/10.1111/ter.12419

Pajang, S., Cubas, N., Letouzey, J., Le Pourhiet, L., Seyedali, S., Fournier, M., et al. (2021). Seismic hazard of the western Makran subduction zone: Insight from mechanical modelling and inferred frictional properties. *Earth and Planetary Science Letters*, **562**, 116789. https://doi.org/10.1016/j.epsl.2021.116789

Peltier, W. R. (2004). Global glacial isostasy and the surface of the ice-age earth: The ICE-5G (VM2) model and GRACE. *Annual Review of Earth and Planetary Sciences*, **32**, 111–149. https://doi.org/10.1146/annurev.earth.32.082503.144359
Ribe, N. M. (2010). Bending mechanics and mode selection in free subduction: A thin-sheet analysis. *Geophysical Journal International*, 180(2), 559–576. https://doi.org/10.1111/j.1365-246X.2009.04460.x

Rosenbaum, G., & Lister, G. S. (2004). Neogene and Quaternary rollback evolution of the Tyrrhenian Sea, the Apennines, and the Sicilian Maghrebides. *Tectonics*, 23(1), 1. https://doi.org/10.1029/2003TC001518

Schellart, W. P. (2008a). Kinematics and flow patterns in deep mantle and upper mantle subduction models: Influence of the mantle depth and slab to mantle viscosity ratio. *Geochemistry, Geophysics, Geosystems*, 9, Q03014. https://doi.org/10.1029/2007GC001656

Schellart, W. P. (2008b). overriding plate shortening and extension above subduction zones: A parametric study to explain formation of the Andes Mountains. *Bulletin of the Geological Society of America*, 120(11–12), 1441–1454. https://doi.org/10.1130/B26360.1

Schellart, W. P., Chen, Z., Strak, V., Duarte, J. C., & Rosas, F. M. (2019). Pacific subduction control on Asian continental deformation including Tibetan extension and eastward extrusion tectonics. *Nature Communications*, 10, 1–15. https://doi.org/10.1038/s41467-019-12337-9

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

Schellart, W. P., Stegman, D. R., & Freeman, J. (2008). Global trench migration velocities and slab migration induced upper mantle volume fluxes: Constraints to find an Earth reference frame based on minimizing viscous dissipation. *Earth-Science Reviews*, 88(1–2), 118–144. https://doi.org/10.1016/j.earscirev.2008.01.005

Schellart, W. P., & Strak, V. (2016). A review of analogue modelling of geodynamic processes: Approaches, scaling, materials and quantification, with an application to subduction experiments. *Journal of Geodynamics*, 100, 7–32. https://doi.org/10.1016/j.jog.2016.03.009

Shemenda, A. (1993). Subduction of the lithosphere and back arc dynamics: Insights from physical modeling. *Journal of Geophysical Research*, 98(93), 16167–16185. https://doi.org/10.1029/93JB01094

Strak, V., & Schellart, W. P. (2016). Control of slab width on subduction-induced upper mantle flow and associated upwellings: Insights from analog models. *Journal of Geophysical Research: Solid Earth*, 121(6), 4641–4654. https://doi.org/10.1002/2015JB012545

Van Der Voo, R., Spakman, W., & Bijwaard, H. (1999). Tethyan subducted slabs under India. *Earth and Planetary Science Letters*, 171(1), 7–20. https://doi.org/10.1016/S0012-821X(99)00131-4

Widiyantoro, S., Kennett, B. L. N., & Van Der Hilst, R. D. (1999). Seismic tomography with P and S data reveals lateral variations in the rigidity of deep slabs. *Earth and Planetary Science Letters*, 173(1–2), 91–100. https://doi.org/10.1016/S0012-821X(99)00216-2

Wortel, M. J. R., & Spakman, W. (2000). Subduction and slab detachment in the Mediterranean-Carpathian region. *Science*, 290(5498), 1910–1917. https://doi.org/10.1126/science.290.5498.1910

Xue, K., Schellart, W. P., & Strak, V. (2020). Effect of plate length on subduction kinematics and slab geometry: Insights from buoyancy-driven analog subduction models. *Journal of Geophysical Research: Solid Earth*, 125(11). https://doi.org/10.1029/2020JB020514

Xue, K., Schellart, W. P., & Strak, V. (2021). overriding plate deformation and topography during slab roll-back and slab rollover: Insights from subduction experiments. Retrieved from https://figshare.com/articles/dataset/OverridingPlate_deformation_and_topography_during_slab_rollback_and_slab_rollover_insights_from_subduction_experiments/14226290

Yang, T., Moresi, L., Zhao, D., Sandiford, D., & Whittaker, J. (2018). Cenozoic lithospheric deformation in Northeast Asia and the rapidly-aging Pacific plate. *Earth and Planetary Science Letters*, 492, 1–11. https://doi.org/10.1016/j.epsl.2018.03.057