Anthony Jourdon

Laetitia Le Pourhiet

Frédéric Mouthereau

Dave May

Modes of propagation of continental break-up and associated oblique rift structures

V-shaped oceanic propagators are widespread around the world. Their geometry combined with magnetic anomalies associated to their opening shows at first order that ridge propagation in the third dimension occurs by pulses. In this study we use 3D thermomechanical numerical models to show how oblique kinematic boundary conditions control both the intracontinental rift development and the oceanic ridge propagation. To do so, we apply a shortening velocity boundary condition in the direction perpendicular to the extension for "strong" and "weak" crustal rheologies. Numerical models results highlight that three ridge propagation modes can occur. For low out of plane velocities (12% to 15% of the extension rate), the ridge propagation is fast (>1.5 cm/yr) and straight. Higher shortening velocities (15% to 17%) lead to a ridge propagation by pulses alternating between fast propagation (~1.5 cm/yr) and stalling phases. Finally, for higher velocities (17% to 20%) a ridge jump propagation mode occurs, localizing a new spreading centre between 100 and 200 km far from the initial ridge. We also show that ridge propagation phases are associated with dip-slip dominated deformation while stalling phases are dominated by strike-slip deformation. These deformation regimes are marked by structures reorientation while kinematic boundary conditions remain constant. We discuss these results in term of plate tectonics reconstructions and regional geological studies. This study is part of the Orogen program, a tripartite research collaboration between CNRS, BRGM, and TOTAL SA. We thank TOTAL SA for the computing time on Pangea supercomputer.

 were made with GeomapApp and post-processing with Paraview.

Introduction

Magnetic anomalies record the past direction and rate of seafloor spreading [START_REF] Vine | Magnetic Anomalies Over Oceanic Ridges[END_REF]. As such they are our best constrain for global plate tectonic reconstructions [START_REF] Müller | Age, spreading rates, and spreading asymmetry of the world's ocean crust[END_REF]2016;[START_REF] Seton | Global continental and ocean basin reconstructions since 200Ma[END_REF]. In perfectly cylindrical condition the orientation of the first magnetic anomaly is parallel to the rifted margin of continents. However, in some areas the first magnetic anomalies intercept the rifted margins of continents at acute angle instead of being parallel to it. These remnants are precious as they allow us to quantify the rate at which continental break-up propagates along the strike of the mid-oceanic ridges [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF]. Indeed, using only basic kinematic considerations, it is clear that if the propagation occurs at the same velocity as the spreading rate, the angle between the margin and the magnetic anomaly should be roughly 45°. Faster propagation leads to the parallelisation of the magnetic anomaly with the margin of the oceanic basin matching with the 2D cylindrical model. Slower propagation produces magnetic anomalies that normalise to the rifted margin, until the limit case of a transform margin (no ridge propagation at all).

In map, these phases of sluggish ridge propagation produce V-shaped margins. V-shaped propagators are found in many places on Earth (Figure 1) and seafloor age map compilation by [START_REF] Müller | Age, spreading rates, and spreading asymmetry of the world's ocean crust[END_REF] is often sufficient to unravel the first order features of the timing of continental break-up at regional or more global scale. Examples displayed in Figure 1 show propagators separated by very long segments where magnetic anomalies are parallel to the margin. This indicates that pulsation in the rate of continental break-up propagation first introduced by [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF] in the Gulf of Aden is a major feature of divergent plate boundaries. Despite the well documented geological record of these changes of rate, and their implications for plate tectonics and kinematic reconstructions or the thermal regime of oblique margins, the dynamics of the propagation at long/geodynamic timescale remains poorly studied or understood.

Kinematic interpretations of V-shaped propagators are based on models established in the late 1970s and early 1980s, which are still in use today. The first model, called the rotational extension [START_REF] Hey | Propagating rifts on midocean ridges[END_REF][START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF][START_REF] Mondy | The role of asthenospheric flow during rift propagation and breakup[END_REF][START_REF] Zwaan | Rift propagation in rotational versus orthogonal extension: Insights from 4D analogue models[END_REF], assumes that seafloor spreading occurs near the Euler pole and is therefore described by very large variations in the opening rate in the direction of the ridge. In this case, the first magnetic anomaly forms parallel to the rifted margin then the linearly increasing rates of opening with distance from Euler pole produce a series concentric V's over time Figure 2c. This model infers that very rapid changes in tectonic regime must occur at the tip of the V-shaped ocean domain. The second kinematic model, published by [START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF], considers that the divergence remains constant in the direction of rift axis but it is variably partitioned between the oceanic accretion and the continental rifting (Figure 2b). This model applies only when magnetic anomalies indicate that the Euler Pole is far enough so that the rate of divergence is constant at regional scale. This is generally the case for a distance greater than 3000 km or 30° from the Euler pole. The Vink model allows the formation of V-shaped propagators when the propagation rate is similar to the oceanic accretion rate. For faster propagation, magnetic anomalies are parallel to the margin (Figure 2a). The Vink model is the basis for [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF] proposal (Figure 2d) who interpreted the change in angle between magnetic anomalies and the continental margins along the direction in the Gulf of Aden as the geological record of pulsating continental break-up propagation. It should be noted here that these two kinematic models are not contradictory but rather complement each other as two end-member cases. However, Vink's and Courtillot's models predictions in terms of magnetic anomalies and associated tectonic structures at the front of the propagator are different and therefore permit the description of different natural examples of V-shaped propagators. The Woodlark Basin only opens 500 km from the Euler Pole. It displays a strong segmentation, variations from extensional to strike-slip structure [START_REF] Holm | Post 8 Ma reconstruction of Papua New Guinea and Solomon Islands: Microplate tectonics in a convergent plate boundary setting[END_REF][START_REF] Mondy | The role of asthenospheric flow during rift propagation and breakup[END_REF] and only the youngest ridge ages at the tip which conform with characteristics of a rotational extension opening basin. While the South China Sea that opened at about 4500 km from its pole of rotation [START_REF] Briais | Updated interpretation of magnetic anomalies and seafloor spreading stages in the south China Sea: Implications for the Tertiary tectonics of Southeast Asia[END_REF][START_REF] Mazur | Displacement along the Red River Fault constrained by extension estimates and plate reconstructions[END_REF][START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] with very little variations in the spreading rate along the strike would better illustrate the Vink kinematic model. Dynamically speaking, on the one hand, the rotational opening model has been reproduced by a number of dynamic simulations based either on analogue materials [START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF][START_REF] Zwaan | Rift propagation in rotational versus orthogonal extension: Insights from 4D analogue models[END_REF] or numerical methods [START_REF] Mondy | The role of asthenospheric flow during rift propagation and breakup[END_REF] while, on the other hand, almost no properly scaled dynamic simulations, including numerical [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF][START_REF] Beaussier | Near-ridge initiation of intraoceanic subduction : Effects of inheritance in 3D numerical models of the Wilson Cycle[END_REF][START_REF] Liao | From continental rifting to sea fl oor spreading : Insight from 3D thermomechanical modeling[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF] or analogue models [START_REF] Ding | Propagated rifting in the Southwest Sub-basin , South China Sea : Insights from analogue modelling[END_REF][START_REF] Zwaan | Insights into the effects of oblique extension on continental rift interaction from 3D analogue and numerical models[END_REF], reproduces the kinematics proposed by [START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF] at lithospheric scale. With these boundary conditions and a realistic viscosity in the asthenospheric mantle, lithospheric necking should propagate approximately ten time faster than plates diverge [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF][START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] resulting invariably in the fast propagation kinematic displayed in Figure 2a. Therefore, to produce a kinematic solution, which resembles [START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF] model (Figure 2b) and a number of geological observations (Figure 1), some additional ingredients must therefore be added to the dynamic simulations to stop or slow down the propagation in order to form a V-shaped oceanic propagator.

The obvious candidate to slow down the propagation of continental break-up is the presence of a strong rheological heterogeneity, so called strong ribbon (Figure 2d). This hypothesis was first proposed by [START_REF] Morgan | Causes and rate-limiting mechanisms of ridge propagation: A fracture mechanics model[END_REF] based on the theory of crack propagation. It was supported by small strain numerical simulations of the propagation of continental rifting published by [START_REF] Van Wijk | Dynamics of continental rift propagation : the endmember modes[END_REF]. The results of this first study are consistent with the results presented by [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF]: the propagation rate of a neck (continental rift) in a frictional layer (upper crust) lying above a viscous layer (lower crust) decreases with the viscosity of the lower layer. Many similar analogue or numerical experiments [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF][START_REF] Brune | Controls of inherited lithospheric heterogeneity on rift linkage: Numerical and analog models of interaction between the Kenyan and Ethiopian rifts across the Turkana depression[END_REF][START_REF] Corti | Tectonophysics Evolution and characteristics of continental rifting : Analog modeling-inspired view and comparison with examples from the East African Rift System[END_REF] have proved that increasing the viscosity of the lower crust resulted in a slower lengthening of grabens in the upper crust. However, at lithosphere scale there is nothing like a high viscosity asthenosphere to delay or slow the propagation of continental break-up. At that scale, the strength of the lower crust has been shown to enhance fast lithospheric necking and tends to increase the rate of propagation of continental break-up [START_REF] Benes | Oblique rifting in the havre trough and its propagation into the continental margin of New Zealand: Comparison with analogue experiments[END_REF][START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF].

Oblique continental rifting controlled by the presence of pre-existing viscous weak zones oblique to the extension direction creates en-échelon rift basins, which are separated by the transfer zones that accommodate strain differences between segments [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Brune | Evolution of stress and fault patterns in oblique rift systems: 3-D numerical lithospheric-scale experiments from rift to breakup[END_REF][START_REF] Duclaux | Rotation , narrowing , and preferential reactivation of brittle structures during oblique rifting[END_REF][START_REF] Mart | Analogue experiments of propagation of oblique rifts[END_REF][START_REF] Mcclay | Analogue modelling of orthogonal and oblique rifting[END_REF][START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF]. These isolated basins form small segments of oceanic accretion before complete lithosphere continental plate break-up is achieved. As they grow in width they also grow in length, propagating from the basin center towards both the extremity. Yet their propagation is limited by the compressional stress build-up that arises from the growth of nearby segment. Le [START_REF] Le Calvez | Experimental designs to model along-strike fault interaction fault interaction[END_REF] proposed a model set-up to study in details the dynamic of linkage between en-échelon segments. This set-up consists of two offset weak zones which grow laterally in map view and link or not depending on spacing, frictional layer thickness and viscosity of the underlying layer. [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF] presented the result for a similar set of numerical experiments and showed the importance of the rate of softening in the frictional layer to produce a transform fault. More recently, Le [START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF] have extended this model to study continental break-up at lithosphere scale (Figure 2e). At that scale, the linkage between segments can be delayed by 10's of millions of years. Meanwhile, two asymmetric V-shaped propagators form, as shown in Figure 2e. Full linkage occurs only after oblique necking of the mantle lithosphere is reached showing similarities with the initial conditions of [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF] study in which an oblique weak zone is imposed. These asymmetric propagators better fit the observations in the Gulf of Aden or in the equatorial Atlantic (Le [START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF] than the original propagator models do (Figure 2a,b), but they do not apply to all geological examples of V-shaped propagators.

To produce oblique extension, pure shear boundary conditions in map view similar to the Alboran Sea set up in Le [START_REF] Pourhiet | Links between long-term and short-term rheology of the lithosphere: Insights from strike-slip fault modelling[END_REF] (i.e. stretching velocity in one horizontal direction and a shortening velocity in the other horizontal direction) are a good alternative to the introduction of oblique heterogeneities in cylindrical boundary conditions. Indeed, strain localisation causes the formation of structures that can rotate internally within the pure shear far field strain and which favours the emergence of instable 4 plates system that stabilizes into stable 3 plates systems [START_REF] Gerya | Nucleation and evolution of ridge-ridge-ridge triple junctions : Thermomechanical model and geometrical theory[END_REF]. It has been demonstrated that pure shear conditions produce V-shaped propagators both during the reorganisation of oceanic plates [START_REF] Gerya | Nucleation and evolution of ridge-ridge-ridge triple junctions : Thermomechanical model and geometrical theory[END_REF] and during continental break-up (Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF]. This results together with the analogue experiments of [START_REF] Ding | Propagated rifting in the Southwest Sub-basin , South China Sea : Insights from analogue modelling[END_REF], who introduced a gradient in gravitational potential energy facing the propagation direction, suggests that applying a small amount of shortening or compression out of the plane of the main driving extension slows down continental rifting propagation rate. One should note that in 3D the transition from transtensional to transpressional regime in pure shear map setting is related to change in surface of the model and that as long as the surface remain in dilatational regime no shortening related structures form (Le [START_REF] Pourhiet | Links between long-term and short-term rheology of the lithosphere: Insights from strike-slip fault modelling[END_REF]. GPS velocities vectors plotted relative to fix Arabian plate in figure 1 clearly show frontal displacement towards actively propagating oceanic ridge exist in active propagating continental rifts of Aden and Red Sea. In the first case, one might suspect the influence of the Afar plume and in the second case the vicinity of the Mediterranean subduction zone. In both cases the frontal shortening is not associated to shortening tectonic structures which could "fossilize" its existence.

Here, we therefore intend to focus our attention in this still largely unexplored direction. We adopt the same boundary conditions as Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] (Figures 2f and3) in order to analyse the sensitivity of continental break-up propagation rate to both the rate of ridge parallel shortening and lower crustal rheology. We first extract from the simulations synthetic propagation rates and expected magnetic anomalies in opening basins that corresponds to data that are widely available, relate them to large scale plate deformations. We then focus on the tectonic structures that accommodates the frontal shortening in the models in order to discuss what are key tectonic observations that might help to deciphers past kinematics and lithospheric dynamics from seismic lines and thermal dataset.

Set Up of the thermo-mechanical simulations

Modelling approach

In order to model the long term deformation of the lithosphere we use pTatin3D [START_REF] May | pTatin3D : High-Performance Methods for Long-Term Lithospheric Dynamics[END_REF][START_REF] May | A scalable , matrix-free multigrid preconditioner for finite element discretizations of heterogeneous Stokes flow[END_REF], a highly scalable, massively parallel implementation of the finite element method which employs an Arbitrary Lagrangian-Eulerian (ALE) discretization together with the material point method to solve the conservation of momentum for an incompressible fluid:

(1)

(2) where η is the non-linear effective viscosity, the strain rate tensor, P the pressure, ρ the density, g the gravity acceleration vector and v the velocity vector of the fluid. This Stokes flow is coupled with the time dependent heat conservation equation:

(3) where T is the temperature, v the velocity vector of the fluid, κ the thermal diffusivity, H the heat source (here only radioactive heat sources are considered) and Cp is the heat capacity. The material density is considered to vary with temperature such as: [START_REF] Sibuet | Geodynamics of the South China Sea[END_REF] where ρ 0 is the initial material density and α the thermal expansion coefficient.

The Stokes problem is solved using a Q 2 -P 1 discretization while the heat equation is discretized on Q 1 elements. pTatin3D uses a free surface boundary condition which can dynamically evolve due to deformation and surface processes that we model with a 2D diffusion equation discretized on Q 1 elements:

(5) where h is the altitude and k the diffusion coefficient. The equation is solved using a first order explicit time integration.

Rheological model

To simulate the long time behaviour of the lithosphere we use temperature and pressure dependent non-linear rheologies. Brittle parts of the lithosphere follow the Drucker-Prager pseudo-plastic yield criterion adapted to continuum mechanics: (6) where C is the cohesion, ϕ the friction coefficient, P the pressure and the second invariant of the strain rate tensor. To simulate the friction drop in deforming faults, we use a simple linear decrease of the friction angle from 30° to 5° with accumulated plastic strain from 0 to 1. However, laboratory experiments show that under high confining pressures (> 1GPa) brittle behaviour changes to plastic behaviour (e.g. [START_REF] Kameyama | Thermal-mechanical effects of low-temperature plasticity (the Peierls mechanism) on the deformation of a viscoelastic shear zone[END_REF][START_REF] Precigout | Strain localisation in the subcontinental mantle -a ductile alternative to the brittle mantle[END_REF]. We account for that change by limiting the Drucker-Prager yield stress to a maximum deviatoric stress of 400 MPa according to the calibration performed by [START_REF] Watremez | Buoyancy and localizing properties of continental mantle lithosphere: Insights from thermomechanical models of the eastern Gulf of Aden[END_REF].

Viscous parts of the lithosphere are simulated using an Arrhenius flow law for dislocation creep: (7) where A, n and Q are material dependant parameters (see Table 1), R is the gas constant and V is the activation volume.

Initial conditions

Our experiments represent a 1200 (in the x direction) x 600 (in the z direction) x 250 (in the y direction) km 3 domain discretized with 512x256x128 elements respectively (Figure 3a). In the y direction, two third of the 128 elements are concentrated in the first third of the 250 km giving approximately 8 km 3 (2x2x2 km 3 ) elements in the first 80 km. Four rheological layers subdivide the models domain to represent the upper and lower crust (each of 20 km thick), the lithosphere mantle (80 km thick) and the asthenosphere (170 km thick). The Moho is located at 40 km depth and the lithosphere-asthenosphere boundary is at 120 km depth (Figure 3a). Lithosphere mantle and asthenosphere are both simulated using a dry olivine flow law [START_REF] Hirth | Rheology of the Upper Mantle and the Mantle Wedge: A View from the Experimentalists[END_REF] (Table 1).

Depending on the model, the crust rheology may vary. The rheology of the lower crust has been proven to be of first order importance for the strain localization and propagation in 3D [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF][START_REF] Pourhiet | Links between long-term and short-term rheology of the lithosphere: Insights from strike-slip fault modelling[END_REF][START_REF] Brune | Controls of inherited lithospheric heterogeneity on rift linkage: Numerical and analog models of interaction between the Kenyan and Ethiopian rifts across the Turkana depression[END_REF][START_REF] Corti | Tectonophysics Evolution and characteristics of continental rifting : Analog modeling-inspired view and comparison with examples from the East African Rift System[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF]2018). Therefore, in this study we use two different crustal lithologies denoted as weak and strong (Figure 3c). The weak crust corresponds to an upper and a lower crust modelled with a quartz rheology [START_REF] Ranalli | Rheological stratification of the lithosphere[END_REF] while the strong crust corresponds to an upper crust modelled with a quartz rheology and the lower crust with an anorthite rheology [START_REF] Rybacki | Dislocation and diffusion creep of synthetic anorthite aggregates[END_REF]. We modelled those two rheological layering with the same boundary and initial conditions in order to assess the role of the lower crust on the evolution of the rifting and continental break-up. The initial temperature field is computed with a steady-state analytical solution (Turcotte & Schubert, 2002) (8) for a T y=0 = 0°C, an incoming mantle heat flux q m = 20 mW.m -2 , a radiogenic heat production H = 1.2x10 -6 W.m -3 , a characteristic radiogenic layer of y p = 40 km and a conductivity of 3.3 W.m -1 .K -1 . This analytical solution gives 610°C at 40 km (depth of the Moho, Figure 3c) and 1300°C at 120 km (depth of the LAB, Figure 3c). Then from 120 km to 250 km we compute a linear increase of temperature for an adiabatic gradient of 0.5°C/km (Figure 3c). Although this second part of the geotherm is not steady state, the cooling by diffusivity is very slow (less than 2°C/Myr for the maximum cooling rate) and it allows to keep reasonable conductivity values in the asthenosphere (3.3 W.m -1 .K -1 ).

Units

The initial radiogenic heat production is set as an exponential decay of heat production with depth according to Turcotte & Schubert (2002) as follow: [START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF] for a surface production H 0 of 1.2x10 -6 W.m -3 .

In order to initiate localisation on one side of the modelling domain, we define a weak zone in which a random amount of initial plastic strain ranging from 0 to 1 is applied to the particles located in a cube of coordinates x=550 to x=650, y=0 to y=100 and z=0 to z=100 instead of 0.0 to 0.03 background value. Since we use plastic strain softening, the initial amount of plastic strain will decrease the initial friction angle in the concerned area as:

Where ϕ is the friction angle, ϕ 0 the initial maximum friction angle (30°), ϕ ∞ the minimal friction angle reached after softening (5°), ε max and ε min are the bounds of maximum and minimum plastic strain amount between which softening is applied and ε p the actual plastic strain computed from equation (6).

Boundary conditions

The boundary conditions to solve the heat conservation equation are T=0°C at the top boundary (y=0 km), T=1440°C at the bottom boundary (y=250km) and null heat fluxes on all vertical boundaries of the domain. To solve the conservation of momentum, we impose a free surface on the top boundary, a free slip boundary conditions on the face of normal z at coordinates x, y, z=0, and Dirichlet boundary conditions on the normal component of velocity on all other vertical sides. An extension velocity v x =5 mm/yr is imposed on both faces of normal x and a shortening velocity v z is imposed on the face of normal z at coordinates x, y, z=600 (Figure 3a,b). This shortening velocity may vary between models (from 1.2 mm/yr to 2.0 mm/yr). It introduces a part of constriction in the model domain that intrinsically modifies the ratio between normal stretching and thinning as compared to a cylindrical model. In this range of shortening rates, propagation rates are indeed of the same order of magnitude (less than 1.5 cm/yr) of the plates velocity imposed in the model (1 cm/yr) and it is for this range of velocities that we expect the highest dependence of deformation pattern to the rheology of the lithosphere (Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF].

Post processing

In order to compare the numerical models results with natural examples and known geological and geophysical data, we computed an equivalent of seafloor magnetic anomalies. Magnetic anomalies measured on the seafloor are used to model seafloor spreading rates and oceanic crust ages in oceans [START_REF] Müller | Age, spreading rates, and spreading asymmetry of the world's ocean crust[END_REF]. In our models, the formation of oceanic crust by partial melting of the mantle and magma crystallization is not implemented. Therefore, we compute an age of mantle exhumation below the 800°C isotherm which can be roughly approximated as the Curie Temperature of ferromagnetic material [START_REF] Rushbrooke | On the Curie points and high temperature susceptibilities of Heisenberg model ferromagnetics[END_REF]. This age corresponds to the time in the model at which the mantle cooled beneath 800°C. In the models we use 0 Myr as the beginning and time goes forward from this initial value while geological time is counted from present (0 Ma) to past. Therefore, the age of the mantle exhumed during the latest stages of the models is greater than the age of the mantle exhumed in the earlier stages. Those ages values need to be interpreted as the greater is also the younger.

Figures 8, 10 and 12 and supplementary figures 11 and 12 display the deformation regime computed from the stress tensor orientation also known as the Regime Stress Ratio [START_REF] Brune | Evolution of stress and fault patterns in oblique rift systems: 3-D numerical lithospheric-scale experiments from rift to breakup[END_REF]Brune & Autin, 2013 and reference therein). The colour intensity is related to the strain rate second invariant value. Figures 4,5,7,8 and 9 show the plastic and/or viscous strain deformation. Those values represent respectively the cumulated finite strain computed with the equations ( 6) and ( 7). Therefore, the plastic strain displays the deformations related to frictional (or brittle) parts of the lithosphere while the viscous strain shows the deformations related to the dislocation creep flow laws in the "ductile" parts of the lithosphere.

Supplementary Figure 1 shows the β factor of the crust computed as (10) where h c t=0 is the initial thickness of the crust at the beginning of the simulation (h c = 40 km, t = 0 Myr) and h c t=n is the thickness of the crust at a given time n. Thus, the β factor value represents the thinning ratio of the crust at a given time.

Figure 6a displays the position of the tip and the ocean propagation with respect to time. The position of the tip is measured as the higher z coordinates of the mantle-continent transition along the strike of the oceanic ridge. Therefore, it represents a particular area at the frontier between a plate boundary (the oceanic ridge) and an intra-continental deformation zone (the rift). As discuss in the numerical results section, this area may evolve with time in a triple junction or remain as a plate boundary to intra-continental deformation transition. The ridge propagation (Figure 6b) is the time derivative of this position computed as: (11) where Vrp is the ridge propagation velocity, z t=n and z t=n-1 are the z coordinates at time t n and t n-1 respectively. Figure 6c shows the tectonic forces evolution through time. Fx and Fz represent the forces normal to the faces of normal x and z respectively. To compute forces we first compute the deviatoric stress tensor as:

(12)
where η is the viscosity and the strain rate tensor. We then compute the traction along faces of normal x and z respectively: ( (13b) where and are the unit vectors in x and z directions respectively. We then compute the tectonic forces per linear meter as:

(14a) and:

(14b)
where Lx and Lz are the length of the domain along x and z axis respectively.

Dynamic of continental break up propagation

Numerical models are aimed at investigating the role of both the crust rheology and the intensity of rift parallel shortening applied in our simulations. We show here that if we impose a sufficient velocity v z (perpendicular to the extension direction) the propagation rate of the continental breakup is not instantaneous and does not exceed 1.5 cm/yr. Moreover, this is independent of the rheology of the crust. For the same amount of shortening velocity applied in the direction perpendicular to the extension, the main difference between strong and weak lower crust models resides in the propagation mode. Furthermore, for large shortening velocities, weak and strong lower crust models tend to show the same behaviour forming oblique rifts (Figures 4 and5). In order to highlight the principal results of this study we first focus on the ridge propagation mode as a function of rift parallel inflow (obliquity).

Rheology of the crust and sensitivity to the shortening velocity

Our numerical modelling results show that the rift-parallel shortening velocity and the crust rheology strongly affect the propagation of the rift and its structure (Figures 4, 5 and6). With a strong lower crust, the minimum shortening velocity necessary to form a strike slip conjugate and isolate a new micro-plate is 17% of the extension velocity (v z =1.7 mm/yr, Figures 4e and5e, model 3). With a weak lower crust, the rift propagation starts to be sensitive to the shortening velocity when the latter reaches 12% of the total extension velocity (Figure 4b and5b, model 1). The influence of this third boundary condition becomes very important for 15% of the extension velocity (Figure 4d and5d, model 2). For both lower crust rheologies (strong and weak), the shortening velocity represents a threshold, across which pure extension gives place to a transtensional and strike-slip deformation regime. A strong lower crust favours the formation of three tectonic plates joining at a Ridge-Ridge-Fault triple junction because it allows for a more localised strain than weaker crust rheologies (Figures 4e,g, h and 5e, g, h). Margins are narrow and few shear zones accommodate the crustal thinning (as shown by the beta factor varying from 1 to 5 in less than 50 km, Supplementary Figure 1). In contrast, a weak lower crust rheology favours a more distributed deformation (Figure 4b, d, f and5b, d, f). Wide rifts are formed and many shear zones develop to accommodate the crustal thinning (beta factor of ~3 along 300 km from the exhumed mantle, Supplementary Figure 1). This distribution prevents the formation of three distinct plates and the effects of the out-of-plane boundary condition are visible in the rift evolution and propagation. Finally, for both weak and strong lower crusts a too large amount of shortening (20% for the strong case, Figures 4h and5h and 17% for the weak case, Figures 4f and5f) results in the formation of an oblique continental rift followed by oblique spreading.

Ridge propagation

The numerical models highlight three types of rifts propagation modes (Figure 4 and Figure 6a). The first mode is a very fast one (from 3 cm/yr to 1.5 cm/yr, Figure 6b) that results in an oceanic opening almost simultaneous along the rift axis. Models with a strong lower crust and a shortening velocity less than 17% of the total extension velocity (1.7 mm/yr) display this evolution (Figure 4 a,c and Figure 6).

The second propagation mode is characterized by a succession of fast propagation and stalling phases (Figure 6). This mode is principally visible in the weak lower crust models for shortening velocities less than 1.7 mm/yr (Figure 4b,d) . Due to the weak lower crust, a long phase of continental rifting occurs (~40 Myr) before ridge propagation takes place (Figure 6d). The propagation velocity starts at about 1.5 cm/yr before slowing down in ~10 Myr to less than 5 mm/yr (Figure 6b). This propagation velocity drop can be related to crustal thickness variations and competition between extension velocity and shortening velocity. Indeed, the initial continental rift phase efficiently thins the crust along 200 km to 300 km (in z direction) where the velocity field is principally perpendicular to the rift. Therefore, once crustal breakup occurs, the propagation along this already thinned crust is fast. However, once the rift tip reaches a less thinned crust the propagation stops or slows down to allow the continental crust domain in front of the tip to thin enough before resuming the propagation (Supplementary Figures 4, 7 and9). During the whole simulations time (~120 Myr) models were able to produce two or three stalling phases, each one lasting less than the previous one. This feature can be explained by the decreasing crustal thickness with time in front of the rift propagator.

Finally, a third propagation mode reveals a "jump" of the spreading centre from the tip of the stalling rift to a new basin. Models with strong and weak lower crust both present this propagation mode for shortening velocities larger than 1.7 mm/yr (Figure 4e, g, h and Figure 6a). As in the two previous propagation modes, the rift starts with a fast propagation (~1.5 cm/yr, Figure 6b) before stalling. But instead of propagating again from where it has stopped, a new spreading centre forms 200 km to 300 km farther (Figure 4e and Figure 6a). The main characteristic of this propagation mode resides in the deformation regime between the onset of the stalling phase and the formation of a new spreading centre. As soon as the along strike oceanic propagation stops, the intra-continental deformation at the front of the ridge switches from a dip-slip dominated regime to a strike-slip dominated. During this deformation regime, the rifting continues but localizes on an offset branch and leads to the formation of two offset oceanic segment and an ultra oblique margin (Figure 4e, g, h).

3.3.Tectonic forces evolution

Figure 6c shows the tectonic forces Fx (equation 14a) and Fz (equation 14b) evolution along boundaries of normal x and z respectively. Along the faces of normal x the Fx force corresponds to the extensional forces. Strong lower crust models display initial extensional forces (from 15 TN/m to 17 TN/m) greater than weak lower crust ones (from 10 TN/m to 12 TN/m) at the beginning of the experiments and through time. The Fx force display a small dependency to the rift parallel shortening magnitude at the beginning of the simulations. However, since the rift parallel shortening influences the strain localization, the force evolution differs between models. At first order, the Fx force decreases with the lithosphere thinning through time for all models. The force decrease is directly related to the thinning and propagation rates of the models: the faster the lithosphere thins, the faster the force decreases. Then, once the continental break-up occurred along the whole length of the modelled domain the force reaches or tend to approach a plateau representing the strength of the lithosphere.

Along the face of normal z a first order observation is that the Fz force is not dependent to the rift parallel shortening velocity but to the strength of the lithosphere. Indeed, the models with a strong lower crust display an initial state between 4 TN/m and 5 TN/m while the models with a weak lower crust show an initial state between at ~3 TN/m. During the deformation of the lithosphere, the Fz force shows few variations related to strain localization events in the models except for models with a strong lower crust and a rift parallel shortening of 1.2 mm/a to 1.5 mm/a. These two models display a fast force increase associated with the fast ridge propagation and a plateau between 8 TN/m and 9 TN/m when the ridge had propagated all along the domain. In absolute values, the Fz force is 2 to 5 times lower than the extensional force Fx.

Characteristic structures associated to propagation of break-up

In order to decipher the geological record of the different propagation modes and at the same time to gain insights on the transition from straight mode of propagation to partitioned mode we have chosen to details the evolution of models with rates of shortening 1.2 mm/a (model 1), 1.5 mm/a (model 2) and 1.7 mm/a (model 3).

4.1.Weak Lower crust, 12% of shortening velocity

The first model (model 1) is characterized by a weak lower crust and a shortening velocity v z =1.2 mm/yr (12% of the total extension velocity). Continental rifting lasts 40 Myr before continental breakup occurs at the southern boundary of the model (Figure 7a, and Figure 8a). During this initial phase, normal faulting is perpendicular to the extension direction. The weak lower crust results in strain distribution over a wide continental rift (Figure 7a, b, cross sections 2 and 3 and Figure 8a,b). At 50 Myr, the continental breakup reaches z=200km and a localized strain band replaces the wide strain distribution in the exhuming mantle between z=0 and z=200 km (Figures 7b, cross section 3 and Figure 8b), while in the rest of the model continental rifting is still active (Sections 1 and 2 of Figure 7b). At 60 Myr, the domain of exhumed mantle reaches z=300 km and the oceanic rift enters a second phase characterized by short episodes of propagation and stalling take place. Propagation first slows down for 10 Myr while the stretching of the continental crust increases in a narrow domain and finally localized as crustal breakup occurs, allowing oceanic rift propagation to continue (Figure 7c,d). During this phase, a minor strike-slip crustal deformation initiates at the edge of the propagator tip (coordinates: 300<z<400; 700<x<800, Figure 7d, cross section 2 and Figure 8d) that accommodates the differential movement between the stalling ridge and ongoing crustal extension away from the ridge. Then, from 80 Myr to 100 Myr the ridge propagation resumes and as the mantle is exhumed, the active deformation in the crust progressively vanishes (Figure 7e-g, Figure 8e-g). The age map of mantle exhumation displays a typical pattern where the oldest ages are decreasing in the direction of propagation (Figure 7g). At the end of the model, the angle between the oceanic rift and the continental margins varies along strike. This angle variation developed during the slowing down of the propagation phase. The direction of normal shear zones changes in time and along the rift axis. The first faults formed during the continental rift phase have been rotated during the first stalling phase (between 50 Myr and 60 Myr, Figure 7c,d). Then, as the rift propagates the shear zones formed in the crust continue to rotate and align with the rifted margins.

4.2.Weak Lower crust, 15% of shortening velocity

The second model (model 2) also presents a weak lower crust but with a velocity v z =1.5 mm/yr (15% of the total extension velocity). Here again continental rifting takes place during approximately 40 Myr before reaching the crustal breakup and mantle exhumation (Figure 8a and Supplementary Figures 5a and6a). Like the previous model a wide rift develops but deformation appears to be more localized. At 50 Myr, the ridge is propagating northward (up to z=200) and normal faulting is reduced (Figure 9b and Figure 10b). At this stage, the propagation starts to slow down and the ridge propagation enters a stalling phase until 80 Myr. During these 30 Myr, the exhumed mantle domain widens, forming stripes of exhumation ages parallel to the direction of extension but perpendicular to the continental margins located at the edge of the propagator (Figure 9b,c,d). This stalling phase results in drastic modification of the rift architecture. The rifted margins located around the tip and normal faults rotate until they form an angle of 90° with the active ridge (Figure 9c,d,e and Figure 10c,d,e). As a result, the rifted margins located at the propagator tip are striking parallel to the extension direction. This rotation is accompanied with the development of a left-lateral strike-slip fault zone in the crust (Figures 9c,d,e and Figure 10c,d,e). Transcurrent deformation first develops at around 50 Myrs as ridge propagation stops (Figure 10b) and evolves as a highly localized crustalscale strike-slip shear zone (Figure 10e). The rotated normal faults of the continental rift reach an obliquity of ~45° with respect to the stretching direction and ultimately form releasing bends aligned in the strike-slip fault direction. During this phase, the continental crust is thinned under a transtensional strain regime, exhuming the lower crust in elongated dome structures parallel to stretching direction as described on the field (e.g. [START_REF] Jolivet | Strain localization during crustal-scale boudinage to form extensional metamorphic domes in the Aegean Sea[END_REF] and reproduced previous modelling studies (e.g. Le [START_REF] Pourhiet | Kinematic interpretation of the 3D shapes of metamorphic core complexes[END_REF][START_REF] Pourhiet | Links between long-term and short-term rheology of the lithosphere: Insights from strike-slip fault modelling[END_REF]. Around 80 Myr, the propagation resumes and the exhumed mantle domain opens obliquely to the stretching direction along the former strike-slip zone (Figures 9e and Figure 10e). At 90 Myr, the ridge propagation slows down once again and a second left-lateral strike-slip structure develops (Figure 9f and Figure 10f). The most remarkable result is that the rift orientation change associated with the development of strike-slip deformation occurs spontaneously without involving changes in plate kinematics (Figure 9 and 10).

4.3.Strong Lower crust, 17% of shortening velocity

The third model (model 3) involves a strong lower crust and a velocity v z =1.7 mm/yr (17% of the total extension velocity). In this model, the continental rifting phase stands for less than 30 Myr before crustal breakup occurs (Figure6d). Crustal deformation is also more localized than in the previous models (Figure 11 and 12). At 30 Myr, the ridge reaches z=200 km (Figures 11a and Figure 12a) and its propagation slows down. This propagation rate change coincides with the formation of two new deformation zones (Figure 12a): a right-lateral strike-slip structure to the left and a left-lateral transtensional shear zone to the right. At 40 Myr, the ridge propagation enters in a stalling phase and stops (Figure 11b). The two branches localize now all the intra-continental deformation at the front of the ridge and a ridge-rift-fault triple junction which rapidly evolves in ridge-ridge-fault forms at the tip of the propagator (Figure 12b). The strike-slip motion in both branches intensifies and becomes predominant (Figure 12b,c) until a new exhumed mantle domain appears at 60 Myr (800<x<900; 400<z<500, Figures 11d, cross-section 1). Between the initial propagator and this new basin, a thin band of continental crust remains before the two exhumed mantle domains join at 70 Myr (Figures 11e, cross section 2). As soon as the new exhumed mantle domain forms, the deformation regime on the right branch changes once again from a strike-slip regime to transtensional extension (Figure 12e). As in the second model (weak lower crust, v z =1.5 mm/yr), the strike-slip deformation in the right branch accommodates the rotation of the continental margin previously located at the edge of the initial propagator. The rifted margin not only rotates but also thins during the transtensional deformation phase, before being torn apart when both exhumed mantle domains join (Figure 11e, cross-section 2). The finite rotation results in a rifted margin oriented parallel to the extension direction, while it was formed perpendicularly to it. This process may be responsible for the initiation of a transform margin. Here again, the lower crust exhumes preferentially in areas where strike-slip deformation is dominant. In the left branch, the right-lateral strike-slip deformation is active from 30 Myr to the end of the simulation (Figure 11 and Figure 12). This localized vertical shear zone progressively offsets the former rifted margin structures (Figure 11, cross section 3) and forms a vertical margin where the transition from continent to exhumed mantle is sharp (Figures 11e,f, cross section 4). Compared with the previous models, the lithosphere-scale strain partitioning is a main characteristic of model 3 and results and the formation of a triple junction and three distinct tectonic plates. As in the model 2, the change in deformation regime and rift orientation is provoked by the stalling phase while kinematic boundary conditions are kept constant.

Discussion

5.1.Geodynamic interpretation of rift parallel inflow

Lithospheric modelling aims at understanding the dynamics that lead to the emergence of new localised plate boundaries at regional scale. These plate boundaries are of 3 types (Fault, Ridge, Trench), which can be studied separately. But plates tectonics usually deal with triple junctions and the velocity field surrounding triple junction is not cylindrical in map, it is rather a combination of pure and simple shear components. The geodynamic pertinence of non cylindrical boundary conditions makes therefore no doubts if one wants to study the emergence of triple junction.

In term of modelling, Le [START_REF] Pourhiet | Links between long-term and short-term rheology of the lithosphere: Insights from strike-slip fault modelling[END_REF] have shown that transtensional tectonic structures (A type domes) emerging from pure shear boundary conditions are similar to those produced by imposing simple shear boundary conditions. Yet, pure shear boundary conditions have the advantage over simple shear boundary conditions because they allow the formation of large strike slip structures without having to impose a step function in the boundary conditions of the models, that is, without prescribing the location of the strike slip faults nor their occurrence at all. They are therefore useful to study the emergence or not of localised strike slip plate boundary rather than oblique rift.

Finally, as for all models run with kinematic boundary conditions, one need to discuss the forces needed to deform the model. Our results and especially the force calculation show some unexpected results which need to be discuss.

For 3D deformation patterns to emerge within cylindrical boundary conditions one must include a form of inheritance in the initial conditions. Models having large area with random noise [START_REF] Duclaux | Rotation , narrowing , and preferential reactivation of brittle structures during oblique rifting[END_REF][START_REF] Naliboff | Development of 3-D Rift Heterogeneity Through Fault Network Evolution[END_REF] demonstrates that the offset of oblique structures in cylindrical conditions is not very large without strong kinematic forcing [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Brune | Modeling suggests that oblique extension facilitates rifting and continental break-up[END_REF][START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF][START_REF] Gerya | Three-dimensional thermomechanical modeling of oceanic spreading initiation and evolution[END_REF][START_REF] Liao | From continental rifting to sea fl oor spreading : Insight from 3D thermomechanical modeling[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF] aassuming that regional deformation is cylindrical. But also through boundary conditions such as oblique extension [START_REF] Brune | Evolution of stress and fault patterns in oblique rift systems: 3-D numerical lithospheric-scale experiments from rift to breakup[END_REF][START_REF] Brune | Modeling suggests that oblique extension facilitates rifting and continental break-up[END_REF][START_REF] Brune | The rift to break-up evolution of the Gulf of Aden: Insights from 3D numerical lithospheric-scale modelling[END_REF], pure shear (e.g. [START_REF] Gerya | Nucleation and evolution of ridge-ridge-ridge triple junctions : Thermomechanical model and geometrical theory[END_REF][START_REF] Pourhiet | Links between long-term and short-term rheology of the lithosphere: Insights from strike-slip fault modelling[END_REF] or combination between pure and simple shear (Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF]. Kinematic boundary conditions involving pure shear showed the ability to reproduce formation of more than two tectonic plates and triple junctions [START_REF] Gerya | Nucleation and evolution of ridge-ridge-ridge triple junctions : Thermomechanical model and geometrical theory[END_REF].

In this study we used Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] kinematic boundary conditions to impose obliquity. As for pure shear boundary conditions we showed that accounting for sufficient rift parallel shortening leads to the formation of a triple junction between three tectonic plates in models with strong lower crust. The coupled rheology of the crust seems to be a first order constraint to produce large scale strike-slip fault in extensional contexts as also shown in analogue models [START_REF] Zwaan | A systematic comparison of experimental set-ups for modelling extensional tectonics[END_REF].

In terms of tectonic forces, the models show that the forces are more dependants to the rheology of the lithosphere than to the rift parallel shortening magnitude. In a weak lower crust lithosphere, the force associated to the rift parallel shortening (Fz) is in the same order of magnitude than the force generated by gravitational potential energy [START_REF] Bird | Patterns of stress and strain rate in southern Africa[END_REF]. For a strong lower crust lithosphere, the gravitational potential energy may not be sufficient to explain such rift parallel shortening and an active tectonic process like a subduction or a collision may be required in the vicinity of the rift to allow rift parallel shortening. However, the force generated by this rift parallel shortening is ten times lower than the extensional force and likely represents the resistance of the strike-slip fault developing in strong lower crust models which is several hundreds of kilometres long. Such a force may be lowered by a more efficient strain softening in the brittle crust.

5.2.Ridge migration and comparison with previous modelling work

Models presented in this study show that a small amount (12% to 20%) of rift parallel shortening leads to alternation of ridge propagation phases and ridge stalling phases. These alternating phases are each associated with strain regime changes. During ridge propagation, the deformation in the continental crust located at the front of the propagator is dip-slip dominated, while during stalling phases the deformation in the continental crust is strike-slip dominated.

According to our models results, the ridge propagation occurs under three different mechanisms (Figure 6 and Figure 14). For low rift-parallel shortening velocities (12% to 15% of the extension rate, Figure 4a,b,c) the ridge propagation is fast (between 1.5 cm/yr and 3.5 cm/yr, Figure 6b) and straight. Therefore, in this case the ridge propagation velocity exceeds the total extension rate (1 cm/yr). This propagation mode has been documented in crustal scale models involving a low viscosity layer beneath a plastic layer [START_REF] Allken | Factors controlling the mode of rift interaction in brittle-ductile coupled systems: A 3D numerical study[END_REF] where the lower viscosities lead to faster propagation. The same behaviour is observed in lithospheric scale models with cylindrical boundary conditions ( Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] where free slip boundary conditions prevent deformation in the rift parallel direction (which is equivalent to a boundary opposing no resistance to propagation) but with the mantle acting as the viscous layer controlling the propagation rate.

Higher rift-parallel shortening velocities (15% to 17%, Figure 4d) lead to ridge propagation by pulses alternating between fast propagation (between 0.5 cm/yr and 1.5 cm/yr, Figure 6b) and stalling phases. This propagation mode obviously requires that the ridge propagation stalls at some point. Until now, this behaviour has been modelled only during soft linkage of en-échelon propagators [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Liao | From continental rifting to sea fl oor spreading : Insight from 3D thermomechanical modeling[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF] and a rift parallel inflow boundary condition model (Le Pourhiet et al., 2018, and this study). The common feature between these models is that during the stalling phase the velocity field at the tip of propagator changes and leads to a strikeslip dominated deformation. This mode of propagation was also proposed by [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF] in the Gulf of Aden to explain the obliquity between seafloor magnetic anomalies and the passive margin. In [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF] model, a hard node in the lithosphere structure would explain why the ridge propagation stopped. All these thermo-mechanical and conceptual models show that to enter in a stalling ridge phase, a non-cylindrical parameter is required. It can be rheological [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF] or kinematic [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Liao | From continental rifting to sea fl oor spreading : Insight from 3D thermomechanical modeling[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF], 2018).

Finally, for higher velocities (17% to 20%, Figure 4e,f,g) a ridge jump propagation mode occurs. This ridge propagation mode also requires a ridge stalling phase, but rather than resuming from the ridge tip, a new spreading centre localizes between 100 and 200 km away from the ridge (Figure 6a and Figure 11d). This propagation mode is particularly well observed in models evolving towards three distinct tectonic plates (Figure 4e, g, h), therefore for highly oblique kinematic conditions. A propagation by jump was also obtained in analogue models with rotational extension coupled with an oblique weak seed [START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF] where only highly oblique week zones lead to jumps of spreading centres. Effectively, a weak zone oblique to the ridge propagation direction has the same effect as a strong zone or an inflow in the rift parallel direction. In the first case it will guide the deformation, while in the second case it will force the deformation to propagate elsewhere.

These different propagation modes are highly dependent of the lithosphere ability to deform in the rift parallel direction, whether it is due to rheological heterogeneities (e.g. [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF][START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF] or velocity field variations (e.g. [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF], 2018).

5.3.Comparison with natural cases

Although the four oceanic basins presented in Figure 1 are not at same time and space scales, they all share similar features between each other and with our models. All of them are well documented examples of ridge propagation characterized by pulses of alternating fast propagation along straight margins and stalling phases which often corresponds to the location of oblique and transform margins. Our numerical simulations show that during propagation phases the deformation occurs orthogonal to the rift axis whereas oblique extension and strike-slip deformation characterises stalling phases. Initiation of strike-slip deformation results in the reorientation of pre-existing rift normal faults. As stretching continues an oblique continental rift forms which may result in the formation of a transform margin or oblique continental break-up propagation (Figure 14). A first order observation is that as in the natural systems described in this study, none of the models involving rift parallel shortening produce shortening structures (as thrusts or folds). However, the non-cylindrical component of the velocity field is well expressed in the GPS velocities of active oblique rift systems (Gulf of Aden and Northern Red Sea basins).

Among the four oceanic basins in Figure 1, the South China Sea is probably the most compelling example of secondary oblique continental break-up propagation following a stalling phase. Figure 1b outlines three V-shaped propagators which indicate that continental break-up propagated from east to west by pulses of rapid propagation, separated by two periods of stalling. The longest period occurred from 32 to 23 Ma where magnetic anomalies in the propagator region form at a 60° angle with the continental margin indicating that the propagation rate is much lower than the spreading rate. Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] proposed that rift parallel inflow could explain this first phase of stalling using a model similar to the one with a weak crust and 1.2 mm/yr rift parallel inflow. They then posit that this long stalling period ends when both the spreading and propagating direction of the mid-oceanic ridge changed by 15° [START_REF] Sibuet | Geodynamics of the South China Sea[END_REF] in the South China Sea (Figure 1b). The origin of this rift parallel inflow has been proposed as the expression of gravitational potential energy gradient (e.g. [START_REF] Ding | Propagated rifting in the Southwest Sub-basin , South China Sea : Insights from analogue modelling[END_REF][START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] between the developing basin and the higher Indochina to the west. [START_REF] Bird | Patterns of stress and strain rate in southern Africa[END_REF] showed that a topographic load can propagate vertical forces of 1 TN/m to 2 TN/m over hundreds of kilometres. Here, with a more thorough exploration of the effect of shortening on the propagation of continental break up, we show that increasing the shortening rate to 1.5 mm/yr reproduces the 15° rotation of the spreading direction reported by [START_REF] Sibuet | Geodynamics of the South China Sea[END_REF] and the propagation of a younger oceanic propagator self consistently. Moreover, associated to the strike-slip deformation, lower crust domes exhume in the continental margin as in the South China Sea. This finding questions whether magnetic anomalies do record a change in far field conditions or simply the localization of the deformation on an oblique continental rift branch. This result may have significant implications on plates reconstruction and geodynamic interpretations. Indeed, changes of structural directions are common during continental rifting. These phases coincide with a distribution of the deformation in the crust in front of the oceanic propagator as in the Atlantic realm [START_REF] Enachescu | Structural Setting and Petroleum Potential of the Orphan Basin , offshore Newfoundland and Labrador[END_REF][START_REF] Mohn | Structural and stratigraphic evolution of the Iberia and Newfoundland hyper-extended rifted margins: A quantitative modeling approach[END_REF][START_REF] Nirrengarten | Kinematic Evolution of the Southern North Atlantic: Implications for the Formation of Hyperextended Rift Systems[END_REF][START_REF] Sibuet | Early Cretaceous motion of Flemish Cap with respect to North America : implications on the formation of Orphan Basin and SE Flemish Cap -Galicia Bank conjugate margins[END_REF], in the Gulf of Aden [START_REF] Fournier | Arabia -Somalia plate kinematics , evolution of the Aden -Owen -Carlsberg triple junction , and opening of the Gulf of Aden[END_REF][START_REF] Leroy | From rifting to oceanic spreading in the Gulf of Aden : A synthesis From rifting to oceanic spreading in the Gulf of Aden : a synthesis[END_REF] or in the northern Red Sea [START_REF] Almalki | Incipient seafloor spreading segments: Insights from the Red Sea[END_REF][START_REF] Bartov | Sinistral movement along the Gulf of Aqaba -its age and relation to the opening of the Red Sea[END_REF][START_REF] Bosworth | Geological Evolution of the Red Sea : Historical Background, Review, and Synthesis[END_REF]2020) where major transform fault zones are segmenting the oceanic basins at locations where stalling phases occurred.

Previous studies proposed that these stalling phases can be related to several processes acting as an opposing force to the rift propagation in the undeformed continental lithosphere like the presence of high strength continental domains or inherited structures in the lithosphere [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF][START_REF] Morgan | Causes and rate-limiting mechanisms of ridge propagation: A fracture mechanics model[END_REF][START_REF] Nirrengarten | Kinematic Evolution of the Southern North Atlantic: Implications for the Formation of Hyperextended Rift Systems[END_REF]. As a result they are generally interpreted as indicating varying kinematic conditions or inheritance. This hypothesis was already questioned by numerical model of rift linkage at oblique continental margin which showed that with increasing strain localisation the tectonic regime changes from initially dip-slip dominated to strike-slip dominated ( Le [START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF]. Here, we show that strain partitioning and change in deformation regime can occur at the tip of oceanic propagators under constant boundary kinematic conditions (Figures 9 and 11). This implies that a relationship seems to exist between major transform margins and the rift stalling phase during ocean opening. In the North Atlantic Ocean we observe that the Newfoundland-Azores-Gibraltar fault zone and the Charlie-Gibbs fault are the two major transform faults separate segment of different spreading ages (Nirrengarten et al., 2018 and references therein). In the same way, the Gulf of Aden (Figure 1a) opened almost synchronously between the Alula-Fartak fault zone and the Shukra el Sheik fault zone but was stalling along the Shukra el Sheik fault zone for 9 Ma to 12 Ma before propagation restarted [START_REF] Fournier | Arabia -Somalia plate kinematics , evolution of the Aden -Owen -Carlsberg triple junction , and opening of the Gulf of Aden[END_REF][START_REF] Leroy | From rifting to oceanic spreading in the Gulf of Aden : A synthesis From rifting to oceanic spreading in the Gulf of Aden : a synthesis[END_REF]. However, in our models we were not able to produce oceanic transform faults segmenting the ridge. Thermo-mechanical models tend to show that magmatism and serpentinization are essential to produce oceanic transform faults (e.g. [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Gerya | Origin and models of oceanic transform faults[END_REF].

The South China Sea and the Gulf of Aden are both good examples to illustrate the rifting of a weak continental lithosphere (Figure 4b, d and Figures 7 and9) due to the geodynamic context in which those two oceanic basins have formed. The Gulf of Aden is located near an important mantle thermal anomaly (the Afar hot spot, e.g. [START_REF] Fournier | Arabia -Somalia plate kinematics , evolution of the Aden -Owen -Carlsberg triple junction , and opening of the Gulf of Aden[END_REF][START_REF] Leroy | From rifting to oceanic spreading in the Gulf of Aden : A synthesis From rifting to oceanic spreading in the Gulf of Aden : a synthesis[END_REF] while the South China sea opening followed a long-term subduction associated with magmatism [START_REF] Pubellier | The building of pericratonic mountain ranges : structural and kinematic constraints applied to GIS-based reconstructions of SE Asia[END_REF][START_REF] Wang | Phanerozoic tectonics of the South China Block: Key observations and controversies[END_REF]. The Gulf of Aden presents the advantage of being an active plate boundary monitored with GPS stations. The GPS velocities with respect to Arabia fixed [START_REF] Arrajehi | Geodetic constraints on present-day motion of the Arabian Plate: Implications for Red Sea and Gulf of Aden rifting[END_REF] display the Eastward motion (i.e. parallel to the ridge) of the continent at the tip of the oceanic propagator and an oblique North-eastward motion for Africa fixed. The oceanic basin and continental margins display stalling phases and propagation phases and the rifting phase is characterized as a wide oblique rift. Figure 13c shows the relative velocity field in the Model 1 (1.2 mm/yr rift parallel inflow) which display the oblique velocity at the tip of the propagator. In the Gulf of Aden, the force that triggers the rift parallel inflow as displayed by GPS velocities could come from the differential gravitational potential energy between the high topography above the Afar hot spot and the opening Gulf of Aden.

However, the strain localization and the evolution of the deformation in a stronger lithosphere (Figure 4e, g, h and Figure 11) results in the development of a localized strike-slip fault and the formation of a micro-plate. These models display similarities with the Northern Red Sea and its two branches, the Gulf of Aqaba to the East and the Gulf of Suez to the West [START_REF] Almalki | Incipient seafloor spreading segments: Insights from the Red Sea[END_REF][START_REF] Bartov | Sinistral movement along the Gulf of Aqaba -its age and relation to the opening of the Red Sea[END_REF][START_REF] Bosworth | Geological Evolution of the Red Sea : Historical Background, Review, and Synthesis[END_REF]2020). In the Northern Red Sea, the Sinai block behaves as a third tectonic micro-plate separated from Arabia by the Levant transform Fault and from Africa by an active oblique rift (with respect to the Red Sea). Although a continental rift is not a complete plate boundary it is the initiation. The GPS velocities with respect to Arabia display the southward motion of the Sinai block along the Levant Fault. Therefore, regarding the absence of shortening structures in the northern Red Sea, the GPS velocities showing the rift parallel displacement and the presence of the transform continental fault of the Levant, we propose that a rotational extension model [START_REF] Molnar | Interactions between propagating rotational rifts and linear rheological heterogeneities: Insights from three-dimensional laboratory experiments[END_REF][START_REF] Smith | Infinitesimal kinematics of rotational rifting with reference to en echelon marginal faults in the Red Sea region[END_REF][START_REF] Zwaan | Rift propagation in rotational versus orthogonal extension: Insights from 4D analogue models[END_REF] may not be the best candidate to explain the dynamics of the Northern Red Sea. However, the formation of the triple junction between a continental transform fault, an active oblique rift and a ridge is well reproduced in our strong lower crust models when the rift parallel inflow is greater than 16% of extension velocity (Figure 14b, model 3). The geodynamic origin of this regional rift parallel inflow arise from the contrast between the fast motion of the Arabian plate driven by subduction in Zagros and Makran and the slow northward motion of Africa in collision with the European plate.

Moreover, the lithosphere-scale strain partitioning between a highly localized strike-slip fault and a very oblique rift leads to the development of two types of transform margins with extremely different kinematic, dynamics and thermal history.

The classical model for transform margins formation involves the development of transform faults during the intra-continental rifting phase [START_REF] Basile | Tectonophysics Transform continental margins -part 1 : Concepts and models[END_REF][START_REF] Francheteau | Marginal Fracture Zones as Structural Framework of Continental Margins in South Atlantic Ocean[END_REF][START_REF] Mascle | Evidence for transform margin evolution from the Ivory Coast-Ghana continental margin[END_REF]. The evolution from rifting to oceanic spreading should preserve the continental structures and develop oceanic transform faults along the previous continental ones. This model has major implications on development of the margin structure and temperature. It implies that oceanic ridge migrates along the transform margin as a triple junction at half the spreading rate (for large offset transform margin it may take several million years) heating the continental margin and its basins (Ruepke et al. 2010). It also implies that the strike-slip faults responsible for the steep structure of the transform margins (Mercier de [START_REF] Mercier De Lépinay | Transform continental margins -Part 2 : A worldwide review[END_REF] are active at the beginning of the rifting phase and that the transform margin was born with characteristics observable today on the fossil transform margins. However, models of this study show different characteristic and evolution for transform margins.

On the one hand, models with high obliquity display the typical kinematic context of a transform margin [START_REF] Basile | Tectonophysics Transform continental margins -part 1 : Concepts and models[END_REF][START_REF] Francheteau | Marginal Fracture Zones as Structural Framework of Continental Margins in South Atlantic Ocean[END_REF]. In this condition one of the oblique continental rift branch (Figure 11, cross sections 1 and 2, Figure 12) reach continental break-up and the system evolves into two en-échellon ridges separated by a zone of soft linkage that evolves into an oblique rift with strike-slip deformation similarly to previous studies [START_REF] Ammann | Oblique continental rifting and long transform fault formation based on 3D thermomechanical numerical modeling[END_REF][START_REF] Liao | From continental rifting to sea fl oor spreading : Insight from 3D thermomechanical modeling[END_REF][START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF]. However, structural properties of these modelled margins correspond to hyper-thinned continental crust but not steep margins along a localized transform fault. On the other hand, simulations with high obliquity and strong lower crust produce steep margins along a pure strike-slip fault (e.g. Figure 11f, corss section 4) that display the structural characteristics of a transform margin as a slope break, a marginal ridge close to the slope break and a strike-slip fault separating the continental crust from the exhumed mantle (or oceanic crust)(Mercier de [START_REF] Mercier De Lépinay | Transform continental margins -Part 2 : A worldwide review[END_REF]. However, these strike-slip margins form in atypical kinematic context for a transform margin. From an initially hyper thinned margin developed during the initial extension phase (until 30Myr), a steep margin develops (Figure 11e, f, cross section 4) as the deformation increases along the strike-slip fault. This strike-slip shear zone puts in contact a 50 km thick continental lithosphere with the asthenosphere. As a result, the formation these steep margins starts 10 Myr after the continental break-up meaning that the margin geometry is acquired after the intra-continental rifting phase. These results highlight that oblique and transform margins can acquire their final structure several million years after continental break-up due to the strike-slip deformation.

Conclusion 3D thermo-mechanical numerical models integrating a rift parallel inflow show that:

-The 3D boundary conditions can be interpreted as an absolute velocity field that can be converted to a relative velocity field comparable with natural cases of triple junctions involving continental transform faults. -At least three propagation modes can exist. A fast propagation resulting in a straight oceanic basin (2 plates), an alternation of propagation and stalling phases (extension opposing velocity of 12% to 15% for weak lower crust) and propagation with a "jump" of the spreading centre (opposing extension velocity over 17% for strong and weak lower crust). -.

-Continental break-up propagation slows down when rift parallel inflow is applied [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF] and we show that during these stalling phases the continental margins and tectonic structures rotate and the deformation regime changes from dip-slip dominated to strike-slip dominated while boundary kinematic conditions are kept constant. -There might be a relationship between stalling phases, transform margin formation and the subsequent location of large oceanic transform fault zones due to margins rotation -Two modes of transform margin can develop, steep margin associated with strike slip faults and hyper-extended ones associated with offset spreading segment similar to Le Pourhiet (2017) -All oblique margins might not emerge from a single geodynamic process and that it could be worthwile in the future to revisit the available data using the soft linkage oblique margin and the strike-slip oblique margin models as a classification

Figure captions

Figure 1: Topographic and bathymetric maps (GMRT, [START_REF] Ryan | Global Multi-Resolution Topography (GMRT) synthesis data set[END_REF] with oceanic crust ages of a) the Gulf of Aden [START_REF] Fournier | Arabia -Somalia plate kinematics , evolution of the Aden -Owen -Carlsberg triple junction , and opening of the Gulf of Aden[END_REF][START_REF] Leroy | From rifting to oceanic spreading in the Gulf of Aden : A synthesis From rifting to oceanic spreading in the Gulf of Aden : a synthesis[END_REF], c) the South China Sea [START_REF] Li | Ages and magnetic structures of the South China Sea constrained by deep tow magnetic surveys and IODP Expedition 349[END_REF][START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF], e) the Red Sea [START_REF] Almalki | Incipient seafloor spreading segments: Insights from the Red Sea[END_REF][START_REF] Bosworth | Geology of Egypt: The Northern Red Sea[END_REF], g) the North Atlantic Ocean [START_REF] Klitgord | Plate kinematics of the central Atlantic[END_REF]. Main transform fault zones are reported with thick black lines. GPS velocities are from [START_REF] Arrajehi | Geodetic constraints on present-day motion of the Arabian Plate: Implications for Red Sea and Gulf of Aden rifting[END_REF][START_REF] Reilinger | Nubia-Arabia-Eurasia plate motions and the dynamics of Mediterranean and Middle East tectonics[END_REF]. Graphs with red curves represent the ridge propagation in b) the Gulf of Aden, d) the South China Sea, f) the Red Sea, h) the North Atlantic Ocean interpreted from the oceanic seafloor ages distribution along strike. Topographic and bathymetric maps made with GeoMapApp (www.geomapapp.org) / CC BY. [START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF], b) slow ridge propagation [START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF], c) Scissor opening [START_REF] Hey | Propagating rifts on midocean ridges[END_REF]. Dynamic models of the formation of oceanic V-shaped propagators and continental rifts d) strong ribbon in the continental crust slowing ridge propagation [START_REF] Courtillot | Propagating rifts and continental breakup[END_REF], e) offset propagators [START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF], f) shortening perpendicular to extension direction ( Le [START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF].

Figure 3: a) Numerical models setup. Black arrows display the initial velocity field in the models domain. b) Velocity boundary conditions to solve the conservation of momentum. c) Yield stress envelopes of the models domain computed for a constant strain rate of 10 -14 s -1 and an initial geotherm represented by the red curve.

Figure 4: Map view of nine models with a strong and weak lower crust and varying v z velocities displaying the plastic strain and the mantle exhumation age. The plastic strain represents all the deformations that occurred under the regime described by the equation 6. The mantle exhumation age is set when the mantle cools beneath the isotherm 800°C.

Figure 5: Serial cross-sections of nine models with a strong and weak lower crust and varying v z velocities representing the plastic strain and the viscous strain. For every models, the cross-sections are regularly spaced every 50 km from z=150 km to z=500 km. The plastic strain represents all the deformations that occurred under the regime described by the equation 6. The viscous strain represents all the deformations that occurred under the regime described by the equation 7. Figure 7: Numerical model with a weak lower crust and a v z shortening velocity of 1.2 mm/yr (model 1). Left panels show a map view evolution of the plastic strain (deformation related to equation 6) in the crust and exhumation ages of the mantle below the 800°C isotherm. Red dashed lines shows the locations of cross-sections displayed in the right panels. On the cross-sections, background colours correspond to the simulated lithologies presented in Figure 3. The deformation drawing is interpreted from Supplementary Figures 2 and3.

Figure 8: Numerical model with a weak lower crust and a v z shortening velocity of 1.2 mm/yr (model 1). The background colours represent the lithologies detailed in Figure 3. The superimposed colours represent the stress inferred deformation regime [START_REF] Brune | Evolution of stress and fault patterns in oblique rift systems: 3-D numerical lithospheric-scale experiments from rift to breakup[END_REF][START_REF] Brune | The rift to break-up evolution of the Gulf of Aden: Insights from 3D numerical lithospheric-scale modelling[END_REF]. The intensity of the colourscale is related to the intensity of the strain rate tensor second invariant.

Figure 9: Numerical model with a weak lower crust and a v z shortening velocity of 1.5 mm/yr (model 2). Left panels show a map view evolution of the plastic strain (deformation related to equation 6) in the crust and exhumation ages of the mantle below the 800°C isotherm. Red dashed lines shows the locations of cross-sections displayed in the right panels. On the cross-sections, background colours correspond to the simulated lithologies presented in Figure 3. The deformation drawing is interpreted from Supplementary Figures 5 and6.

Figure 10: Numerical model with a weak lower crust and a v z shortening velocity of 1.5 mm/yr (model 2). The background colours represent the lithologies detailed in Figure 3. The superimposed colours represent the stress inferred deformation regime [START_REF] Brune | Evolution of stress and fault patterns in oblique rift systems: 3-D numerical lithospheric-scale experiments from rift to breakup[END_REF][START_REF] Brune | The rift to break-up evolution of the Gulf of Aden: Insights from 3D numerical lithospheric-scale modelling[END_REF]. The intensity of the colourscale is related to the intensity of the strain rate tensor second invariant.

Figure 11: Numerical model with a strong lower crust and a v z shortening velocity of 1.7 mm/yr (model 3). Left panels show a map view evolution of the plastic strain (deformation related to equation 6) in the crust and exhumation ages of the mantle below the 800°C isotherm. Red dashed lines shows the locations of cross-sections displayed in the right panels. On the cross-sections, background colours correspond to the simulated lithologies presented in Figure 3. The deformation drawing is interpreted from Supplementary Figures 8 and9.

Figure 12: Numerical model with a strong lower crust and a v z shortening velocity of 1.7 mm/yr (model 3). The background colours represent the lithologies detailed in Figure 3. The superimposed colours represent the stress inferred deformation regime [START_REF] Brune | Evolution of stress and fault patterns in oblique rift systems: 3-D numerical lithospheric-scale experiments from rift to breakup[END_REF][START_REF] Brune | The rift to break-up evolution of the Gulf of Aden: Insights from 3D numerical lithospheric-scale modelling[END_REF]. The intensity of the colourscale is related to the intensity of the strain rate tensor second invariant.

Figure 2 :

 2 Figure 2: Kinematic models of ridge propagation and continental rift geometry. a) Fast ridge propagation[START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF], b) slow ridge propagation[START_REF] Vink | Continental rifting and the implications for plate tectonic reconstructions[END_REF], c) Scissor opening[START_REF] Hey | Propagating rifts on midocean ridges[END_REF]. Dynamic models of the formation of oceanic V-shaped propagators and continental rifts d) strong ribbon in the continental crust slowing ridge propagation[START_REF] Courtillot | Propagating rifts and continental breakup[END_REF], e) offset propagators[START_REF] Pourhiet | A genetic link between transform and hyper-extended margins[END_REF], f) shortening perpendicular to extension direction ( Le[START_REF] Pourhiet | Continental break-up of the South China Sea stalled by far-field compression[END_REF].

Figure 6 :

 6 Figure 6: a) Curves for every models representing the position of the tip along the z axis (in km) with respect to time (in Myr). b) Curves for every models representing the velocity of the rift propagation (in mm/yr) with respect to time (in Myr). The velocity of the rift propagation is computed as the time derivative of the position of the tip. c) Evolution of the tectonic forces through time. Fx (left panel) corresponds to the forces applied in the direction normal to the face of normal x (equation 14a) and Fz

Figure 13 :

 13 Figure 13: Map view of the horizontal velocity field extracted from model 3 (v z = 1.7 mm/a, strong lower crust) at 50 Myr. a) Absolute horizontal velocity field. b) Relative horizontal velocity field (plate A fixed). c) Relative horizontal velocity field (plate C fixed). d) Relative horizontal velocity field (plate B fixed) Figure 14: Synthetic sketch representing the three propagation modes highlighted by models of this study. Propagation modes are related to the degree of obliquity and strain partitioning. Relative velocities vector directions are extracted from the model with a strong lower crust and a shortening velocity v z =1.7 mm/yr.

Acknowledgment

Strong / Coupled

Lower Crust Weak / Decoupled Lower Crust

Strong -Coupled

Weak -Decoupled

Vz (mm/yr)

Strong -Coupled

Weak -Decoupled

Strong -Coupled

Weak -Decoupled

Continental