# The initiation of subduction by crustal extension at a continental margin

F. Lévy, C. Jaupart

► To cite this version:

F. Lévy, C. Jaupart. The initiation of subduction by crustal extension at a continental margin. Geophysical Journal International, 2012, 188, pp.779-797. 10.1111/j.1365-246X.2011.05303.x. insu-03583364

HAL Id: insu-03583364

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

Submitted on 22 Feb 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_1.png" /></figure>
<figure><img src="image_2.png" /></figure>

# The initiation of subduction by crustal extension at a continental margin

F. Lévy\* and C. Jaupart

Institut de Physique du Globe de Paris, Sorbonne Paris Cité, Université Paris Diderot CNRS (UMR 7154), Paris, France.

E-mail: florence.levy@saint-venant-lab.fr

Accepted 2011 November 14. Received 2011 November 14; in original form 2011 May 22

## SUMMARY

We investigate how subduction may be triggered by continental crust extension at a continental margin. The large topography contrast between continental and oceanic domains drives the spreading of continental crust over oceanic basement. Subduction requires the oceanic plate to get submerged in mantle, so that negative buoyancy forces may take over and drive further descent. This is promoted by two mechanisms. Loading by continental crust bends the oceanic plate downwards. Extension in the continental domain induces crustal thinning, which acts to raise mantle above the oceanic plate. In this model, the width of the continental region undergoing extension is an important control parameter. The main physical controls are illustrated by laboratory experiments and simple theory for elastic flexure coupled to viscous crustal spreading. Three governing dimensionless parameters are identified. One involves the poorly constrained oceanic plate buoyancy. We find that the oceanic plate can be thrust to depths larger than 40 km even if it is buoyant, enabling metamorphic reactions and density increase in the oceanic crust. Another parameter is the ratio between the width of the continental extension region and the flexural parameter for the oceanic plate. Initiating subduction is easier if the continent thins over a short lateral distance or if the oceanic plate is strong. The third important parameter is the ratio of oceanic plate thickness to initial continental crust thickness, such that a weak plate and a thick crust do not favour subduction. Thus, the change from a passive to an active margin depends on the local characteristics of the continental crust and is not determined solely by the age and properties of the oceanic lithosphere. It is shown that the spreading of continental crust induces uplift of the margin as the adjacent seafloor subsides. Evidence for the emplacement of continental crust over oceanic basement at passive margins is reviewed.

**Key words:** Subduction zone processes; Continental margins: convergent; Continental tectonics: extensional; Dynamics of lithosphere and mantle; Lithospheric flexure; Mechanics, theory, and modelling.

## 1 INTRODUCTION

Subduction is a key geodynamic process, yet some of its most basic features have not been explained satisfactorily. From the standpoint of thermal convection, the Earth is expected to be in a stagnant lid regime, such that a cold and highly viscous upper boundary layer remains stable and does not participate in convective motions (Ogawa et al. 1991; Davaille & Jaupart 1993; Solomatov & Moresi 2000). From the standpoint of subduction zones as they exist on Earth today, one striking observation is that they are almost entirely absent from the Atlantic Ocean despite the presence of very old seafloor in the North. The most challenging problem about subduction is, perhaps, that oceanic lithosphere may not be denser than the

asthenosphere (Hynes 2005; Afonso et al. 2007). The average density of an oceanic plate depends on the compositions of crust and depleted mantle and on the thermal structure, which evolves as the plate ages and cools down. Using a comprehensive data set, Afonso et al. (2007) have suggested that the oceanic plate remains buoyant with respect to the asthenosphere at all ages, as proposed earlier by Hynes (2005). Density calculations may not be accurate enough yet for definitive statements on plate buoyancy but there can be no doubt that the negative thermal buoyancy is nearly balanced by a positive chemical buoyancy. Thus, subduction may be triggered not by an intrinsic instability process but by some mechanism external to the oceanic plate.

Previous studies have established that self-sustaining subduction driven by negative buoyancy can only occur if the oceanic plate is first pushed into the asthenosphere by some external process (McKenzie 1977; Cloetingh et al. 1989). This is obvious if

<footnote>*Now at: Saint-Venant Laboratory for Hydraulics, Chatou, France.</footnote>
the oceanic plate is buoyant because a phase change is required to augment its density. One candidate is the gabbro-to-eclogite phase change, which occurs at depths of at least 40 km (Peacock 2003). Even if the oceanic plate is denser than the asthenosphere, it has elastic strength, so that it must be bent downwards to some finite depth before negative buoyancy may dominate the force balance, as shown by McKenzie (1977). Ridge push is not large enough for this, implying that other forces must operate. Many different models have been proposed but none has gained widespread acceptance. Loading by sediments can only work for large sediment thicknesses acting on a weak plate with small rigidity, which may be achieved on thin (Cloetingh et al. 1982) or hydrated lithosphere (Regenauer-Lieb et al. 2001). This predicts that subduction gets triggered on seafloor that is either very young or older than 100 Ma. Observations, however, indicate that subduction may begin on oceanic basement of all ages (Jarrard 1986). Other authors have argued that present-day trenches are indirect consequences of older mature subduction zones (Mueller & Phillips 1991; Kemp & Stevenson 1996). Burov & Cloetingh (2010) have recently proposed that a mantle plume impinging on continental lithosphere acts to thrust downwards cold lithospheric mantle slabs. By construction, however, this model requires line plumes to generate downwellings over large horizontal distances and does not address oceanic subduction. The initiation process that has been pursued most extensively involves compression acting on pre-existing zones of weakness in oceanic basement (Toth & Gurnis 1998; Gurnis et al. 2004). In the calculations of Gurnis et al. (2004), the force that allows thrusting along the weak zone is not required to match known plate-driving forces, such as ridge push. This is reasonable for intra-oceanic settings, because one deals with part of a much larger plate boundary, such that the large trench pull exerted by an extensive system of mature subduction zones overwhelms resistance on a small plate segment. This model accounts for many features of young or incipient subduction zones in the Pacific domain (Gurnis et al. 2004). It does not address subduction initiation in primordial ocean basins, however, and leaves us with challenging questions on the fate of recent ocean basins that are devoid of large subduction zones, such as the Atlantic.

As emphasized by McKenzie (1977), one place on Earth where forces are large enough to initiate subduction is at an ocean–continent boundary. Continental crust is thicker and less dense than its oceanic counterpart, and hence, stands at a higher elevation than adjacent seafloor. In such a setting, lateral variations of density and elevation generate horizontal stresses in the 50–100 MPa range, which are larger than the strength in tension of the upper continental crust (Brace & Kohlstedt 1980) and are, therefore, sufficient to induce crustal extension towards the adjacent oceanic domain. These stresses are also large enough for the decoupling of oceanic and continental domains, so that the oceanic plate can be bent downwards into the asthenosphere. Processes at an ocean–continent boundary have been studied using laboratory experiments (Mart et al. 2005; Goren et al. 2008) and numerical calculations (Nikolaeva et al. 2010). In the laboratory studies of Mart et al. (2005) and Goren et al. (2008), self-sustaining subduction was not observed for reasons that are not entirely clear. This may be due to the purely viscous rheologies that were used and we shall see indeed that elastic behaviour is an important ingredient. The numerical model of Nikolaeva et al. (2010) deals with a much more elaborate physical system involving different types of materials in continental and oceanic domains and predicts that self-sustaining subduction can occur in certain conditions. In this model, two types of processes operate: spreading of continental crust onto seafloor and thrusting of continental lithospheric mantle

over oceanic lithosphere. The first is topography-driven, whereas the second one is due to the buoyancy of continental lithospheric mantle. According to Nikolaeva et al. (2010), subduction may not begin in earnest if continental lithosphere has not been thinned to less than about 120 km over a large horizontal distance. Such conditions seem far removed from what is known about continents (Yuan & Romanowicz 2010; Lévy & Jaupart 2011). The conclusions of Nikolaeva et al. (2010) depend on a large number of input parameters, flow laws and initial conditions and it is worth revisiting the mechanism invoked using a different approach. We have been studying this mechanism independently (Lévy 2009; Lévy & Jaupart 2009) and summarize our main results in this paper.

In this study, we aim for a simple diagnostic analysis of necessary conditions for subduction initiation at a passive margin. The interplay between key variables is made explicit through dimensionless numbers. We focus on initial stages where crustal spreading is the dominant process, which allows us to study a reduced physical system involving a minimum number of input parameters. We demonstrate that a key variable is the width of the region of continental extension. We incorporate several factors that Nikolaeva et al. (2010) did not account for, such as oceanic lithosphere that may not be denser than the asthenosphere. In addition, we have made a series of laboratory experiments that illustrate the bending of an elastic plate as material spreads over it and its unbending in later phases. We finally discuss observations that support the spreading of continental crust over oceanic basement, a process which is expected on physical grounds but which has rarely been recognized as such in passive margins. In recent years, a number of detailed geophysical surveys have made it clear that passive margins may not be as passive as they were thought to be and that some of them are affected by extension after the onset of seafloor spreading (Shillington et al. 2006; Peron-Pinvidic et al. 2007; Pallister et al. 2010). For clarity purposes, the discussion of the relevant evidence is deferred until the end of the paper.

The paper is organized as follows. Laboratory experiments are described in the first part and allow some understanding of the key processes at work. In the theoretical part, we derive the governing equations that are required to describe crustal flow and flexure of the oceanic lithosphere. Dimensional analysis leads to a small set of dimensionless numbers that characterize the dynamics and the various regimes in compact form. We then obtain numerical solutions and investigate requirements for the passive-to-active transition at a continental margin. We review evidence for extension at continental margins and emplacement of continental crust over oceanic basement. We discuss the implications of our model for both the initiation of subduction and the dynamics of passive margins.

## 2 LABORATORY EXPERIMENTS

In the laboratory experiments and theory that are described in this paper, continental crust spreads over an oceanic plate, which behaves elastically. This marked contrast in deformation behaviour is due to the large differences of driving force and deformation amplitude that exist: The oceanic lithosphere bends slowly over small vertical distances, whereas large topography gradients drive continental crust extension over long horizontal distances. Such behaviour can be observed in the elaborate numerical calculations of Nikolaeva et al. (2010), which account for the full range of rheological behaviours. We use a viscous rheology to describe the spreading of continental crust but show that it does not limit the validity of the model. We have also sought observations that support the basic
<table><caption>**Table 1.** Properties of materials used in the laboratory experiments.</caption><thead><tr><th colspan="4">PVC elastic sheets</th></tr></thead><tbody><tr><td>v</td><td>Poisson's ratio</td><td>0.4</td><td>—</td></tr><tr><td>E</td><td>Young's modulus</td><td>1.23</td><td>GPa</td></tr><tr><td>ρ</td><td>Density</td><td>1293</td><td>kg m<sup>-3</sup></td></tr><tr><td>h<sub>p</sub></td><td>Thickness</td><td>300 and 690</td><td>μm</td></tr><tr><td>W</td><td>Width</td><td>19 to 19.8</td><td>cm</td></tr><tr><td>L</td><td>Length</td><td>52 to 58</td><td>cm</td></tr><tr><td colspan="4">Natrosol and salt solutions</td></tr><tr><td>ρ<sub>1,2</sub></td><td>Density</td><td>1003–1309</td><td>kg m<sup>-3</sup></td></tr><tr><td>η<sub>1,2</sub></td><td>Viscosity</td><td>0.1–175</td><td>Pas</td></tr></tbody></table>

tenets of the mechanism proposed for subduction initiation. On the one hand, it is known that bending of the oceanic lithosphere at subduction zones can indeed be accounted for by an elastic model, as demonstrated convincingly by Levitt & Sandwell (1995). On the other hand, evidence for the spreading of continental crust over oceanic basement is difficult to find in the literature because it has rarely been sought. We, thus, devote part of the discussion section to several salient features of continental margins.

The purpose of the laboratory experiments is to illustrate a few features of the mechanism proposed for the initiation of subduction. At the start of an experiment, an elastic plate is resting on a viscous fluid, with one edge clamped and the three others free. This represents the oceanic plate that lies next to a continental domain made of a buoyant fluid layer lying on the same lower fluid than in the oceanic domain. The plate gets overridden by the buoyant fluid and gets deflected downwards. The plate may eventually be pushed into the lower fluid, so that it is subjected to a buoyancy force whose sign and magnitude depend on the density contrast. In the experiments, the plate is denser than the lower fluid and, hence, may sink under its own weight in a situation analogous to subduction. We shall show later, however, that a dense plate does not always get submerged in the lower fluid and that a buoyant plate may do so in some cases. Thus, subduction requires specific conditions that

are more involved than restrictions on the sign and magnitude of oceanic plate buoyancy.

## 2.1 Experimental set-up

To model the elastic oceanic lithosphere, we used thin polyvinyl chloride (PVC) sheets allowing deformations that can be observed and measured at the scale of the laboratory. We worked with two thicknesses (300 and 700 μm) and, hence, two rigidities. Continental crust and mantle were represented as viscous fluids. We used aqueous solutions of salt and a cellulose compound called Natrosol (Tait & Jaupart 1989), which set the density and viscosity of the fluids, respectively. We were able to vary independently these two physical properties and investigate a large range of conditions. The physical properties of both the elastic plates and the viscous fluids were measured in our laboratory (see Table 1).

We report on three different sets of experiments using the basic set-up illustrated in Fig. 1. In the first set, a fixed volume of fluid is released at the edge of the plate. Thus, a constant volume of buoyant fluid spreads over the plate so that the same total mass gets distributed over an increasing distance as time increases. In the second set, the buoyant fluid overriding the oceanic plate is fed from the adjacent continental domain and its volume increases with time. The third set of experiments is aimed at the behaviour of a two-layer crust.

## 2.2 Spreading of a constant volume of fluid

At time $t = 0$, a constant volume of buoyant fluid is poured on the edge of the plate from an adjacent reservoir. As shown in Fig. 2, the upper, buoyant, fluid initially spreads over the plate and bends it downwards until the lower fluid begins to submerge the plate. We follow the behaviour of the plate and of the fluid by tracking the depth of plate tip and the leading edges of the upper and lower fluids as they spread horizontally (Fig. 3).

Three phases are observed. The volume of buoyant fluid released on the plate acts as an end load which initially induces a large downward deflection. In the second phase, as the buoyant fluid

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Experimental set-up for laboratory experiments. Two different working fluids with different physical properties and elastic sheets of known properties are used. In one set of experiments (Section 2.2), a slightly different set-up is used. A fixed volume of buoyant fluid is released at one end of the plate and spreads over the plate. In a second set of experiments (Section 2.3), using the set-up shown here, a lock initially separates an oceanic-like domain with an elastic plate resting on dense red fluid and a continental-like domain with buoyant viscous yellow fluid on top of the same red fluid. The lock is lifted at time $t = 0$, allowing the buoyant fluid to undergo extension in the continental domain and spread over the elastic plate. The tip of the plate has a thick front to prevent leakage of small amounts of buoyant yellow fluid below the plate in the first few seconds of an experiment. In nature, the oceanic plate is thick and does not allow such leakage.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Experiment with a fixed volume of buoyant fluid released on the elastic plate. The dashed black line delineates the elastic plate. Arrows indicate the direction of motion of plate tip. These three photographs illustrate three different phases for plate flexure: rapid sinking due to loading (phase A), uplift as the load spreads over a large distance (phase B) and slow sinking due to the negative buoyancy of the plate (phase C). Fluid properties are as follows: $\rho_1 = 1006 \text{ kg m}^{-3}$, $\rho_2 = 1022 \text{ kg m}^{-3}$, $\eta_1 = 17.06 \text{ Pa s}$ and $\eta_2 = 2.45 \text{ Pa s}$. The plate is 300 $\mu\text{m}$ thick.</figcaption></figure>

spreads over an increasing distance, the plate unbends and its tip comes back up. This uplift is an elastic response to the changing distribution of a load of constant volume and we shall see that this is an important phenomenon. Finally, the plate is overthrust by the lower fluid and, once it is pushed into this fluid over some threshold length, sinks under its own weight. In these experiments, self-sustaining subduction only occurs in a late phase, and the initial phases are entirely controlled by the upper fluid, such that the plate responds passively to an evolving load. These experiments also illustrate that the amount of upper fluid that is available is a key parameter.

## 2.3 Spreading of buoyant fluid fed from an adjacent region

In these experiments, the oceanic and continental domains are initially separated by a lock (Fig. 1). On the left, the plate is on top of the lower fluid, as in the previous experiment. On the right, there is a layered system with a layer of buoyant fluid resting on lower fluid. Hydrostatic equilibrium between the two compartments implies that the top of the buoyant fluid lies at a higher elevation than the plate. At time $t = 0$, the lock is opened and the buoyant fluid spreads onto the plate due to gravity. One typical experiment is illustrated in Fig. 4. As extension proceeds in the continental domain, the volume of fluid that spreads over the elastic plate increases, which pushes to plate tip deeper into the lower fluid. In this experiment, lower fluid eventually overthrusts the plate, which then sinks under its own weight. As discussed in the introduction, self-sustaining subduction requires that the oceanic plate gets submerged in mantle material over a finite depth interval. This can be achieved in two different ways, by pushing the oceanic plate into the mantle or by raising the mantle over the plate. Both are driven by crustal extension: the former is due to the spreading of continental crust, whereas the latter is a consequence of crustal thinning.

The experiment shown in Fig. 4 was used to assess the validity of the quantitative physical model that is developed in the following section.

## 2.4 A stratified continental domain with two buoyant layers

The experimental set-up is comparable to the previous one, but the buoyant material is now made of two fluids with the most buoyant one at the top, to mimic a stratified continental crust with different upper and lower crustal layers. The fluid that stands for the lower crust is less viscous and slightly denser than that for the upper crust. Fig. 5 illustrates the time evolution of a typical experiment,

which is again characterized by bending of the elastic plate leading to subduction. One should note that the distal end of the stretched continental crust is made solely of upper crustal material resting on oceanic basement.

## 2.5 Experimental conclusions

These experiments show that elastic deformation of the plate may be responsible for unexpected behaviour, such as an initial phase of subsidence followed by one of uplift as crustal spreading proceeds over an increasing distance. They all lead to the sinking of the elastic plate in the lower fluid under its own weight for the limited range of parameter values that could be investigated in the laboratory. We now investigate the full range of behaviours using theory.

# 3 GOVERNING EQUATIONS

## 3.1 Model set-up

The basic set-up is that of Fig. 6. The oceanic plate is modelled as an elastic plate resting over viscous mantle. It lies next to a continental domain with continental crust of thickness $H'_{Lo}$. In the continent,

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Measurements of flexure and spreading for the fixed-volume experiment of Fig. 2. Depth of the plate tip (bottom panel), and positions of the continental (top panel, dashed line) and mantle (top panel, solid line) fronts as shown as a function of time. Labels A, B and C correspond to the three photographs of Fig. 2.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Experiment where crust-like fluid undergoes extension in the continental domain and spreads over the adjacent elastic plate. The plate is delineated by the dashed line. Only one phase of flexure is observed in this experiment, with the plate sinking continuously. The parameter values of the fluids are $ρ_1 = 1000 \text{ kg m}^{-3}$, $ρ_2 = 1180 \text{ kg m}^{-3}$, $η_1 = 5.23 \text{ Pa s}$ and $η_2 = 3.98 \text{ Pa s}$. The plate is $300 \text{ µm}$ thick.</figcaption></figure>

extension occurs over a finite horizontal distance $X_L$. This distance is dictated by either the width of the continent or the size of the thermally weakened region at the continental edge. Weakening may be due to tectonic and magmatic activity that went on before the onset of seafloor spreading or to a deep thermal perturbation, such as a mantle plume. We shall show that, $X_L$ is an important parameter and, hence, shall vary it within a large range.

We have aimed our physical set-up at a few key features and at the minimum number of parameters. We did not reach for a comprehensive model because it would involve many physical properties and variables that remain poorly constrained. As regards mantle rheology, for instance, not only is the creep mechanism debated but the water content and the various parameters in the creep laws are only known within large uncertainty ranges (Korenaga & Karato 2008). Even if these parameters were known perfectly, one would still need to specify the absolute potential temperature of the Earth's mantle, which is not determined to better than about ±50 °C. This leads to one order of magnitude uncertainty on mantle viscosity. By focusing on the initial stages of crustal spreading, we can dispense with specifying the mantle rheology. As shown by the experiments, initial stages involve crustal thinning in the continent, crustal spreading over oceanic basement and plate flexure due to loading. The asthenosphere beneath the oceanic and continental domains does deform and flow in response to these processes, but it does so over large vertical distances, so that the stresses involved can be neglected in the momentum balance. In essence, the asthenosphere

responds passively to the slow superficial changes that are brought by the spreading crust. In later stages, the oceanic plate may get overridden by mantle, which heralds the beginning of subduction. In this configuration, a thin mantle wedge spreads on top of the oceanic plate and the stresses associated with this mantle flow are no longer negligible. We, therefore, do not deal with this. We define subduction to be such that the oceanic plate gets submerged in mantle material, so that buoyancy effects within the mantle enter the force balance on the oceanic plate. This definition allows an unambiguous separation between subduction and no subduction. That the oceanic plate begins to get submerged by mantle material does not guarantee self-sustaining subduction, as discussed above, but it is a pre-requisite. We return to this question later.

Crustal rheologies and temperatures are other grey areas of geological knowledge, with additional difficulties due to the large lateral variations of thermal structure that are associated with changes of radiogenic heat production (Lévy & Jaupart 2011). We discuss the impact of crustal rheology on the model results in the discussion section. We show that the crustal rheology only sets the timescale for the spreading and, hence, can be factored in easily.

As will be shown below, even when it is reduced to a few essentials, our model involves a large number of variables and physical properties. To illustrate the collective impact of the most important parameters and to evaluate the combinations that are most favourable to subduction, we have extracted the relevant dimensionless numbers out of the governing equations.
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Experiment with two buoyant fluid layers mimicking upper and lower crust. Note that the lower liquid (analogous to the lower crust) does not spread over a large distance and is missing from the distal region. The fluid properties are $\rho_{1u} = 1010 \text{ kg m}^{-3}$, $\rho_{1l} = 1090 \text{ kg m}^{-3}$, $\rho_2 = 1200 \text{ kg m}^{-3}$, $\eta_{1u} = 11.16 \text{ Pa s}$, $\eta_{1l} = 1.67 \text{ Pa s}$ and $\eta_2 = 3.66 \text{ Pa s}$. The plate is 300 $\mu\text{m}$ thick.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Diagram showing the variables that are used for a model of crustal spreading on oceanic basement at a continental margin.</figcaption></figure>

## 3.2 Flexure of an elastic plate

Most geological models of elastic plate bending rely on the assumption of small displacements (e.g. Turcotte & Schubert 1982). The vertical deflection of the plate, $w$, is the solution of a fourth-order

differential equation,

$$ \frac{d^2}{dx^2} \left( D \frac{d^2 w}{dx^2} \right) + F \frac{d^2 w}{dx^2} - q(x) = 0, \quad (1) $$
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Flexure of an elastic plate bent under its own weight. The plate properties are $L = 20$ cm, $h_p = 300$ µm, $\rho = 1296$ kg m$^{-3}$, $E = 1.2$ GPa and $\nu = 0.4$. Stars indicate measurements made with a laser positioning device. Theoretical predictions made with the small-displacement approximation (dashed line) and with beam equations for large deflections (continuous line) are shown. The small-displacement solution overestimates the vertical deflection and the plate length. Furthermore, it does not allow for the opening of a gap at the tip of the downgoing plate.</figcaption></figure>

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

where $x$ is the horizontal coordinate, $D$ is the flexural rigidity of the plate, $F$ is a horizontal force per unit thickness applied at the edge of the plate and $q(x)$ is a vertical downward normal stress. For our present purposes, this model is not appropriate. To illustrate why, we compared the solution for a plate bending under its own weight in air with experimental measurements. In this case, no horizontal force is applied to the plate ($F = 0$), and the vertical load is the weight of the plate,

$$q(x) = -\rho g h_p, \qquad (2)$$

where $\rho$ is the plate density, $h_p$ is the plate thickness and $g$ is the acceleration of gravity. The boundary conditions are $w = 0$ and $dw/dx = 0$ at the clamped edge of the plate at $x = 0$. No bending moment and no stress are applied at the free tip of the plate located at $x = L$, so that we specify that $d^2w/dx^2 = 0$ and $d^3w/dx^3 = 0$ there. We obtained solutions for the PVC sheets of our experiments. Fig. 7 shows the calculated and measured deflections for a 20-cm-long plate. The solution clearly overestimates the deflection. Moreover, as the free tip of the plate is supposed to remain at a fixed horizontal distance equal to $L$ from the clamped edge, this theory implies an unrealistic increase of plate length and does not allow for the opening of a large gap at that end. Both are important factors in our problem and must be accounted for (Fig. 4).

We use beam theory for the finite amplitude flexure of an elastic plate and use local coordinates along the deformed plate: curvilinear abscissa $s$ along the plate and deflection angle $\alpha$, defined as the angle between the horizontal and the tangent vector to the plate $\vec{t}$ (see Fig. 8). The force and momentum balances lead to a system of four differential equations [see Timoshenko & Young (1945), pp. 424–427],

<figure><img src="image_11.png" /><figcaption>**Figure 8.** Definition of the variables involved in the beam equations for flexure with large displacements.</figcaption></figure>

$$\left\{ \begin{aligned} \frac{dp}{ds} &= \tau \frac{d\alpha}{ds} & -\rho g \sin \alpha \\ \frac{d\tau}{ds} &= -p \frac{d\alpha}{ds} & -\rho g \cos \alpha \\ \frac{d^2\alpha}{ds^2} &= \frac{12 \tau}{h_p^2 A} \\ \frac{dw}{ds} &= \sin \alpha \end{aligned} \right. , \qquad (3)$$

where $p$ is the mean axial stress (i.e. stress oriented along vector $\vec{t}$, averaged over the plate thickness) and $\tau$ is the mean shearing stress (i.e. oriented along the normal vector $\vec{n}$, averaged over the plate thickness; see Fig. 8). $A$ is an elastic coefficient defined as,

$$A = \frac{(1 - \nu)}{(1 - 2\nu)(1 + \nu)} E,$$

where $\nu$ is Poisson's ratio and $E$ is Young's modulus.

The boundary conditions are $w = 0$ and $\alpha = 0$ at $s = 0$, and
$d\alpha/ds = 0$ (no moment applied), $\tau = 0$ and $p = 0$ at $s = L$. The solution obtained for the same sheet as before is shown on Fig. 7. Unlike the small-displacement solution, the numerical results fit the experimental data very well. Perfect agreement could be achieved with a slight modification of either the density or the elastic coefficient for PVC within bounds of uncertainty. Note the large gap that opens up at the end of the plate, which plays an important role in the process of crustal spreading and thinning.

### 3.3 Spreading of a viscous fluid

Continental crust is extremely heterogeneous so that it is difficult to derive its rheological behaviour from laboratory data on individual rock types. Large-scale deformation proceeds through two different mechanisms. In the upper crust, deformation is distributed over fault arrays and involves block structures (King et al. 1994; Thatcher 1995). This can be accounted for in active regions, where the main faults have been mapped and kinematics are well understood. Deeper crustal levels deform continuously, which allows implementation of dynamic models where the deformation rate is solved for as a function of the driving forces involved. Using the approximation of a deforming continuum, which is either viscous or plastic for the whole crust has met with considerable success and will be employed here (Bird & Piper 1980; England & McKenzie 1982; Sonder & England 1986). For our present purposes, the exact rheology is not critical because we focus on the elastic behaviour of the oceanic plate. The rheology determines the timescale and geometrical characteristics of the spreading. As regards time, we shall show that the effective viscosity of the crust only serves to set the timescale and shall review the relevant constraints in the discussion section. Geometrical aspects, such as lateral thickness variations of the crust that spreads over the elastic plate, only play a minor role. For instance, modelling the crust as a granular medium would lead to a triangular-shaped leading edge and we have verified that a triangular load has the same influence on plate behaviour as a viscous gravity current with a rounded tip. We, therefore, adopt a bulk viscous rheology for the crust and describe extension and spreading in two dimensions using Navier-Stokes equations. For crust extending over a laterally extensive region, it is appropriate to use the lubrication approximation for crustal flow (McKenzie et al. 2000),

$$ \left\{ \begin{array}{l} \frac{\partial P}{\partial x} = \eta_1 \frac{\partial^2 u}{\partial z^2}, \\ \\ \frac{\partial P}{\partial z} = -\rho_1 g \end{array} \right. , \qquad (4) $$

where $z$ is vertical coordinate, $P$ is pressure, $u$ is horizontal velocity, $\eta_1$ is the fluid viscosity and $\rho_1$ is density. We also use the continuity equation to derive an equation for crustal thickness $H$. The equations were solved numerically using an implicit method. The validity of the numerical code was assessed through a test comparison with the spreading of a 2-D gravity current over a horizontal plane with a constant volume. In this case, the time evolution of the fluid thickness $H$ is described by

$$ \frac{\partial H}{\partial t} = \frac{\rho_1 g}{3 \eta_1} \frac{\partial}{\partial x} \left( H^3 \frac{\partial H}{\partial x} \right). \qquad (5) $$

Scaling analysis indicates that thickness $H$ and leading edge position vary as $t^{-1/5}$ and $t^{1/5}$, respectively. We verified that the code did lead to such scalings.

For spreading over a bent substratum with deflection $w(x)$, eq.

(5) becomes

$$ \frac{\partial h}{\partial t} = \frac{\rho_1 g}{3 \eta_1} \frac{\partial}{\partial x} \left( H^3 \frac{\partial h}{\partial x} \right), \qquad (6) $$

with $H(x, t)$ is the total thickness of fluid and $h(x, t) = H(x, t) + w(x)$, where $w(x) \le 0$ stands for elevation above a reference level.

### 3.4 Coupling of elastic flexure and viscous spreading

We consider an elastic plate lying over a fluid of density $\rho_2$ which gets deflected by the spreading of a fluid of density $\rho_1$ and viscosity $\eta_1$. Eq. (3), describing the flexure of a plate under its own weight, must be modified to account for the pressures $P_l$ and $P_u$ that are exerted at the base and top of the plate, respectively. In a thin spreading fluid, the pressure field is locally identical to the hydrostatic field. For the lower fluid, which extends over a very large thickness, viscous stresses are negligible and the pressure field can also be written in the hydrostatic approximation. The elastic beam eq. (3) need only be modified in the expression for shearing stress $\tau$,

$$ \frac{d\tau}{ds} = -p \frac{d\alpha}{ds} - \rho g \cos \alpha + \frac{1}{h_p} (P_l - P_u). \qquad (7) $$

The pressure distributions at the top and bottom of the plate change as the upper fluid spreads. We use a shooting method to solve the plate equations. The plate flexure profile at time $t$ is used to calculate the spreading of the fluid between $t$ and $t + dt$, using eq. (6). Plate flexure is then determined for the new pressure distribution (at $t + dt$) using eq. (7). The convergence of the solutions was established by changing the time step in a systematic manner.

For geological cases (Fig. 6), fluids 1 and 2 correspond to continental crust and mantle, respectively. The elastic plate, which corresponds to the oceanic lithosphere, initially lies under a layer of water of thickness $h_o$ and density $\rho_o$. The pressure acting on the lower surface of the plate is, therefore,

$$ P_l = \rho_o g h_o - \rho_2 g w + \rho g h_p, \qquad (8) $$

and the pressure acting on the upper surface of the plate is,

$$ \begin{cases} P_u = \rho_1 g H & \text{for } h > h_o \\ P_u = \rho_1 g H + \rho_o g (h_o - h) & \text{for } h < h_o, \end{cases} \qquad (9) $$

where $H$ is the thickness of continental crust overthrusting the plate and $h$ is elevation relative to the undeflected seafloor. $h_p$ is the elastic thickness of the oceanic plate; $h_p$ and $h_o$ are assumed to remain constant.

Away from the oceanic plate, the continental crust and mantle are in hydrostatic equilibrium. There, thickness $H'_l$ of the continental crust (Fig. 6) is specified by the following equation:

$$ H'_l = \frac{\rho_2}{\rho_2 - \rho_1} (h_L + h_p) - \frac{\rho_o}{\rho_2 - \rho_1} h_o - \frac{\rho}{\rho_2 - \rho_1} h_p, \qquad (10) $$

where $h_L$ is the elevation of the crust located beyond the plate end. Away from the oceanic plate, crustal thinning affects a large thickness so that viscous stresses are small compared to those in the thinner crustal regions. This implies that lateral thickness variations are small and allows the useful simplification that $H'_l$ and $h'_L$ are both uniform laterally.

The total volume of continental crust is conserved, which implies the following condition on the flux $\Phi$ of continental material at the
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Calculations of plate flexure for the parameter values of the laboratory experiment in Fig. 4. Top: initial condition for the calculation. Bottom: plate flexure due to the spreading of the viscous fluid 9 s after the start of the experiment. Calculated displacements and fluid thicknesses are very close to the experimental observations in Fig. 4.</figcaption></figure>

continent/ocean boundary,

$$ \left\{ \begin{aligned} \Phi &= -\frac{\rho_1 g}{3 \eta_1} H^3 \frac{\partial h}{\partial x} \\ \Phi &= -\frac{d}{dt} \{ H'_{L} [(L - x_k) + X_L] \}, \end{aligned} \right. \qquad (11) $$

where $L - x_k$ denotes the width of the opening that has been made available by plate bending. $H'_L$ decreases as the continental crust spreads on the oceanic plate. At each time step, $H'_L$ is determined by the constraints of hydrostatic equilibrium (eq. 10) and volume conservation (eq. 11).

We assessed the validity of the model equations through a comparison with data from the experiment of Fig. 4. Fig. 9 shows calculations of plate deflection and upper fluid thickness, which can be compared to the first experimental snapshot at the top of Fig. 4. Over the short time span of the experiments, the initial time is ill-defined because it is not possible to lift the lock instantaneously. This did not allow a precise comparison between experimental data and theory, but we did verify that the calculations reproduce the observed characteristics and amplitude of the deformation.

## 3.5 Two different regimes

Figs 10 and 11 illustrate the two different behaviours that have been obtained. To describe the solutions, we focus on three particular variables, the location of the tip of the spreading crust, $x_o$, the depth of the plate tip, $w_L$, and the depth difference between this tip and the continental Moho, $\Delta h$. For subduction to occur, the oceanic plate must be submerged in mantle material, which implies that $\Delta h$ must be negative with the sign convention adopted here. For the two cases that are described in this section, all parameters are identical save for the plate buoyancy, which is either positive or negative (Table 2). Times that are indicated in Figs 10 and 11 have been obtained for a bulk crustal viscosity of $5 \times 10^{20}$ Pa s. We show later (Section 4.2) that time scales linearly with the crustal viscosity, so that one may easily evaluate the impact of other viscosity values.

In both calculations, we find that plate gets pushed downwards as crustal extension proceeds. With time, spreading of the increas-

ingly thin crust decelerates and the downward bending of the plate slows down. In many cases, the plate eventually comes back up, a behaviour which was illustrated by the first set of laboratory experiments described above. Simultaneously, thinning of the adjacent continental crust continues. $\Delta h$ is the difference between the depth of the plate tip and the base of the adjacent continental crust, such that it is negative when the oceanic plate is submerged in the asthenosphere. The sign of $\Delta h$ provides a straightforward criterion to discriminate between two regimes. The first regime is such that $\Delta h$ tends to a constant positive limit noted, $\Delta h_f$, and is illustrated in Fig. 10 for a buoyant oceanic plate. In this case, $\Delta h$ remains positive at all times and the plate tip never gets submerged in the asthenosphere, implying that subduction is prevented. In the other regime, $\Delta h < 0$ and the oceanic plate can proceed to self-sustaining subduction. Fig. 11 shows one example of this behaviour for a negatively buoyant plate. We observe that $\Delta h$ becomes negative as the plate dips below the continental crust. In this case, mantle material overrides the oceanic plate and the negative plate buoyancy drives further subduction, as observed in the laboratory. As stated earlier, we do not consider time-dependent calculations beyond that point because we focus on conditions for the initiation of subduction. The two regimes can be conveniently separated by the simple threshold of $\Delta h_f = 0$. The two different regimes have been illustrated for plates that were positively and negatively buoyant, but we shall see that subduction can, in fact, be initiated for buoyant oceanic lithosphere. We now investigate a large range of parameters and properties and extract the important dimensionless numbers that characterize the solutions.

# 4 CONDITIONS FOR THE INITIATION OF SUBDUCTION ON EARTH

## 4.1 Geological parameter values

The large number of parameters and physical properties that are involved in this simple physical model are recapitulated in Table 3. Some of them are not well constrained and dimensional analysis is useful in providing general results that can be applied to a large range
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Calculations of margin evolution with the parameter values listed in Table 2. Top left: location of the leading edge of spreading crust as a function of time. Bottom left: depth of the plate tip ($w_L$) and depth difference between the plate tip and the base of the continental crust ($\Delta h$) as a function of time. Red dots mark times corresponding to profiles shown on the right-hand side figures. The plate tip never sinks below the continental crust ($\Delta h$ remains positive): the plate never gets submerged in mantle material and no subduction can occur. Right: plate flexure as three different times showing that subduction does not occur.</figcaption></figure>

of parameter values. For example, as mentioned earlier, the density contrast between oceanic lithosphere and asthenosphere cannot be measured directly and is calculated using thermodynamic models for melt depletion and physical properties for a large number of mineral phases. Results rely on the assumption of thermodynamic equilibrium as well as on an assumed water content in the upper mantle and, hence, cannot be taken as perfectly exact. We, therefore, allow for both positively and negatively buoyant oceanic lithospheres. The intrinsic density contrast due to melt depletion in mantle rocks cannot exceed $50 \text{ kg m}^{-3}$, a value which is only achieved in Archean continental lithosphere (Poudjom Djomani et al. 2001; Schutt & Lesher 2006). The oceanic plate density increases as the plate cools down and thickens, but we consider plate buoyancy and flexural rigidity as two independent parameters for clarity purposes.

Viscosities in the $10^{18}-10^{22} \text{ Pa s}$ range are used to describe the bulk continental crust (e.g. Block & Royden 1990; Kruse et al. 1991; Kaufman & Royden 1994; McKenzie et al. 2000; Peltier & Drummond 2008). We will show that this viscosity only sets the timescale of the process and does not affect the loading and bending of the plate, that is, whether subduction can be initiated or not. It is, therefore, not crucial to specify this viscosity precisely and a straightforward scaling procedure allows one to evaluate the impact of other viscosity values.

Values for the initial thickness of continental crust, $H'_{Lo}$ and for the size of the extension zone, $X_L$, are taken within rather large ranges. At passive margins, crustal thickness varies laterally with distance to the ocean-continent boundary. These variations may be attributed, in part, to rifting processes before ocean opening and, in part, to later phases of crustal spreading driven by topography, as studied here. This may be ascertained by considering an ocean basin in the making, such as the Red Sea, for example (Al-Damegh et al. 2005). There, one observes that the magnitude of extension in the young margins is about half as large as those of the older passive margins of the Atlantic Ocean. Thus, the initial crustal thickness $H'_{Lo}$ may be close to the value estimated in stable regions, that is, between 35 and 45 km (Christensen & Mooney 1995; Mooney et al. 1998). The width of the extension zone, $X_L$, depends on thermal events that were active before seafloor spreading and may vary amongst geological provinces. We use seismic studies of existing passive margins to set the range of values for this important parameter. For instance, the thickness of the continental crust decreases from 35 km to essentially zero over a distance of about 500 km in the Norwegian margin (Scheck-Wenderoth et al. 2007). In other regions, continental crust is stretched over smaller distances, for instance over about 180 km in the Nova Scotia margin (Funck et al. 2004).
<figure><img src="image_14.png" /><figcaption>**Figure 11.** Same as Fig. 10, but for a different plate density (Table 2). The elastic plate eventually sinks below the continental crust and subduction occurs. Curves are shown in thin lines once the plate tip reaches into the mantle, because the overthrusting of the plate by mantle material is not modelled accurately.</figcaption></figure>

<table><caption>**Table 2.** Parameter values and $\Delta h_f$ for the two examples shown in Figs 10 and 11.</caption><thead><tr><th></th><th>ν<br>(-)</th><th>ρ<br>(kg m<sup>-3</sup>)</th><th>L<br>(km)</th><th>E<br>(GPa)</th><th>h<sub>p</sub><br>(km)</th><th>δ<sub>p</sub><br>(km)</th><th>η<br>(Pa s)</th><th>X<sub>L</sub><br>(km)</th><th>H'<sub>Lo</sub><br>(km)</th><th>ρ<sub>1</sub><br>(kg m<sup>-3</sup>)</th><th>ρ<sub>2</sub><br>(kg m<sup>-3</sup>)</th><th>Δh<sub>f</sub><br>(km)</th></tr></thead><tbody><tr><td>Fig. 10</td><td>0.25</td><td>3000</td><td>800</td><td>30</td><td>30</td><td>81</td><td>5.10<sup>20</sup></td><td>75</td><td>42</td><td>2700</td><td>3300</td><td>13</td></tr><tr><td>Fig. 11</td><td>-</td><td>3350</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-</td><td>-4</td></tr></tbody></table>

## 4.2 Dimensional analysis

To identify the combinations of parameters that control the margin and plate behaviours, we proceed with a dimensional analysis and write all the equations in dimensionless form:

(i) equations for elastic flexure (eq. 7),

(ii) the momentum equation for crustal spreading (eq. 6),

(iii) the conservation of the total volume of continental crust (eq. 10),

(iv) hydrostatic equilibrium (eq. 11).

Plate bending is achieved over horizontal distances that depend on the flexural parameter, $\delta_p$, which has the dimensions of length and provides the relevant length-scale,

$$ \delta_p = \left( \frac{Ah_p^3}{\Delta\rho g} \right)^{\frac{1}{4}}, \qquad (12) $$

where $\Delta\rho$ is the density difference between plate and water ($\Delta\rho =$

<table><caption>**Table 3.** Parameters involved in the calculations and ranges of values for geological cases.</caption><tbody><tr><td colspan="3" style="text-align: center;">Mantle</td></tr><tr><td>ρ<sub>2</sub></td><td>Density</td><td>~3300 kg m<sup>-3</sup></td></tr><tr><td colspan="3" style="text-align: center;">Oceanic lithosphere</td></tr><tr><td>ρ</td><td>Density</td><td>3250–3350 kg m<sup>-3</sup>?</td></tr><tr><td>L</td><td>Length</td><td>0–1000 km</td></tr><tr><td>ν</td><td>Poisson's ratio</td><td>0.23–0.33</td></tr><tr><td>E</td><td>Young's modulus</td><td>50–100 GPa</td></tr><tr><td>h<sub>p</sub></td><td>Elastic thickness</td><td>5–70 km</td></tr><tr><td colspan="3" style="text-align: center;">Continental crust</td></tr><tr><td>ρ<sub>1</sub></td><td>Density</td><td>2700–2900 kg m<sup>-3</sup></td></tr><tr><td>η<sub>1</sub></td><td>Viscosity</td><td>10<sup>20</sup>–10<sup>22</sup> Pa s</td></tr><tr><td>H'<sub>Lo</sub></td><td>Thickness</td><td>35–40 km</td></tr><tr><td>X<sub>L</sub></td><td>Thinning width</td><td>100–? km</td></tr></tbody></table>

$\rho_2 - \rho_o$). We deduce from the flexure equations a stress scale $[\tau]$,

$$ [\tau] = A \frac{h_p^2}{\delta_p^2}. \qquad (13) $$
The dimensionless equations for plate flexure are thus,

$$
\left\{
\begin{aligned}
\frac{\mathrm{d}p^*}{\mathrm{d}s^*} &= \tau^* \frac{\mathrm{d}\alpha}{\mathrm{d}s^*} - \frac{\rho}{\Delta\rho} \frac{h_p}{\delta_p} \sin\alpha \\
\frac{\mathrm{d}\tau^*}{\mathrm{d}s^*} &= -p^* \frac{\mathrm{d}\alpha}{\mathrm{d}s^*} - \frac{\rho}{\Delta\rho} \frac{h_p}{\delta_p} \cos\alpha + \frac{\delta_p}{h_p} (P_l^* - P_u^*) \\
\frac{\mathrm{d}^2\alpha}{\mathrm{d}s^{*2}} &= 12\tau^* \\
\frac{\mathrm{d}w^*}{\mathrm{d}s^*} &= \sin\alpha
\end{aligned}
\right.
, \quad (14)
$$

where asterisks indicate dimensionless variables.

To describe the loading of the plate in dimensionless form (eqs 8 and 9), we scale the plate deflection, w, and crustal thickness, H, with the initial crustal thickness, $H'_{Lo}$. We obtain the following equations for pressure:

$$P_l^* - P_u^* = \frac{h_p}{\delta_p} \left[ \frac{\rho}{\Delta\rho} \left( \frac{h_p}{\delta_p} \right) - \frac{H'_{Lo}}{\delta_p} w^* - \frac{\rho_1 - \rho_o}{\Delta\rho} \frac{H'_{Lo}}{\delta_p} H^* \right], \quad (15)$$

in the oceanic domain ($h_o > h_L$) and

$$P_l^* - P_u^* = \frac{h_p}{\delta_p} \left[ \frac{\rho_o}{\Delta\rho} \left( \frac{h_o}{\delta_p} \right) + \frac{\rho}{\Delta\rho} \left( \frac{h_p}{\delta_p} \right) - \frac{\rho_2}{\Delta\rho} \frac{H'_{Lo}}{\delta_p} w^* - \frac{\rho_1}{\Delta\rho} \frac{H'_{Lo}}{\delta_p} H^* \right], \quad (16)$$

in the continental domain ($h_o < h_L$).

Timescale [t] is deduced from the crustal spreading equation:

$$[t] = \frac{3\eta'_1\delta_p^2}{\rho_1 g H_{Lo}^{3/2}} \quad (17)$$

and the differential equation for crustal thickness H becomes

$$\frac{\partial H^*}{\partial t^*} = \frac{\partial}{\partial x^*} \left( H^{*3} \frac{\partial h^*}{\partial x^*} \right). \quad (18)$$

The other equations are

$$H^{*3} \frac{\partial h^*}{\partial x^*} = \frac{X_L}{\delta_p} \frac{\mathrm{d}H_L^{*3}}{\mathrm{d}t^*}, \quad (19)$$

and

$$H_L^{*3} = \frac{\rho_2}{\rho_2 - \rho_1} \left( h_L^* + \frac{\rho_2 - \rho}{\rho_2} \frac{h_p}{H_{Lo}'} - \frac{\rho_o}{\rho_2} \frac{h_o}{H_{Lo}'} \right). \quad (20)$$

We note that the initial continental elevation does not appear explicitly because it is set by the crustal thickness.

$$h_L = \frac{\rho_2 - \rho_1}{\rho_2} H_L' + \frac{\rho - \rho_2}{\rho_2} h_p + \frac{\rho_o}{\rho_2} h_o. \quad (21)$$

We also note that the bulk continental viscosity only appears through timescale [t] and does not enter the dimensionless equations. We have verified that, with all other parameters fixed, $\Delta h_L^*$ stabilizes at the same value regardless of the chosen crustal viscosity.

All in all, 10 dimensionless numbers appear in this large set of eqs (14)–(20), including a host of density ratios and their differences. Fortunately, most of them are not important because they do not affect the quantitative results significantly when they are varied within their respective ranges. As regards density, for instance, only differences are significant and likely to change by significant amounts. The densities of continental crust, water or sublithospheric mantle cannot deviate much from the chosen values. In the end, we find that only three dimensionless numbers are truly important: the dimensionless plate buoyancy, $(\rho_2 - \rho)/\rho_2$, a thickness ratio, $h_p/H'_{Lo}$ and a width ratio, $X_L/\delta_p$.

We have seen that the two regimes for the oceanic plate can be characterized by the value of $\Delta h_f$, which is the depth difference between the plate tip and the base of the continental crust. In dimensionless form, this depth difference can be written as a function of the three key numbers,

$$\frac{\Delta h_f}{H'_{Lo}} = f \left( \frac{\rho_2 - \rho}{\rho_2}, \frac{h_p}{H'_{Lo}}, \frac{X_L}{\delta_p} \right). \quad (22)$$

For the parameter values and physical properties that have been discussed earlier, appropriate ranges for the three important dimensionless numbers are as follows:

$$
\left\{
\begin{array}{l}
-0.015 \le \frac{\rho_2 - \rho}{\rho_2} \le 0.015 \\
0.125 \le \frac{h_p}{H'_{Lo}} \le 1.75. \\
0.1 \le \frac{X_L}{\delta_p} \le 10
\end{array}
\right.
\quad (23)
$$

We discuss the influence of each one of them in the following section.

## 4.3 Three dimensionless numbers governing oceanic plate behaviour

### 4.3.1 Dimensionless plate buoyancy, $(\rho_2 - \rho)/\rho_2$

The influence of this dimensionless number on plate behaviour is quite obvious: the denser the plate is, the easier subduction is. We find that final position of the plate tip, which is measured by the value of $\Delta h_f^*$, increases linearly with plate buoyancy. Two sets of results with the same value of the thickness ratio, $h_p/H'_{Lo}$, and two different values of the distance ratio, $X_L/\delta_p$, are illustrated in Fig. 12. For a large distance ratio of 4.32, subduction initiation is predicted only for density ratios that smaller than −0.024, that is,

<figure><img src="image_15.png" /><figcaption>**Figure 12.** Limit value (dimensionless) reached by $\Delta h$, the depth difference between the plate tip and the base of the continental crust, as a function of dimensionless plate buoyancy. For the sign convention adopted in the calculations, $\Delta h$ becomes negative when the plate tip reaches into the mantle. Thus, subduction does not occur if $\Delta h$ tends to a positive limit value. These results show that negative plate buoyancy enhances the likelihood of subduction, as expected. Note, however, that subduction of negatively buoyant plates is prevented in some cases. Furthermore, buoyant plates can be thrust below continental crust in some other cases.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>**Figure 13.** Limit value reached by dimensionless depth difference $\Delta h$ as a function of thickness ratio $h_p/H'_{Lo}$. Subduction requires that $\Delta h$ becomes negative (such that the plate tip sinks below the continental crust). A thick plate or a thin continental crust are favourable to subduction.</figcaption></figure>

for oceanic lithosphere that is at least 80 kg m⁻³ denser than the underlying mantle. For a smaller distance ratio of 1.67, subduction can occur even for a buoyant plate [for $(\rho_2 - \rho)/\rho_2 \le +0.014$. In this case, however, self-sustaining subduction will proceed only if the plate is forced to a depth that allows metamorphic reactions within the oceanic crust. These calculations show that the density of oceanic lithosphere may, in fact, not be an important control on subduction.

### 4.3.2 Thickness ratio, $h_p/H'_{Lo}$

This number affects the depth that can be reached by the plate tip. For a small plate thickness, loading does not induce flexure over a large lateral distance and the force balance in distal regions reduces to a hydrostatic equilibrium with buoyant crustal material over oceanic basement. In this case, the plate is flat away from the spreading crust and its tip lies at a depth which is dictated by the thickness of crustal material that has spread over it. In contrast, a thick elastic plate gets bent over a large lateral distance, which acts to push the tip to large depths. The thickness ratio also involves the initial crustal thickness. What matters here is not the magnitude of the crustal load on the oceanic plate, but the depth that the plate tip must reach to be submerged in mantle material. With an initially thin crust, and hence, a large thickness ratio, the plate does not need to be bent strongly to reach into the mantle.

We, thus, expect that the likelihood of subduction gets enhanced by decreasing the thickness ratio. We find indeed that the position of the plate tip relative to the base of the crust, $\Delta h'_t$, decreases linearly as $h_p/H'_{Lo}$ increases (Fig. 13). In the Earth, elastic thickness increases with age and we conclude that subduction is easier for old lithosphere.

### 4.3.3 Width ratio, $X_L/\delta_p$

The size of the continental extension region comes into play because it affects the amount of crustal thinning that occurs due to spreading on the oceanic plate. The width ratio scales the size of the extension region to the flexural parameter. As for the previous dimensionless

<figure><img src="image_17.png" /><figcaption>**Figure 14.** Limit value reached by dimensionless depth difference $\Delta h$ as a function of distance ratio $X_L/\delta_p$. Subduction requires that $\Delta h$ becomes negative (such that the plate tip sinks below the continental crust). A narrow extension zone (small values of $X_L$) is favourable to subduction.</figcaption></figure>

number, the ratio involves two competing effects on the depth of the plate tip. The influence of the flexural parameter is similar to that of the elastic thickness because the two parameters are linked to one another (eq. 12). A thick plate and a large flexural parameter, $\delta_p$, allow large downward deflections. Even a thin plate can get submerged in mantle material, however, if the continental extension zone is narrow. In this case, large amounts of crustal thinning in the continental domain act to raise continental mantle over the plate even at small degrees of flexure.

We indeed find that, $\Delta h'_t$, decreases with decreasing width ratio (Fig. 14). There is a rather sharp transition between behaviours at small and large values of this ratio, indicating that it must be within a restricted range for subduction to occur.

## 4.4 Summary: favourable conditions for the initiation of subduction

The conditions that make subduction possible can be summarized in two diagrams (Fig. 15). One diagram involves the dimensionless plate buoyancy and the distance ratio, for a fixed thickness ratio. In the other diagram, the distance ratio is fixed. These diagrams show of course that a dense oceanic plate is more likely to subduct than a buoyant one, but they emphasize that the properties of the continental extension zone are also important. A thin continental crust that extends over a small width provides a favourable setting for subduction. The diagrams also show that, all else being equal, a strong oceanic plate is more likely to subduct than a weak one, which may be due to age, through its effect on elastic thickness, or to dehydration, through its effect on elastic properties.

Finally, we show in Fig. 16 the position of the plate tip with respect to the base of the continental crust, $\Delta h$, as a function of time for several types of passive margins (i.e. with different widths $X_L$ and initial crustal thicknesses $H'_{Lo}$). All parameters are chosen within their geological ranges discussed above (Table 4). Given the uncertainty on plate buoyancy, we have assumed that it is zero, such that it does not contribute to plate bending. One sees that subduction may start at various times and various stages of extension, and may take more than 150 Ma.
<figure><img src="image_18.png" /><figcaption>**Figure 15.** Diagrams separating parameter domains with or without subduction. Left: the thickness ratio is fixed. Right: the distance ratio is fixed.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Figure 16.** Depth difference $\Delta h$ between the plate tip and the base of the continental crust (scaled to the initial crustal thickness), as a function of time for various cases. Subduction begins when $\Delta h$ becomes negative and occurs at different times depending on the continental margin properties. Parameter values are given in Table 4. The width of the extension zone plays an important role in the evolution of a continental margin.</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>**Figure 17.** Calculated flexure of a buoyant plate ($((\rho_2 - \rho)/\rho_2 = +0.006)$ deflected by highly extended continental crust and mantle. The plate is 50 km thick, has a Young's Modulus of 100 GPa and is 20 kg m$^{-3}$ less dense than the mantle. The thin light grey line stands for the sea surface and stretched continental crust appears a thick black line.</figcaption></figure>

## 4.5 Fate of buoyant lithosphere

We have shown that even buoyant oceanic lithosphere can be flexed into the mantle. In this case, however, this is not sufficient for self-sustaining subduction, which requires a phase change that increases the density of oceanic crust. If this involves eclogitization, for ex-

ample, the minimum depth is about 40 km (Peacock 2003). The calculations presented earlier were stopped before plate submer- sion into the mantle because this would introduce yet another set of variables. For the purposes of evaluating conditions that are re- quired for the plate tip to reach a certain depth, one can use simple plate bending calculations with no time-dependent spreading, how- ever. Fig. 17 illustrates one case where a buoyant oceanic plate has been pushed to a depth of 45 km. The plate is 50 km thick, has a Young's Modulus of 100 GPa and is 20 kg m$^{-3}$ less dense than the underlying mantle. The load on the plate is mostly made of mantle with highly stretched continental crust above it.

<table><caption>**Table 4.** Parameter values for Fig. 16.</caption><thead><tr><th></th><th colspan="3">Mantle</th></tr></thead><tbody><tr><td>ρ<sub>2</sub> (kg m<sup>-3</sup>)</td><td colspan="3">3300</td></tr><tr><td colspan="4">Oceanic lithosphere</td></tr><tr><td>ρ (kg m<sup>-3</sup>)</td><td colspan="3">3300</td></tr><tr><td>ν</td><td colspan="3">0.25</td></tr><tr><td>E (GPa)</td><td colspan="3">100</td></tr><tr><td colspan="4">Age-dependent parameters</td></tr><tr><td>age (Ma)</td><td>20</td><td>50</td><td>130</td></tr><tr><td>h<sub>p</sub> (km)</td><td>25</td><td>35</td><td>50</td></tr><tr><td>h<sub>o</sub> (km)</td><td>4.50</td><td>5.25</td><td>5.75</td></tr><tr><td colspan="4">Continental crust</td></tr><tr><td>ρ<sub>1</sub> (kg m<sup>-3</sup>)</td><td colspan="3">2700</td></tr><tr><td>η<sub>1</sub> (Pa s)</td><td colspan="3">5.10<sup>20</sup></td></tr></tbody></table>

## 5 DISCUSSION

Subduction initiation at passive margins has been suggested in several cases, including the North Afghanistan platform in the Devonian and the southern Iranian plate during the late Triassic (Bradley 2008), but the driving mechanism has not been identi- fied. As illustrated by our laboratory experiments and theoretical model, spreading of continental crust on oceanic basement may do the job. With all the parameters and rheological requirements that are involved, however, it may not be efficient in practice and it is important to review supporting evidence at passive margins. We first discuss physical constraints on the feasibility of this process. We, then, review some relevant features of passive margins, dealing first with structural characteristics and timing of events, and then with the late uplift that affect many margins.
## 5.1 Feasibility of the mechanism for subduction initiation: physical considerations

Feasibility depends on two factors. One is the magnitude of plate bending that can be achieved and the other is the time that it takes for crustal spreading to occur. As shown in this paper, one should not focus on bending only: What really matters is that the plate must get submerged in mantle. This can be achieved through bending, of course, but also through crustal thinning, because it reduces the depth that the plate tip must reach. The magnitude of bending depends on the volume of continental crust that can be set in motion, because it determines the magnitude of the load on the plate. Increasing the width of the continental extension zone increases the volume of crust available for spreading but it also reduces the amount of crustal thinning that can occur. The interplay between the processes of crustal spreading and thinning is subtle and our results show that, in the end, subduction initiation is favoured by relatively narrow extension zones. The time taken by crustal spreading depends on the crustal flow law for crustal rocks as well as temperature. For a given rheological law, crustal temperatures must be above a certain threshold for crustal spreading to be rapid enough. According to Nikolaeva *et al.* (2010), the threshold Moho temperature for a wet quartzite flow law is about 650 °C, close to melting conditions in hydrated crustal rocks. Our approach relies on an effective viscosity value, which lumps together the effects of flow law parameters and temperature. Thus, instead of a threshold temperature, we determine a threshold viscosity. With the viscosity value chosen, $5 \times 10^{20}$ Pa s, crustal spreading and plate bending proceed in a few tens of million years. The exact time taken for subduction initiation depends on the other parameters. For this time to be less than 200 Myr, the bulk crustal viscosity must not exceed a value of about $2 \times 10^{21}$ Pa s. For reference, viscosity values that have been proposed for lower crustal flow in intracontinental extension zones are smaller than $5 \times 10^{19}$ Pa s (Block & Royden 1990; McKenzie *et al.* 2000).

## 5.2 Crustal extension at passive margins

We rely mostly on the extensive number of surveys that have been made on the eastern Canadian margin, including Newfoundland and Nova Scotia, and on the conjugate Iberian and Irish margins. In a profile perpendicular to the continental edge, these margins can be split into three different domains. At the landward end of the profile, one finds extended continental crust that thins away from the continent, going from typical thicknesses of 30–35 km to a few kilometres (Funck *et al.* 2004). Oceanic crust lies at the other end and is associated with well-developed linear magnetic anomalies allowing dating of the seafloor. Between these two well-defined units and structures, one finds transitional basement that defies simple explanations. Such basement is made of very thin (typically ≈ 2–3 km) continental crust lying on top of serpentinized peridotite (Shillington *et al.* 2006; Gerlings *et al.* 2011). Sediments recovered by drilling indicate that such basement was emplaced in deep water (Shillington *et al.* 2008). According to Peron-Pinvidic *et al.* (2007), two major and separate phases of extension can be identified. A first phase of extension and rifting ended at about 135 Ma with the onset of seafloor spreading and was followed by a separate phase around 110 Ma. Results from ODP sites suggest that transitional basement with thin continental crust was formed around that time. Interestingly, deep sediments of the Newfoundland margin were intruded by sills dated at about 100 Ma (Peron-Pinvidic *et al.* 2010),

showing that the margin was active long after continental breakup. There is no consensus on the mechanisms that are responsible for the late extension phase. We suggest that seafloor spreading was followed by topography-driven spreading of continental crust onto young oceanic basement. One clue may be provided by the presence of serpentinized mantle below the thin continental crust of transitional basement (Shillington *et al.* 2006; Gerlings *et al.* 2011). In a simple model, extension and rifting affect continental lithosphere as a coherent unit, with continental mantle that remains beneath continental crust. This model, however, does not explain how such mantle got serpentinized, which requires the invasion of water from above. Pérez-Gussinyé & Reston (2001) have suggested that this proceeds through extension-induced faults, but it is not clear how this can generate the laterally extensive pervasive serpentinization that is documented by the seismic data. In our model, continental crust spreads over oceanic basement which was once in direct contact with sea water.

The exact sequence and nature of events that have affected the Newfoundland and Iberia margins are difficult to decipher because the evidence lies below sea level. Recent on-land observations in the Red Sea rift support the proposed scenario. There, crustal extension and magma intrusion are currently affecting the eastern passive margin of the rift, 200 km away from the active oceanic ridge and 4 Myr after the onset of seafloor spreading (Pallister *et al.* 2010).

## 5.3 Late uplift of passive margins

In the model of this paper, passive margins are not fully formed in the initial rifting phase before ocean opening and continue to evolve as seafloor spreading proceeds. One poorly understood behaviour of continental margins is uplift that occurs over a few tens of million years with an amplitude of several hundred metres. These uplift phases clearly post-date seafloor spreading and affect passive margins as the adjacent seafloor subsides. They are observed all around the North Atlantic (Rohrman & van der Beek 1996; Japsen & Chalmers 2000) as well as in the South Atlantic, on the Angola–Congo and Brazil margins (Gallagher *et al.* 1994; Lavier *et al.* 2001). Several processes have been proposed but fail to explain the amplitude of the uplift. In the North Atlantic, the Iceland hotspot has been invoked (Vågnes & Amundsen 1993; Rohrman & van der Beek 1996; Stuevold & Eldholm 1996; Clift *et al.* 1998), but this does not generate a large-enough uplift and runs into chronological inconsistencies. Further, this accounts for neither the elongated shape of the uplifted regions parallel to the coastline nor the continuing subsidence of adjacent seafloor. The hotspot hypothesis has also been proposed for the west African margins with similar difficulties (Lavier *et al.* 2001; Séranne & Anka 2005). Alternative explanations include postglacial rebound for the North Atlantic (Riis & Fjeldskaar 1992) and magmatic underplating (Gallagher *et al.* 1994), with again poor agreement with the amount of uplift that is observed. All these models rely on processes that are specific to individual areas, and yet the uplift phenomenon affects almost all the Atlantic margins. Thermal relaxation following rifting and oceanic opening does lead to a small amount of uplift, but it is not sufficient (Leroy *et al.* 2008). We now show that crustal spreading over adjacent oceanic seafloor does lead to uplift of the margin shoulders.

The model for crustal spreading led to equations for crustal thickness and horizontal velocity, which can be used to derive an equation for the vertical velocity using the continuity equation. Denoting
<figure><img src="image_21.png" /><figcaption>**Figure 18.** Time evolution of a passive margin (top), and vertical velocity at the surface of the spreading continental crust. An upward velocity is observed at the edge of the continental margin, whereas the adjacent oceanic plate subsides. Parameter values are summarized in Table 5.</figcaption></figure>

vertical velocity by, v, we obtain,

$$v(h) = v(w) + \frac{\rho_1 g}{\eta_1} (h-w)^2 \left[ \frac{1}{3} \frac{\partial^2 h}{\partial x^2} (h-w) + \frac{1}{2} \left( \frac{\partial h}{\partial x} \right)^2 - \frac{\partial h}{\partial x} \frac{\partial w}{\partial x} \right], \quad (24)$$

where $v(h)$ is vertical velocity at the top of continental crust and $v(w)$ is vertical velocity at $z = w$, that is, at the top of the oceanic plate. The other variables have been defined in Section 3. Fig. 18 shows the vertical velocity at the top of the crust as a function of distance and at different times for parameter values that are listed in Table 5. Near the continental edge, this velocity is positive, corresponding to uplift, whereas the oceanic plate is subsiding at the same time, which is consistent with the geological observations. The amount of uplift depends on the volume of crustal material that undergoes extension and can be as large as 1 km. In the particular example shown here, uplift velocities are on the order of a few hundred metres per million years, which is larger than observed. A slight increase of bulk crustal viscosity would bring the rate down to observed values. A detailed parametric study of the uplift phenomenon would be outside the scope of this study but these results do show that crustal spreading on adjacent seafloor leads to passive margin uplift.

## 5.4 Subduction in the Atlantic

We have found that a narrow continental extension region is favourable to subduction. There are currently two active subduction zones in the Atlantic, in the Caribbean and the Scotia-South

Sandwich regions (Fig. 19). Both are located where the adjacent continent is narrow and some features of the South Sandwich subduction are worth recalling. The Scotia Sea is surrounded by continental fragments that have undergone extension and that were part of a continuous connection between South America and the Antarc- tica Peninsula 40 Ma (Barker 2001). Subduction, which began at about that time, was thus initiated next to a thin continental bridge. Another interesting feature is that the trench has been retreating towards the east over the whole 40 Ma interval (Barker 2001). This has been attributed to subduction roll-back due to the increasing flexure of the oceanic plate sinking under is own weight (Royden & Husson 2006). This mechanism, however, cannot be invoked in early stages of subduction because it requires a plate that has already been thrust over at least 150 km into the mantle (McKenzie 1977). In our model, trench retreat is automatically induced even when the oceanic plate has not started to subduct in earnest (Fig. 4). It is not the oceanic plate that retreats but continental crust that overrides oceanic basement and displaces the trench location.

The other Atlantic subduction zone, which occurs beneath the Caribbean, shares many characteristics with the South Sandwich one. Both are comparable in form and dimension and are surrounded by continental fragments (Barker 2001; James 2009a). James (2009b) has suggested that the Caribbean plate is made of highly stretched continental crust lying below a basaltic layer resulting from extension-induced decompression melting. Interestingly, parts of the Caribbean basement are made of serpentinized mantle. A definitive test of the mechanism proposed could well be provided by detailed geophysical surveys of the Atlantic subduction zones and their margins.

## 6 CONCLUSIONS

Recent studies have suggested that oceanic lithosphere may never become denser than the underlying mantle (Hynes 2005; Afonso et al. 2007). Should this conclusion hold, it would make it very clear that subduction is not due to an intrinsic property of oceanic lithosphere and requires an external process. Subduction could only develop if the plate gets pushed into the mantle at depths that are sufficiently large for metamorphic reactions to kick in and change its buoyancy.

We have studied how continental extension at a passive margin can induce flexure of the oceanic plate. We show that this mechanism depends on the properties of both the oceanic lithosphere and adjacent continental crust, and draw attention to the width of the crustal extension zone. Subduction is easier for an old oceanic plate, that is, with a high average density and a large elastic thickness. A minimum rigidity is required for the plate to be flexed to sufficiently large depths into the mantle. An initially thin continental crust deforming over a region of limited width favours subduction. The size of the crustal extension region may well be the key parameter controlling the location of the first subduction zone at the margin of a new ocean. This size is determined by thermal processes that operated before the onset of seafloor spreading, which may involve hotspot activity. As a consequence, subduction may start at different ages in the same ocean basin.

Our study suggests that passive margins continue to evolve in important ways after the onset of seafloor spreading and draws attention to late uplift phases that have affected many of them. It is becoming increasingly clear that passive margins are not as passive as they were once thought to be and that post-breakup processes deserve further study.
<table><caption>**Table 5.** Main parameter values and dimensionless numbers for Fig. 18, showing margin uplift.</caption><thead><tr><th>E (GPa)</th><th>h<sub>p</sub> (km)</th><th>δ<sub>p</sub> (km)</th><th>ρ (kg m<sup>-3</sup>)</th><th>η (Pa s)</th><th>X<sub>L</sub> (km)</th><th>H'<sub>Lo</sub> (km)</th><th>h<sub>p</sub>/H'<sub>Lo</sub></th><th>X<sub>L</sub>/δ<sub>p</sub></th><th>ρ<sub>2</sub>-ρ/ρ<sub>2</sub></th></tr></thead><tbody><tr><td>100</td><td>25</td><td>95</td><td>3300</td><td>5.10<sup>21</sup></td><td>600</td><td>35</td><td>0.7</td><td>6.3</td><td>0.0</td></tr></tbody></table>

<figure><img src="image_22.png" /><figcaption>**Figure 19.** The two active subductions of the Atlantic ocean (thick red lines). Zones of widespread extension with scattered continental blocks are encased in white frames.</figcaption></figure>

## ACKNOWLEDGMENTS

We thank two anonymous reviewers for their useful comments, criticisms and suggestions which helped us improve the manuscript.

## REFERENCES

Afonso, J.C., Ranalli, G. & Fernández, M., 2007. Density structure and buoyancy of the oceanic lithosphere revisited, *Geophys. Res. Lett.*, **34**(10), doi:10.1029/2007GL029515.

Al Damegh, K., Sandvol, E. & Barazangi, M., 2005. Crustal structure of the Arabian plate: new constraints from the analysis of teleseismic receiver functions, *Earth planet. Sci. Lett.*, **231**, 177–196.

Barker, P.F., 2001. Scotia Sea regional tectonic evolution: implications for mantle flow and palaeocirculation, *Earth-Sci. Rev.*, **55**, 1–39.

Bird, P. & Piper, K., 1980. Plane-stress finite-element models of tectonic flow in southern California, *Phys. Earth planet. Inter.*, **21**, 158–175.

Block, L. & Royden, L., 1990. Core complex geometries and regional scale flow in the lower crust, *Tectonics*, **9**, 557–567.
Brace, W.F. & Kohlstedt, D.L., 1980. Limits on lithospheric stress imposed by laboratory experiments, *J. geophys. Res.*, **85**, 6248–6252.

Bradley, D.C., 2008. Passive margins through earth history, *Earth-Sci. Rev.*, **91**, 1–4.

Burov, E. & Cloetingh, S., 2010. Plume-like upper-mantle instabilities drive subduction initiation, *Geophys. Res. Lett.*, **37**, L03309, doi:10.1029/2009GL041535.

Christensen, N.I. & Mooney, W.D., 1995. Seismic velocity structure and composition of the continental crust: a global view, *J. geophys. Res.*, **100**, 9761–9788.

Clift, P.D., Carter, A. & Hurford, A.J., 1998. The erosional and uplift history of NE Atlantic passive margins: constraints on a passing plume, *J. geol. Soc. Lond.*, **155**, 787–800.

Cloetingh, S.A.P.L., Wortel, M.J.R. & Vlaar, N.J., 1982. Evolution of passive continental margins and initiation of subduction zones, *Nature*, **297**, 139–142.

Cloetingh, S., Wortel, R. & Vlaar, N.J., 1989. On the initiation of subduction zones, *Pure appl. Geophys.*, **129**(1-2), 7–25.

Davaille, A. & Jaupart, C., 1993. Transient high-Rayleigh-number thermal convection with large viscosity variations, *J. Fluid Mech.*, **253**, 141–166.

England, P. & McKenzie, D., 1982. A thin viscous sheet model for continental deformation, *Geophys. J. Int.*, **70**, 295–321.

Funck, T., Jackson, H., Louden, K., Dehler, S. & Wu, Y., 2004. Crustal structure of the northern Nova Scotia rifted continental margin (eastern Canada), *J. geophys. Res.*, **109**, B09 102, doi:10.1029/2004JB003008.

Gallagher, K., Hawkesworth, C.J. & Mantovani, M.S.M., 1994. The denudation history of the onshore continental margin of SE Brazil inferred from apatite fission track data, *J. geophys. Res.*, **99**(B9), 18 117–18 145.

Gerlings, J., Louden, K. & Jackson, H., 2011. Crustal structure of the Flemish Cap Continental Margin (eastern Canada): an analysis of a seismic refraction profile, *Geophys. J. Int.*, **185**, 30–48.

Goren, L., Aharonov, E., Mulugeta, G., Koyi, H.A. & Mart, Y., 2008. Ductile deformation of passive margins: a new mechanism for subduction initiation, *J. geophys. Res.*, **113**, doi:10.1029/2005JB004179.

Gurnis, M., Hall, C. & Lavier, L., 2004. Evolving force balance during incipient subduction, *Geochem. Geophys. Geosystems.*, **5**, 7001, doi:10.1029/2003GC000681.

Hynes, A., 2005. Buoyancy of the oceanic lithosphere and subduction initiation, *Int. Geol. Rev.*, **47**, 938–951.

James, K.H., 2009a. In situ origin of the Caribbean: discussion of data, in *The Origin and Evolution of the Caribbean Plate*, Vol. 328, pp. 77–125, eds James, K.H., Lorente, M. & Pindell, J., Geol. Soc. London, London, UK.

James, K.H., 2009b. Evolution of Middle America and the in situ Caribbean Plate model, in *The Origin and Evolution of the Caribbean Plate*, Vol. 328, pp. 127–138, eds James, K.H., Lorente, M.A. & Pindell, J.L., Geol. Soc. London, London, UK.

Japsen, P. & Chalmers, J.A., 2000. Neogene uplift and tectonics around the North Atlantic: overview, *Glob. Planet. Change*, **24**, 165–173.

Jarrard, R.D., 1986. Relations among subduction parameters, *Rev. Geophys.*, **24**, 217–284.

Kaufman, P.S. & Royden, L.H., 1994. Lower crustal flow in an extensional setting: constraints from Halloran Hills region, eastern Mojave Desert, California, *J. geophys. Res.*, **99**, 15 723–15 739.

Kemp, D.V. & Stevenson, D.J., 1996. A tensile, flexural model for the initiation of subduction, *Geophys. J. Int.*, **125**, 73–93.

King, G., Oppenheimer, D. & Amelung, F., 1994. Block versus continuum deformation in the Western United States, *Earth planet. Sci. Lett.*, **128**, 55–64.

Korenaga, J. & Karato, S. I., 2008. A new analysis of experimental data on olivine rheology, *J. geophys. Res.*, **113**, D133B-06.

Kruse, S., McNutt, M., Phipps-Morgan, J. & Royden, L., 1991. Lithospheric extension near Lake Mead, Nevada—a model for ductile flow in the lower crust, *J. geophys. Res.*, **96**, 4435–4456.

Lavier, L.L., Steckler, M.S. & Brigaud, F., 2001. Climatic and tectonic control on the Cenozoic evolution of the West African margin, *Mar. Geol.*, **178**, 63–80.

Leroy, M., Gueydan, F. & Dauteuil, O., 2008. Uplift and strength evolution of passive margins inferred from 2-D conductive modelling, *Geophys. J. Int.*, **172**, 464–476.

Levitt, D.A. & Sandwell, D.T., 1995. Lithospheric bending at subduction zones based on depth soundings and satellite gravity, *J. geophys. Res.*, **100**, 379–400.

Lévy, F., 2009. Two aspects of the influence of continents on the mantle: subduction initiation-thermal effect of continental roots, *PhD thesis*, University Paris 7 - Diderot, France.

Lévy, F. & Jaupart, C., 2009. The initiation of subduction by crustal extension at continental margins, *EOS Trans. AGU Fall Meet. Suppl.*, **90**(52), doi:10.1029/2010JB007726.

Lévy, F. & Jaupart, C., 2011. Temperature and rheological properties of the mantle beneath the North American craton from an analysis of heat flux and seismic data, *J geophys. Res.*, **116**, doi:10.1029/2007JB005100.

Mart, Y., Aharonov, E., Mulugeta, G., Ryan, W., Tentler, T. & Goren, L., 2005. Analogue modelling of the initiation of subduction, *Geophys. J. Int.*, **160**, 1081–1091.

McKenzie, D., 1977. The initiation of trenches: a finite amplitude instability, in *Island Arcs, Deep Sea Trenches and Back-Arc Basins*, Maurice Ewing Series, Vol. 1, pp. 57–61, eds Talwani, M. & Pitman, W., American Geophysical Union, Washington, DC.

McKenzie, D., Nimmo, F., Jackson, J.A., Gans, P.B. & Miller, E.L., 2000. Characteristics and consequences of flow in the lower crust, *J. geophys. Res.*, **105**, 11 029–11 046.

Mooney, W.D., Laske, G. & Guy Masters, T., 1998. CRUST 5.1: a global crustal model at 5° × 5°, *J. geophys. Res.*, **103**, 727–748.

Mueller, S. & Phillips, R.J., 1991. On the initiation of subduction, *J. geophys. Res.*, **96**, 651–665.

Nikolaeva, K., Gerya, T.V. & Marques, F., 2010. Subduction initiation at passive margins: numerical modeling, *J. geophys. Res.*, **115**, doi:10.1029/2009JB006549.

Ogawa, M., Schubert, G. & Zebib, A., 1991. Numerical simulations of three-dimensional thermal convection in a fluid with strongly temperature-dependent viscosity, *J. Fluid Mech.*, **233**, 299–328.

Pallister, J.S *et al.*, 2010. Broad accommodation of rift-related extension recorded by dyke intrusion in Saudi-Arabia, *Nature Geosci.*, **3**, 705–712.

Peacock, S.M., 2003. Thermal structure and metamorphic evolution of subducting slabs, *Geophys. Monogr.*, **138**, 7–22.

Peltier, W.R. & Drummond, R., 2008. Rheological stratification of the lithosphere: a direct inference based upon the geodetically observed pattern of the glacial isostatic adjustment of the North American continent, *Geophys. Res. Lett.*, **35**, 16 314–16 319.

Pérez-Gussinyé, M. & Reston, T.J., 2001. Rheological evolution during extension at nonvolcanic rifted margins: onset of serpentinization and development of detachments leading to continental breakup, *J. geophys. Res.*, **106**, 3961–3975.

Peron-Pinvidic, G., Manatschal, G., Minshull, T. & Sawyer, D., 2007. Tectonosedimentary evolution of the deep Iberia-Newfoundland margins: evidence for a complex breakup history, *Tectonics*, **26**, doi:10.1029/2006TC001970.

Peron-Pinvidic, G., Shillington, D.J. & Tucholke, B., 2010. Characterization of sills associated with the U reflection on the Newfoundland margin: evidence for widespread early post-rift magmatism on a magma-poor rifted margin, *Geophys. J. Int.*, **182**, 113–136.

Poudjom Djomani, Y.H., O'Reilly, S.Y., Griffin, W.L. & Morgan, P., 2001. The density structure of subcontinental lithosphere through time, *Earth planet. Sci. Lett.*, **184**(3-4), 605–621.

Regenauer-Lieb, K., Yuen, D.A. & Branlund, J., 2001. The initiation of subduction: criticality by addition of water?, *Science*, **294**, 578–581.

Riis, F. & Fjeldskaar, W., 1992. On the magnitude of the Late Tertiary and Quaternary erosion and its significance for the uplift of Scandinavia and the Barents Sea, *Nor. Pet. Soc. Spec. Publ.*, **1**, 163–185.

Rohrman, M. & van der Beek, P., 1996. Cenozoic postrift domal uplift of North Atlantic margins: an asthenospheric diapirism model, *Geology*, **24**(10), 901–904.

Royden, L.H. & Husson, L., 2006. Trench motion, slab geometry and viscous stresses in subduction systems, *Geophys. J. Int.*, **167**, 881–905.
Scheck-Wenderoth, M., Raum, T., Faleide, J., Mjelde, R. & Horsfield, B., 2007. The transition from the continent to the ocean: a deeper view on the Norwegian margin, *J. Geol. Soc.*, **164**, 855–868.

Schutt, D.L. & Lesher, C.E., 2006. Effects of melt depletion on the density and seismic velocity of garnet and spinel herzolite, *J. geophys. Res.*, **111**, doi:10.1029/2003JB002950.

Séranne, M. & Anka, Z., 2005. South Atlantic continental margins of Africa: a comparison of the tectonic vs. climate interplay on the evolution of equatorial west Africa and SW Africa margins, *J. Afr. Earth Sci.*, **43**, 283–300.

Shillington, D.J., Holbrook, W., Van Avendonk, H., Tucholke, B., Hopper, J., Louden, K., Larsen, H. & Nunes, G., 2006. Evidence for asymmetric nonvolcanic rifting and slow incipient oceanic accretion from seismic reflection data on the Newfoundland margin, *J. geophys. Res.*, **111**, 5408, doi:10.1029/2005JB003981.

Shillington, D.J., White, N., Minshull, T.A., Edwards, G.R.H., Jones, S., Edwards, R.A. & Scott, C.L., 2008. Cenozoic evolution of the eastern Black Sea: a test of depth-dependent stretching models, *Earth planet. Sci. Lett.*, **265**, 360–378.

Solomatov, V.S. & Moresi, L.-N., 2000. Scaling of time-dependent stagnant lid convection: application to small-scale convection on Earth and other terrestrial planets, *J. geophys. Res.*, **105**, 21 795–21 818.

Sonder, L.J. & England, P., 1986. Vertical averages of rheology of the continental lithosphere: relation to thin sheet parameters, *Earth planet. Sci. Lett.*, **77**, 81–90.

Stuevold, L. & Eldholm, O., 1996. Cenozoic uplift of Fennoscandia inferred from a study of the mid-Norwegian margin, *Glob. Planet. Change*, **12**, 359–386.

Tait, S. & Jaupart, C., 1989. Compositional convection in viscous melts, *Nature*, **338**, 571–574.

Thatcher, W., 1995. Microplate versus continuum descriptions of active tectonic deformation, *J. geophys. Res.*, **100**, 3885–3894.

Timoshenko, S. & Young, D., 1945. *Theory of Structures.*, 1st edn, McGraw Hill, New York, NY.

Toth, J. & Gurnis, M., 1998. Dynamics of subduction initiation at preexisting fault zones, *J. geophys. Res.*, **103**, 18 053–18 068.

Turcotte, G. & Schubert, G., 1982. *Geodynamics: Application of Continuum Physics to Geological Problems*, John Wiley, New York, NY.

Vågnes, E. & Amundsen, H.E.F., 1993. Late Cenozoic uplift and volcanism on Spitsbergen: caused by mantle convection?, *Geology*, **21**(3), 251–254.

Yuan, H. & Romanowicz, B., 2010. Lithospheric layering in the North American continent, *Nature*, **466**, 1063–1069.