Valentin Mikhailov
email: valentin@ipgp.jussieu.fr

Randell Stephenson
email: r.stephenson@abdn.ac.uk

Michel Diament
email: diament@ipgp.jussieu.fr

Schmidt

MODELLING OF COMPRESSION AND EXTENSION OF THE CONTINENTAL LITHOSPHERE: TOWARDS REHABILITATION OF THE NECKING-LEVEL MODEL

Keywords: lithosphere extension/compression, necking level, yield strength envelope, isostasy

We present a dynamic model of continental lithosphere deformation under extension or compression, focusing on the role of an effective mechanical parameter called -necking level‖ or -necking depth‖, a widely used concept in basin modelling studies. Though it has generally been assumed that -necking depth‖ depends strongly upon the rheological structure of the lithosphere (especially the depth distribution of its strong layers), such a dependency has never been demonstrated. Our model, which accommodates small deformations of a thin inhomogeneous plate induced by in-plane as well as by mantle boundary forces (applied to the model sides and base, respectively), shows that -necking depth‖ is a function of the horizontal position and depends mainly on the relative thicknesses and strengths of the rigid layers in the uppermost crust and below the Moho. Using different yield strength envelopes we demonstrate that the final structure of the lithosphere formed as a result of deformation and its consequent isostatic adjustment can be closely approximated by a model with a flat necking level. In the process of extension and compression of the continental lithosphere all boundaries, including the topographic surface and the Moho, deform. As a result, the total disturbance of the isostatic equilibrium state (specified as a load) is only a part of the topographic weight. Estimates of the correct load can be made using the depth to the necking level inferred from lithosphere structure, composition and thermal state. The final topography of lithospheric interfaces depends on both necking depth and effective flexural rigidity. Any attempt to estimate simultaneously strain distribution, necking depth and effective flexural rigidity, however, represents an ill-posed problem and is not possible without reliance upon strong independent assumptions constraining lithosphere structure.

Introduction

The topography of continentsboth uplifted orogens and downwarped sedimentary basins results from a series of processes: tectonic deformations, erosion, sedimentation, isostatic response of the lithosphere, etc. and, naturally, there is a long history of their quantitative modelling. Numerical models typically incorporate complexities such as the rheological effects of pressure, temperature and strain rate, compositional variations, inherited structures and so forth. In such models the lithosphere is considered as a rheologically layered entity with lateral heterogeneity [e.g. [START_REF] Mitrovica | Tilting of the continental interiors by the dynamical effects of subduction[END_REF]Braun and Beaumont, 1989a;[START_REF] Bassi | Factors controlling the style of continental rifting": insights from numerical modelling[END_REF]Govers and Wortel, 1999 and many others]. These kinds of complex models have not fully replaced simple kinematic ones postulating a relationship between horizontal shortening or lengthening and the amplitude of vertical displacement of lithospheric boundaries [e.g. [START_REF] Mckenzie | Some remarks on the development of sedimentary basins[END_REF][START_REF] Weissel | Flexural uplift of rift flanks due to mechanical unloading of the lithosphere during extension[END_REF][START_REF] Kooi | Lithospheric necking and regional isostasy at extensional basins: part 1, Subsidence and gravity modeling with an application to the Gulf of Lions margin (SE France)[END_REF][START_REF] Van Der Beek | Mechanism of extensional basin formation and vertical motions at rift flanks: constraints from tectonic modeling and fission track thermochronology[END_REF][START_REF] Spadini | Thermomechanical modelling of the Tyrrhenian sea: lithospheric necking and kinematics of rifting[END_REF][START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF]. Indeed, although kinematic models often do not have a fully appropriate physical background, they are still used in quantitative modelling studies, for example in sedimentary basin analysis. These models are especially convenient because they permit formulation of effective algorithms for solving inverse problems, such as reconstruction of the history of sedimentary basin formation using data on infill thickness, age and lithology and on the structure of underlying crust. Using horizontally varying strain these models are capable of reproducing complex tectonic structures forming in laterally inhomogeneous lithosphere and compare well different geophysical and geological data.

One of the earliest ways in which the lithosphere (or crust) was idealised in order to facilitate modelling was a thin elastic plate. It was used for the first time by Vening Meinesz [1931] and later by [START_REF] Gunn | Quantitative aspects of juxtaposed ocean deeps, mountain chains and volcanic ridges[END_REF] as a proxy for modelling the response to loading. This simple model, reintroduced by [START_REF] Walcott | Flexural rigidity, thickness and viscosity of the lithosphere[END_REF], is still widely in use today (e.g. [START_REF] Jordan | Gravity anomalies, flexure and the elastic thickness structure of the India-Eurasia collisional system[END_REF][START_REF] Leever | Late orogenic vertical movements in the Carpathian Bend Zoneseismic constraints on the transition zone from orogen to foredeep[END_REF][START_REF] Pérez-Gussinyé | Effective elastic thickness of South America and its implications for intracontinental deformation[END_REF]. The key parameter characterising the A c c e p t e d M a n u s c r i p t 4 rheology (or strength) of the lithosphere and, hence, its isostatic response to external loads is then the effective elastic thickness T e (simply linked to the effective flexural rigidity of the plate). Of course, T e is an -effective‖ parameter that serves as a proxy for what undoubtedly is a more complex lithospheric rheological structure. Several studies were made in which T e was compared to how more complex, more realistic, numerical rheological models of lithosphere responded, in order to understand better the controls on T e . [START_REF] Burov | The effective elastic thickness (Te) of continental lithosphere: what does it really mean? (Constraints from rheology, topography, and gravity)[END_REF]1996], for example, found that non-elastic effects and possible decoupling between layers in the lithosphere could play an important role in how a layered rheology was expressed as effective T e value, especially in continental lithosphere.

Thin elastic plate theory has provided a very successful way of modelling regional isostatic response, in which only vertical forces (surface load and buoyancy) are balanced, although sometimes horizontal (in-plane) forces have also been incorporated into the mechanical equilibrium equation [e.g. [START_REF] Cloetingh | On a tectonic mechanism for regional sea level variations: Earth Planet[END_REF][START_REF] Stephenson | Isostatic response of the lithosphere with in-plane stress: application to central Australia[END_REF]. When dealing with active tectonic deformation in extensional or compressional processes, of course, the horizontal dynamics must be considered as well. Braun and Beaumont (1989b) found that the results of their numerical modelling, incorporating complex rheologies, could be interpreted as if extension occurred in two stages. The first stage included the deformation by external forces in the absence of isostatic rebound and the second one was the isostatic response of the lithosphere to that deformation. They hypothesized that within the lithosphere there exists a flat level that does not move vertically during the first stage. For this level they introduced the term -necking depth‖or z n . Necking depth was suggested to depend on depth variations of lithospheric strength.

Later, the necking depth was used as a mechanical parameter in numerous basin modelling studies (involving backstripping of sedimentary loads and forward modelling of thermal subsidence; cf. [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF], in which z n controlled the geometry of the rifted A c c e p t e d M a n u s c r i p t 5 lithosphere at the end of the active extensional period while T e controlled its response to subsequent sedimentary and thermal loads. [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF] considered z n , like T e , to be a proxy for the actual, more complex, rheology of the lithosphere and, accordingly, that insights into the rheology of the lithosphere could be gained by looking at possible systematic relationships between T e and z n (given other -known‖ parameters such as crustal and lithosphere thicknesses, heat flow, and tectonic age).

The results of this endeavour were somewhat ambiguous. On the basis of numerical modelling, [START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF] concluded that there is no one-to-one relationship between depth to the necking level and strength distribution within the lithosphere.

However, the use of the necking level appeared to be successful for detailed analyses of the structure and evolution of numerous sedimentary basins (for bibliography see [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF]) and in explaining the formation of rift shoulders [e.g. Braun and Beaumont, 1989b and many others]. Indeed, the necking-level model approximates well the numerical results obtained adopting complex lithospheric rheologies [START_REF] Mitrovica | Tilting of the continental interiors by the dynamical effects of subduction[END_REF]Braun and Beaumont, 1989a;[START_REF] Bassi | Factors controlling the style of continental rifting": insights from numerical modelling[END_REF]and others]. Accordingly, it is of interest to analyse further the meaning of -necking level‖ in the context of more realistic lithospheric models, and to investigate its possible relationship to other tectonic parameters such as crustal thickness, heat flow or tectonic age, and, indeed, T e . Because of the inherent ambiguities in the complex numerical approaches to this problem [cf. [START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF], we choose to apply an analytical approach, reducing the number of free parameters and considering only several archetypical yield strength envelopes [Ranalli and Murphy, 1987] characterizing different lithospheric structures.

In the following part of the paper we focus on a dynamic description of continental lithosphere deformation under extension or compression. We have adopted the two stage scheme of Braun and Beaumont [1989b], mentioned above (deformation and consequent A c c e p t e d M a n u s c r i p t 6 isostatic rebound) and investigate the position of the level at which the first stage vertical displacement is equal to zero. Even though this level exists for both compressional and extensional deformations, we retain Braun and Beaumont's term -necking level‖ although they initially introduced it to describe extensional necking of the lithosphere only. We do simply this because the term is already well established in the geodynamics literature.

We present the general analysis of the problem and give an analytical solution for small deformations of a thin inhomogeneous elastic plate by forces applied to its side boundaries or to its bottom. This solution can be presented in the form of deformation around a necking surface ) ( x z n which, in general, is not flat. In the last part of the paper we present and discuss some results of our numerical calculations and show that a model with a flat equivalent

necking depth const z n 
provides a very close approximation to the exact solution after the isostatic adjustment stage. We then discuss which data are required to estimate the depth to the necking level and effective elastic rigidity of the lithosphere simultaneously.

Kinematic necking-level model.

Let us now briefly recall the approach frequently used in modelling lithospheric deformation by intraplate forces [e.g. [START_REF] Kooi | Lithospheric necking and regional isostasy at extensional basins: part 1, Subsidence and gravity modeling with an application to the Gulf of Lions margin (SE France)[END_REF][START_REF] Van Der Beek | Mechanism of extensional basin formation and vertical motions at rift flanks: constraints from tectonic modeling and fission track thermochronology[END_REF]Spadini et al., 1995]. The lithosphere is assumed to be composed of one or several layers, such that: (a)

U -the horizontal component of velocity or displacement vector within each layer does not depend upon the vertical coordinate z, and (b) the lithosphere deforms about some necking

level const z n 
, which does not move vertically during intraplate deformation in the absence of isostatic rebound. In the particular case of one layered homogeneous lithosphere, the horizontal (U) and vertical (W) components of velocity or displacement vector are related by the equation:

A c c e p t e d M a n u s c r i p t 7 W x z t z z U x t x n ( , , ) ( ) ( , )      (1)
(for the sake of simplicity, hereafter we consider only 2D problems with horizontal Ox axis and Oz axis directed upward). The deformation about the necking level is a first stage of the previously mentioned two-stage scheme. The second stage being the restoration of the isostatic balance which could have been disturbed during the first stage. It is implicitly assumed that the position of z n depends upon the depth distribution of mechanical properties of the lithosphere, which in its turn is mostly determined by temperature profile and lithology [START_REF] Ranalli | Nonlinear flexure and equivalent elastic thickness of the lithosphere[END_REF]Murphy, 1987, Afonso and[START_REF] Afonso | Crustal and mantle strengths in continental lithosphere: is the jelly sandwich model obsolete?[END_REF]. [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF] estimated z n to vary from 4 to 35 km and made an attempt to relate its depth to the thermal age of the lithosphere, its thickness and strain rate at the stage of extension. Considering these results, [START_REF] Fernandez | The role of rheology in extensional basin formation modeling[END_REF] concluded that -the relationship between strength envelopes and kinematic level of necking is more complicated than previously thought‖. They classified the necking-level based model as a kinematic one with rheological constraints.

Moreover, kinematic considerations presented in Appendix A show that this relationship cannot be found by a purely kinematic approach and calls either for additional assumptions (for example, of local isostatic equilibrium giving a 2D analogue of the [1978] model or the presence of a flat necking level, etc.) or for consideration of forces and strain-stress relationships, i.e. for a dynamic approach.

McKenzie

3. Dynamic model.

Statement of the problem

The rheology of the medium (i.e. the equations coupling the stress and the strain and/or strain rate tensors) must be assigned when constructing a dynamic model. The rheology of the lithosphere depends upon several factors including rock compositions, temperature, pressure,

A c c e p t e d M a n u s c r i p t 8 stress distribution and magnitude, strain rate [e.g. [START_REF] Goetze | Stress and temperature in the bending lithosphere as constrained by experimental rock mechanics[END_REF]Ranalli and Murphy, 1987;Kohlsted et al., 1995;[START_REF] Afonso | Crustal and mantle strengths in continental lithosphere: is the jelly sandwich model obsolete?[END_REF]. In recent years, various yield strength envelopes accounting for these factors were proposed. They mostly vary by the assumption of what has been adopted as a constant: the strain rate [e.g. Ranalli and Murphy;1987;[START_REF] Carter | Flow properties of continental lithosphere[END_REF][START_REF] Cloetingh | Mechanical structure, in: A Continental Revealed: the European Geotraverse[END_REF], external forces [START_REF] Kusznir | The distribution of stress with depth in the lithosphere: thermorheological and geodynamical constraints[END_REF] or some combination [Ershov and Stepheson, 2006]. Fig. 1 shows five typical yield strength envelopes for different thermal regime and lithospheric structure accounting for heat generation in the upper crust and for a strain rate    10 14  s 1 . Following Turcotte and Shubert [2002] the steady-state continental lithosphere isotherm is assigned: Ranalli and Murphy [1987] and «Normal» and «Alpine» are close to those of [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF].

    k z q T z T m s / ) ( k h Q h z r rad r / )) / exp( 1 ( 2     ,
To construct these envelopes, we used the linear Coulomb frictional law [START_REF] Sibson | Frictional constraints on thrust, wrench and normal faults[END_REF] which is mostly constrained by experimental data obtained for conditions corresponding to the upper crust. For higher pressure and temperature it overestimates the strength of rocks [see experimental data compiled in [START_REF] Ranalli | Nonlinear flexure and equivalent elastic thickness of the lithosphere[END_REF]Murphy, 1987 and[START_REF] Kohlstedt | Strength of the lithosphere: Constraints imposed by laboratory experiments[END_REF]. To account for this effect, we limited the strength of the upper mantle to 1 GPa. We use these According to numerous geological and geophysical studies, discontinuities generally accommodate compression or extension of rigid layers, hence elastic or plastic deformations are comparatively small. However, when a rigid layer contains numerous discontinuities, an equivalent continuous medium can be used for a generalized description of the macroscale deformations, although the choice of the stress-strain relationships is still open. Below we use linear stress-strain relationships to investigate the problem of lithosphere extensioncompression. For a Newtonian compressible stratified medium the problem can be investigated using the equations obtained in [START_REF] Mikhailov | Dynamics of the Earth' outer shell evolution under extension and compression[END_REF][START_REF] Mikhailov | Modelling of extension and compression of the lithosphere by intraplate forces[END_REF].

Simple constitutive laws have been successfully applied to model many geodynamic

processes. An impressive collection of such results can be found, for example, in Turcotte and Schubert [2002]. In our case, a simple model permits investigation of the main characteristics, and their dependence on rheology, of extensional and compressional structures of the lithosphere. In particular, the topography of the necking surface in our model is close to that of [START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF], which was found using a dynamic numerical modelling techniques and more realistic rheologies. Thus, the complex rheological behaviours, transient effects, and so on, incorporated in many numerical models can be considered to produce second order effects. Our results simply demonstrate the validity of the necking level model, which has indeed been successfully used for modelling many sedimentary basins [see [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF] for an overview].

Effective elastic model.

Let us use a heterogeneous pure elastic model to consider relatively fast deformations of the lithosphere by intraplate forces. Such a model provides an analogue (effective) description [START_REF] Lomakin | Theory of elasticity of inhomogeneous media[END_REF][START_REF] Ranalli | Nonlinear flexure and equivalent elastic thickness of the lithosphere[END_REF][START_REF] Burov | The effective elastic thickness (Te) of continental lithosphere: what does it really mean? (Constraints from rheology, topography, and gravity)[END_REF]1996]. To investigate the mechanical response to the loading due to extension or compression one may consider the problem of isostatic equilibrium of a thin elastic plate floating on a non-viscous substratum [START_REF] Turcotte | Geodynamics[END_REF]. Since deformations by intraplate forces are supposed to be irreversible and non-elastic, they can be considered separately from the associated isostatic response. As before, these processes will be referred as stage (a) and stage (b). For both stages, the purely elastic model only provides an approximation of the real, more complicated, processes, its parameters are effective and cannot be determined in laboratory (for stage (b)

this question was studied in detail by Burov andDiament [1995, 1996) depends only on the relative distribution of the effective Young's modulus with depth and to characterize the depths of compliant and rigid layers we used the yield stress envelopes (Fig. 1). In problem (b), the behaviour of the lithosphere can be specified by the flexural rigidity estimated from a combined analysis of gravity and topography.

In our model inhomogeneous mechanical properties of the lithosphere have to be assigned at the onset of deformation. Such inhomogeneity can be intrinsic (e.g. vertical rheological stratification of the lithosphere) or formed after intraplate forces are applied, as a result of faulting or rock damage. Subsequently, the mechanical properties of the lithosphere change only as a result of material displacements; accordingly, all transient effects resulting, for example, from time varying temperatures are ignored.

Let us first consider problem (a). We adopt the following notations: E x z ( , ) -apparent Young's modulus,  -Poisson's ratio assumed to be equal to 0.5 in order to obtain formulas ) and the asthenosphere ( 0 a z ) were flat.

In Appendix B we solve the problem of a deformation of the lithosphere by intraplate forces considering only small deformations of a thin plate. This implies that we have assumed (1) that L -equal to the characteristic scale of variations of E x z ( , ) and U x z ( , ) along horizontal axis Ox -is much larger than H -the characteristic scale of variations of E x z ( , ) and the vertical component W x z ( , ) along vertical axis Oz (this yields a small parameter

   H L 1); and
(2) that -magnitudes of displacements are also small, i.e. characteristic scales for the horizontal and vertical components of displacements are of the order of

u L 0 2   , w H 0 2   .
We then proceed to use dimensionless values by introducing characteristic scales and expand all components of displacement and stress in power series of the squared small parameter (  2 ).

Assuming further that the distribution of apparent Young's modulus can be presented

as ) z ( F ) x ( E ) z , x ( E 1 1  where ) z ( F 1
is a non-dimensional function that can be associated with the yield strength profile, we arrive at an analytical solution, which, in the first expansion terms, is analogous to a 2D state of plane stress [e.g. Turcotte and Schubert, 2002, section 3.5). But to find the displacement fields and position of the necking surface, it is necessary to analyse the equations for the second order terms that strongly depend on vertical stratification of the model.

From the obtained equations, it follows that the first term of the expansion of the horizontal component of displacement U 1 does not depend on the vertical coordinate z, i.e. , which can be presented in a form close to equation

(1):

x U x z z x f x U z W n     1 1 1 1 )) ( ( ) (       , (2)
where z x n ( ) is the position of the necking surface.

In Appendix B we obtain equation for z x n ( ) by considering second order terms of the expansion (equation B9). To express this equation in a more useful form for our analysis, let us introduce the depth:

dz ) z , x ( E dz z ) z , x ( E ) x ( d 0 l 0 a 0 l 0 a z z z z 1      . ( 3
)
By analogy with the centre of mass this depth can be identified as the centre of rigidity of every vertical section. Then, introducing n-th central moment ( 2  n

) in respect to the centre of rigidity as:

dz ) z , x ( E dz )) x ( d z ( ) z , x ( E ) x ( c 0 l 0 a 0 l 0 a z z z z n 1 n       (4)
one transforms equation for the necking surface (B9) to the form:

0 ) d z ( x 2 c x ) d z ( c 1 n 2 1 x 2 1 x 2 3 2 1 x 1 n 2 2             . ( 5
)
where

x U 1 1 x    
is the main component of the strain in x-direction.

During the first stage the model has to be somehow fixed relative to the vertical axis Oz by assigning an asymptotic value to The final result (after isostatic rebound) does not depend on this asymptotic value. Indeed, any shift of the asymptote only produces an additional constant in the total perturbation of the A c c e p t e d M a n u s c r i p t 13 initial isostatic balance (see below); thus, during the isostatic adjustment stage the lithosphere will return to the same place independent of its flexural rigidity. We place the asymptotic value at the centre of rigidity 1 d whose depth is different for different yield strength envelopes used in this study.

From equation ( 5) it follows that in the general case when  x const 1  , the necking surface



and then to the same value again) one finds that ) ( x z n is different at the first and the second steps. All this is in good agreement with the results of the numerical modelling of [START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF]. Note that, in the particular case of a horizontally

homogeneous plate, ) ( ) , ( z E z x E  ; thus (equation (B6c) of the Appendix B)  x a 1   const .
In this specific case, when a  0 (extension), equation ( 5) has the solution:

  1 n d z const.
When a  0 (compression), in addition to

  1 n d z
const, a periodic solution corresponding to a loss in stability exists as well (if the compressional force is big enough). In this particular case the necking level does not exist.

Our solution allows calculating the necking surface depth (equation 5) as well as the position of the top and the base of the lithosphere when the distribution of strength (apparent Young's modulus) and intraplate force

x ) x ( U ) x ( e 3 4 F 1 pl   
, where

dz ) z , x ( E ) x ( e ) x ( z ) x ( z l a   ,
(equation B6c of Appendix B) are known.

For applications of the model to actual observations it is important that, instead of F pl and the distribution of Let us consider now the isostatic response stage (problem (b)). If at t  0 the lithosphere was in local isostatic equilibrium, then the disturbance of the equilibrium state due to a deformation by external forces (referred below as a load) at time t is

q x t g x z t dz z x t z z x t z x t x z dz z z a a a a l a l ( , ) [ ( , , ) ( ( , ) ) ( , ) ( , ) ( , , ) ]           0 0 0 0 , ( 6
)
where g is the gravity acceleration,  a is the asthenosphere density, 0

l z , 0 a z , ) ( x z l and ) ( x z a
are the topography of the top of the lithosphere and the asthenosphere before and after deformation. Using the definition of the free mantle (floating) level (A2) and considering that, under assumed conditions,  x1 is independent from the coordinate z, one obtains:

) 1 /( ) ) ( ( ) , ( 1 1 x x fm n a z x z g t x q       . ( 7
)
This means that the load due to an extensional or compressional deformation varies in direct proportion with the spacing between the necking and the floating (free mantle) levels and does not depend on the density distribution within the lithosphere.

To describe the isostatic response we use a model of a thin homogeneous elastic plate [e.g.

Turcotte and Schubert

, 2002]:

) x ( q g dx d D a 4 4 s      (8),
where s D is the flexural rigidity of the lithosphere, (x) -the magnitude of elastic flexure under the load q x ( ) calculated from equation (7). For dependence of the flexural rigidity on the rheology of the lithosphere see Burov andDiament [1995, 1996]. Now we present and discuss some results of our model. As mentioned above, in order to arrive at our analytical solution we assumed for stage (a) that the distribution of the apparent Young's modulus at t=0 can be approximated as According to equations ( 4)-( 5), the depth to the necking surface as well as the disturbance of the isostatic equilibrium (specified by value of load q) depend on the relative depth distribution of the mechanical properties of the lithosphere, but not on their absolute values. The yield stress envelopes (Fig. 1) characterize the distribution of rigid and compliant layers and can be used to assign thicknesses and relative strengths of rigid layers at the top of the crust and below the Moho (i.e. function ) z ( F 1

) z ( F ) x ( E ) z , x ( E 1 1  where ) z ( F 1 is a non-dimensional
) within the lithosphere. We computed the depth to the necking surface shown in Fig. 2 for the five types of lithosphere (Fig. 1) with the horizontal component of strain tensor specified as

) x b exp( a 2 0 0 1 x      x 0 , (9)
where parameter a 0 governs the amplitude of extension or compression; It is important to emphasise that the complex geometry of Note that at the borders of the extensional basin, all of the lithosphere moves upward and, accordingly, the necking surface ) x ( z n for the -Alpine‖ lithosphere is situated above the top of the model (curve A1 on Fig. 2). The same behaviour was found in [START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF] even on their figures when necking level came out of the limits of the crust, it was fixed at the nearest crustal boundary.

Thus we conclude that even if the topography of the necking surface significantly changes close to the edges of the tectonic structure (where  x approaches zero on Fig. 2 A,B), the equivalent constant necking level provides a close final solution everywhere. For all tested cases, the difference in topography of all boundaries calculated for ) x ( z n and eqv n z after isostatic rebound never exceeds 1% of displacement in the centre of the structure.

The depth to the necking level mainly depends on the relative thickness of the rigid layers in the upper crust and below the Moho.

According to formulae (3)-( 4), when a depth distribution of mechanical properties in the lithosphere has odd symmetry about the centre of rigidity It is easy to calculate the position of the centre of rigidity when the lithosphere consists of n-layers each having a constant strength

i F ) z ( F  , when i 1 i z z z    . It is equal to:    i str i i mid i T z F z d 1 ,
where

1 i i i z z z     is the thickness of the i-th layer, 2 / ) z z ( z 1 i i mid i   
is the depth to its central part and

   k k k str z F T
is total (integrated) strength of the lithosphere. In other words, the depth to 1 d and consequently to egv n z , is equal to the weighted average of all layer's mid i z .

In particular, when the lithosphere contains two equally rigid layers with a negligibly small strength for all other layers the necking level is situated just in between them. If one of these layers is thicker and stronger than the others, the necking level shifts closer to the stronger layer. This result coincides with the guess of Spadini et al [1995] and with the numerical calculations of [START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF].

The depth to the necking surface slightly depends on the amount of strain. If parameter 0 a in equation ( 9) varies within two orders of magnitude, the depth to the necking surface 3)-( 4)).

Parameter b 0 in equation ( 9) controls the horizontal gradients of strain 

x , which through equation ( B6c) is related to the horizontal gradients of the vertically averaged apparent Young's modulus It is important to note that, although the equivalent necking depth can be located anywhere within these depth intervals since these intervals are quite narrow (see Fig. 1), the necking level depth can be used to characterize the rheological properties of the lithosphere.Fig. 4 shows the shape of the lithosphere upper and lower boundaries, as well as crustal thickness and depth to the necking surface computed for the -Normal‖ lithosphere model (Fig. 1). As a result of extension and consequent isostatic rebound, uplifts (shoulders) are formed at the both sides of an extensional basin. The height of shoulders depends on horizontal gradient of mechanical properties at the periphery of tectonic structure. For the same maximum extension ratio, the lower gradient results in the wider structure and both the depth to the necking level and the amplitude of the elastic flexure depend on the width of this gradient zone.

The simultaneous estimation of the strain distribution

) ( x x 
, depth to the equivalent necking level z n eqv and effective flexural rigidity of the lithosphere D s is an ill-posed inverse problem.

Let us now investigate how to determine simultaneously the strain distribution, depth to the equivalent necking level and effective elastic rigidity. If it is assumed that the shape of the upper and lower crustal boundaries ) ( x z l and ) (x z Moho after deformation are known from geophysical data, then we get:

) 1 ( ) ( ) ( 0 0 x x eqv n l l l z z z x z       ; ) 1 ( ) ( ) ( 0 , 0 , x x eqv n Moho Moho Moho z z z x z       (10),
where functions with index -0‖ denote the topography before deformation and . Equations ( 7)-( 8) reveal that this problem has a unique solution. However, the problem is ill-posed. Solid lines in Fig. 5 In this case z n eqv =21 km, while for the first solution z n eqv =50 km. In both cases z x l ( ) and the Moho topography are almost identical and would not be distinguishable with seismic imagery.

This example suggests that there is a trade-off between z n eqv and D s and, as a result they cannot be simultaneously estimated. To resolve the problem, independent data must be used. Probably this explains why no clear dependence of the necking level depth on thermal age and thickness of the lithosphere was previously inferred [START_REF] Cloetingh | Role of pre-rift rheology in kinematics of extensional basin formation: constraints from thermomechanical models of Mediterranean and intracratonic basins[END_REF].

The depth to z n

eqv can be estimated from seismic, seismological or geothermal data.

Let us now investigate how to estimate z n eqv . For this, we assume that the lithosphere has two uniform rigid layers: in the upper crust and below the Moho with thicknesses H 1 and H 2 respectively. We note A, the strength ratio (expressed as effective Young's modulus) of the two rigid layers. Fig. 6 shows the computed depth to z n eqv versus the thickness of the rigid layer below the Moho for three different values of A (line 1-3 on Fig. 6). In addition, Fig. 4  seem unrealistic for the lithosphere of any tectonically active region (see Fig. 1).

Estimates of the relative thicknesses of the rigid layers H

1 and 2 H and of their relative strengths A depend on the assumed geotherm and the composition and rheological parameters of lithospheric rocks.

The load resulting from lithospheric extension or compression does not equal the topographic weight;

To estimate the flexural rigidity D s from equations ( 7)-( 8) one has to specify the load that causes the flexure. Where D s is estimated from gravity data it is generally assumed that this load is equal to the weight of topography after tectonic deformation but before isostatic adjustment ( ) ( x z l

) [see e.g. [START_REF] Forsyth | Subsurface loading and estimates of the flexural rigidity of continental lithosphere[END_REF]. According to equations ( 7), ( 10) the load can be expressed as:

q x ( )  )] ) ( ( [ 0 l l c z x z g K   , where )] ( ) ( [ 0 0 , 0 eqv n l Moho l z z z z K    / ) z z ( eqv n fm   ) z z ( 0 , Moho fm  , 
c is the average density of the crust, g is the acceleration of gravity, z fm is the free mantle level, z l ,0 and z Moho,0 are the topographies of the top of the crust and of the Moho before deformation (supposed to be flat for the sake of simplicity).

Thus, the load is equal to the weight of the topography only in the specific case z n eqv = z Moho,0 .

A c c e p t e d M a n u s c r i p t 23 Fig. 1 shows the depth interval in which the equivalent necking level is situated. For all types of strength profiles except -Shield‖ this interval is situated within the crust; hence, 1  K .

Accordingly, all estimates made under the assumption K=1 underestimate the flexural rigidity. This is consistent with [START_REF] Forsyth | Subsurface loading and estimates of the flexural rigidity of continental lithosphere[END_REF] who explained abnormally low values of D s estimated for orogenic belts using topography only (through isostatic admittance) in terms of neglected subsurface loads.

When the necking surface is situated within the crust, then, during extension or compression, the topographies of the upper surface and of the Moho and other density interfaces increase in amplitude simultaneously but in opposite directions, thus partly isostatically compensating each other. When necking depth coincides with free mantle depth (which can occur for the lithosphere having -Hot‖ or -Alpine‖ yield strength profiles), the weight of the topography is completely compensated by displacement of the subsurface density interfaces. In this case, the tectonic structure keeps its initial local isostatic equilibrium state independent of the flexural rigidity of the lithosphere. [START_REF] Fernandez | The role of rheology in extensional basin formation modeling[END_REF] named this depth as -neutral necking‖.

Fig. 7 demonstrates the topography of lithospheric interfaces formed as a result of compression of the lithosphere having the -Normal‖ yield strength profile (Fig. 1). The zoomed in area shows flexural subsidence at the periphery of the compressional belt. In this specific case the load ( 7) is 20% less than the weight of the topography before isostatic adjustment. This presumably changes the morphology of admittance and coherence functions routinely used to estimate the flexural rigidity of the lithosphere in such areas. It may also be relevant to the geodynamics of foreland basin formation [e.g. [START_REF] Mikhailov | Foredeep basins: the main features and model of formation[END_REF].

Our study suggests that for structures formed by extension or compression (for example post-rift basins and foredeeps) the ratio between subsurface and surface loads (equal A c c e p t e d M a n u s c r i p t 24 to K-1) can be found a priori using the estimate of the depth to the necking level we discussed above.

Conclusions.

We investigated the deformation of the continental lithosphere using a model incorporating small deformations of a thin inhomogeneous plate by boundary forces: in-plane forces applied to its side boundaries and/or mantle forces applied to its base. Our results

demonstrate that the necking level model is valid at least for a first order description of the geodynamics of regional tectonic structures such as sedimentary basins, continental rifts, orogenic belts and passive continental margins. According to this model:

1. The vertical component of displacement W can be expressed as

x n )) x ( z z ( W    
; where x  is the horizontal component of strain, which does not depend on the vertical coordinate z; and z x n ( ) is the necking surface, where there is no vertical displacement before isostatic response takes place.

2. Beneath the area of the most significant deformation, the necking surface is almost flat. Thus z x n ( ) can be replaced by a constant z n eqv . The depth to this constant equivalent necking level depends on the strength distribution within the lithosphere.

3. The depth to z n eqv depends mainly on the relative thickness and strength of the rigid layers in the uppermost crust and below the Moho.

4. The simultaneous estimation of distribution of strain, equivalent depth to the necking level and effective flexural rigidity is an ill-posed problem. It may be solved when eqv n z is estimated a priori from independent data on the structure, composition and thermal state of the lithosphere (e.g., yield strength diagrams), seismic or seismological data.

A c c e p t e d M a n u s c r i p t 25 5. The load on the lithosphere in extensional or compressional areas must include both topographic weight and subsurface loads. In the absence of -geological loads‖ due to lateral heterogeneity in the lithosphere the total load can be evaluated using estimates of the depth to the necking level.

A c c e p t e d M a n u s c r i p t

where α is the thermal expansion coefficient. We assume also that the following initial conditions are known at t=0:

spatial distribution of temperature - ) 0 , , ( z x T
, density which corresponds to this temperature -) 0 , , ( z x  , and initial topographies of the top of the sedimentary cover -) 0 , (x z s , the lithosphere -) 0 , (x z l and the asthenosphere -) 0 , (x z a . To reduce the length of the formulas we consider initial topographies to be constants:

z z x s s 0 0  ( , ), z l 0  z x l ( , ) 0 and z a0  ) 0 , x ( z a
. Let us finally assume that before being deformed, the lithosphere was in a state of local isostatic equilibrium: ) x z t dz z z x t Here z fm corresponds to the so-called floating or free mantle level. Equation (A2) can be used as a mathematical definition of z fm . Its physical meaning is as follows: it is a level to which the top of the asthenosphere would rise if rocks of every vertical column were compressed to the density of the asthenosphere a  .

  ( , , ) ( ( , )), ( , ) ( ,
The topography of any material boundary z x t p ( , ) can be determined from the equation:

    z t U x t z x W x z t p p p    ( , ) ( , , ) , (A3)
assuming that at each instant of time this boundary marks the same material points. For the top of the sedimentary cover ( z z s  ) the term ( , )

x t should be added to the right side of the equation (A3) to account for material brought by sedimentation or removed by erosion [e.g. [START_REF] Mikhailov | Mathematical model of the processes of evolution of structures formed as a result of vertical movements[END_REF].

The temperature distribution (T) is given by:

dT dt T Q x z t rad     ( ) ( , , )   , ( A4
)
where ( , )

x z is thermal diffusivity and Q rad is radioactive heat generation per mass rock unit.

Combining equations (A1) and (A4) with the equation of continuity yields:

0 ) , , ( ) , ( ) , , ( ) , , (            z t z x W x t x U t z x dt t z x d       . ( A5
)
The vertical component of the velocity vector -W x z t ( , , ) can be expressed as:

) t , x ( f ) t , x ( dz ) Q ) T ( ( ) T T ( 1 x U ) z z ( ) t , z , x ( W a s z ) x ( z rad 0 fm a                      (A6)
where f x t ( , ) is an unknown function that depends on a number of factors including the mode of deformation and physical properties of the constituent lithospheric rocks. Thus, in For many tectonic processes one can assume, as [START_REF] Mckenzie | Some remarks on the development of sedimentary basins[END_REF] did, that the duration of lithospheric deformation by external forces is considerably shorter than the time required to re-establish the thermal equilibrium in the lithosphere. If so, two stages of structure formation can be recognized. Rather short stage I, with a duration of several My, includes two substages introduced by [Braun and Beaumont, 1989b]: I.1 -deformation by external forces and I.2isostatic rebound (if the equilibrium had been disturbed at stage 1 a ). Thus, for stage I equation (A6) takes a following form:

) , ( ) ( ) , , ( 1 t x f x U z z t z x W fm       (A7.1)
Stage II, that lasts up to 100 My, besides restoring thermal equilibrium, comprises sedimentation and erosion and the isostatic movements associated with them: Another way to determine this function is to consider the actual mechanism of the lithosphere deformation.

) , ( ) , ( ) ) ( ( ) ( 1 ) , , ( 2 ) ( 0 t x f t x dz Q T T T t z x W a s z x z rad a                 . (A7.

Appendix B:

Effective elastic model.

For the elastic problem, we used the following notations: E x z ( , ) -apparent (effective)

Young's modulus,  -Poisson's ratio assumed to be equal to 1/2 to obtain formulas similar to

(1 or A6), U x z ( , ) and W x z ( , ) -the horizontal and vertical components of the displacement vector. As above, we consider that at t  0 the lithosphere was in local isostatic equilibrium To find a solution of any elastic problem, functions of displacements have to obey Hooke's Law, which under our assumptions can be written as follows:

x x z x E    ) , ( 3
2 ) (   ,   x z   , E x z xz xz ( , ) 3    ,
where:

   x U x  ,    z W z  ,      xz U z W x  
are the components of strain tensor,

 x ,  z ,  xz -components of stress tensor,      ( ) /
x z 2 -mean normal stress, Equilibrium equations must also be included:

    x xz x z   0 ,     xz z x z   0 . (B1)
To describe lithospheric deformation by in-plane forces we set the following boundary conditions:

(a) free -surface condition at the top at z z x l  ( ) and at the base at z z x a  ( ) of the lithosphere: (B3)

  x xz nx nz cos( ) cos( )   0 ,   xz z nx nz cos( ) cos( )   0. ( B2
The problem was solved suggesting that deformations are small and the plate is thin.

Then, using the assumptions listed in the section 3.2, introducing dimensionless values as follows:

x L   , z H   , U x z L u ( , ) ( , )     2 , W x z Hw ( , ) ( , )     2 , ) , ( ) , ( 0   e E z x E  , A c c e p t e d M a n u s c r i p t 31  ) z , x (  ) , ( s E 0 2   
and replacing components of displacement by the components of strain tensor, one can rearrange equations (B1)-(B2) as:

                       2 2 2 2 3 3 0 3 2 3 0 ( ) ( ) ( ) ( ) e u s e u w e u w s e u                (B1´)     u w   0 (B2´)
with the boundary conditions accounting for

     w x z 2 
at the top and the bottom of the lithosphere being:

               0 w ) w u ( 3 e ) u e 3 2 s ( 0 w ) s u e 3 2 ( ) w u ( 3 e 2 2 4 2                            (B3´)
Let us expand the components of displacement and the mean normal stress in a power series of the squared small parameter ( 2  ). For example for the non-dimensional mean normal stress s:

... ) , ( ) , ( ) , ( ) , ( 3 4 2 2 1               s s s s

Substitute the expansion in terms of the small parameter in (

A  1 )-( A  3 ) and set equal to each other the terms of the same power of  . As a result for the zero and second power of  one obtains: Boundary conditions at     l ( ) and     a ( ) are as follows:

0 ) u e ( 1       (B4a) 0 ) u 3 e ( ) u e 3 2 s ( 1 1 1             (B4b) 0 )
s e u 1 1 2 3 0     (B5a),   u 1 0  , (B5b)     u w 2 1 0   , (B5c) s e u 2 2 2 3 0     , (B5d)
The solution of the system (B4) with the boundary conditions (B5) takes the form:

u u 1 1  ( )  , (B6a) s e u 1 1 2 3    , ( B6b
) const u e      1 ) ( (B6c)     u w 2 1 0   (B6d) s e u 2 2 2 3    (B6e)
where

      d E e l a z z   ) ( ) ( ) , ( ) (
. Equations (B6 a-c) are valid for arbitrary function ) gives the equation:

w u f u n 1 1 1 1                ( ) ( ( )) , ( B7
u u f f 2 2 2 1 2 1 2 2 ( , ) ( )             .
Equation (B7) is presented in the form containing   n ( ) , which specifies the position of the necking surface. The functions f 1 and f 2 can be determined from the following equations for the zero and the first moments ( n  0 1 , ) which under the adopted boundary conditions can be written in the form: [Ranalli and Murphy, 1987]. The crust has a quartzite rheology, except in the model of -Normal‖ lithosphere where the lower part of the crust has a diorite rheology (shown in light gray). The underlying mantle has an olivine rheology. Thicknesses of the crust and the lithosphere are listed in Table 1. M a n u s c r i p t

  characterize the relative position and thickness of -rigid‖ and -compliant‖ layers within the different types of the lithosphere.

  used below to specify large scale irreversible deformations including sliding along faults. The effective elastic models are widely used to describe deformation in non-elastic media [e.g.

  ]). Moreover, parameters in problems (a) and (b) describe different physical mechanisms of deformation and as a result have different values. Fortunately, as it is shown below, the solution of problem (a)

  are the horizontal and vertical components of the displacement vector. For the sake of simplicity, we assume that at t  0 the lithosphere was in local isostatic equilibrium (equation A3 in Appendix A) and the initial topography of the top of the lithosphere (



  ( ) . As a result, the vertical component of displacement is related to the horizontal strain through an unknown function ) ( 1 x f

  difficult to assign a priori), the horizontal component of the strain tensor ) ( 1 x x  can be assigned. This function is related to the final topographies of the lithospheric interfaces and can be found by modelling real data. Moreover, when assigning A c c e p t e to the deformation by in-plane forces, our solution also describes the deformation by force applied to the base of the lithosphere (by mantle drag). See Appendix B for more details.

  the tectonic structures as well as (through equation B6c) the horizontal gradients of mechanical properties of the lithosphere at the periphery of the tectonic structure.The coordinate origin is at the centre of the model. Fig. 2A and 2B show the position of the necking surface for narrow (b 0 =3.2 km -2 ) and wide (b 0 =0.8 km -2 ) areas of extension, respectively, when the asymptotic value for z x n ( ) is assigned as the centre of rigidity d1. (The A c c e p t e d M a n u s c r i p t 16 choice of asymptotic value does not affect the final result after isostatic rebound, as explained above). All curves z x n ( ) are nearly flat beneath the area of main deformations. Their forms become more complicated only at the periphery of the tectonic structure where 1 x is close to its asymptotic value x 0 (the same behaviour was obtained by[START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF]. As a consequence, the solution of the elastic problem with the necking surface ) can be closely approximated by a solution with an equivalent flat necking level value affects both the vertical component of displacement W 1 (equation 2) and the load (7). Therefore, the depth of this equivalent level is close to but not equal to the depth of the necking surface at the centre of the considered structure.

  the periphery of the area of main deformation almost completely vanish after isostatic adjustment. Therefore, the model with a flat necking level approximates considerably better the final geometry after isostatic rebound than after the extension of the first stage. This can be understood in terms of the main deformation occurring in a weaker zone (pre-existing or formed during deformation), while the surrounding non-deformed areas are stronger. Deformation in the transition zone, where the mechanical properties of the lithosphere rapidly change, is actually more complicated than predicted by the model with a flat necking level. The solid curve on Fig.3Ashows the topography of the top of the lithosphere having an -Alpine‖ yield strength profile and extended around the exact necking level calculated from equation (5). The graph marked by filled circles corresponds to a flat necking level at 19.8 km depth, which provides the best fit for the final result after isostatic rebound. Basin shoulders, during extension, relative to ) x ( z n move upward and, therefore, it A c c e p t e d M a n u s c r i p t 17 is more difficult to match the results of the first stage (before isostatic rebound) by the model with a flat necking level. The more complicated necking ) x ( z n manifested in the geometry of all other interfaces produces a more complicated perturbation of isostatic equilibrium in comparison to the model with flat necking level (Figure 3 B). Finally, after isostatic rebound, the topographies of the top of the lithosphere as well as all other boundaries (not shown) become almost identical (Figure 3 C).

  this particular case the necking surface is flat and situated at the depth 1 n d z  . Hence, the topography of the necking surface depends on asymmetry of the function ) rigidity 1 d , i.e. it mainly depends on the relative thickness and -strength‖ of rigid layers in the upper crust and below the Moho and only slightly depends on the mechanical properties of the middle part of the crust. As a result, the necking levels for a one-layered and a two-layered crustal model are very close, especially within the region of main deformation.

Fig. 2 A

 2 Fig.2 A, B presents the shape of the necking surface for the five different yield strength

in

  of main deformation moves less than 10%. Probably, the strong dependence of the depth to the necking surface with the extensional ratio)[START_REF] Govers | Some remarks on the relation of vertical motions of the lithosphere during extension and the necking depth parameter inferred from kinematic modelling studies[END_REF] can be explained by the difficulty to numerically integrate a finite element solution for a very small extensional ratio. In their calculations the necking surface was more stable and changed in the same range of 10%.4.3The depth to the equivalent necking level egv n z depends on the horizontal gradients of mechanical properties of the lithosphere and falls in the range from

  Fig.2Bshows the necking surface for the same previously studied five models, but for

  show the topographies of the top of the crust and of the Moho in the area of extension (the horizontal component of the strain x  is shown at the top of the figure) for the -Shield‖ model (Fig. 1) when D s is equal to 10 22 N m  . A very close estimate with slightly different strain (shown by crosses) was obtained for D s   5 10 22 N m  .

  6

  curves: curve 4 corresponds to the necking depth computed as the distance between the middle of the upper and lower rigid layers

  (A4) and (A6) can be used to determine the vertical component of the velocity vector and the topography of the boundaries if the following functions are available:(a) horizontal component of velocity vector U x t ( , ) , (b) rate of sedimentation and/or erosion ( , ) x t , (c) thermal diffusivity ( , )x z , (d) initial conditions, including temperature

  of the material reloading by sedimentation / erosion and into movements within the lithosphere and asthenosphere related to their intrinsic stress relaxation [e.g.[START_REF] Mikhailov | Dynamics of the Earth' outer shell evolution under extension and compression[END_REF]. Consideration of 2D or 3D dynamics at the second stage is beyond the scope of this paper. Below we consider only the first stage of the processthe initial deformation by external forces (stage I.1) and the consequent isostatic rebound (stage I.2).

  A7.1), let us suppose that at stage I the deformation rate was small and resulting tectonic structures wide enough to neglect the contribution of elasticity in the isostatic balance. Hence, the condition of the local isostatic equilibrium (equation 3 in the main text) can be used as the supplementary equation at any t  0 .

  and the initial topographies of all density interfaces including the top of the lithosphere (

  ) taking into account the following relationship, which follows from the equation of continuity (B2´):

Fig. 1 .

 1 Figure captions

Fig. 2 Fig. 3 .

 23 Fig.2. a. Topography of the necking surface for a relatively narrow area of extension (the

Fig. 4 .

 4 Fig.4. Topography of the top ( Z lith ) and the bottom ( a Z ) of the lithosphere and the crustal

Fig. 5 .

 5 Fig.5. Close solutions for two different values of effective elastic rigidity of the lithosphere.

Fig. 6 .Fig. 7 .

 67 Fig.6. Depth to the equivalent necking level versus thickness of the rigid layer below the

Figure

  Figure

Figure

  Figure

Figure

  Figure

  function, which can be associated with the yield strength profile. In this case the depth to the centre of rigidity 1 d and all central moments n c are constant.

	4.1. The necking surface z x n ( ) is almost flat in the area of main deformation; thus	z x n ( ) can
	be replaced by an equivalent flat necking level z n eqv .

  For sedimentary basins an alternative approach is to use subsidence curves to estimate the extensional ratio  

	not depend on the depth to z n eqv , thus the function		x	1 (		x	)
							  1	[McKenzie, 1978].
	) x a n u s c r i p t ( x x  is the horizontal component of the strain tensor. The crustal thickness ) ( ) ( ) ( z x z x H l Moho   does A c c e p t e d M 21 ( ) ( ) x x x Thus, since )  is known (probably with some error), the problem reduces to a ( x x simultaneous estimate of z n eqv and the flexural rigidity D s from the topography of z x ) (x l ( ) and z Moho



can be estimated from crustal thickness before and after deformation.

Table 1

 1 Parameters of the lithosphere used to calculate strength diagrams (Fig.1).

	Table
		Thickness of	Thickness of the	Heat flow at the base of
		the crust, (km)	lithosphere, (km)	the lithosphere,	mW	m	2
	"Normal"	40	100	41.9
	"Alpine"	40	75	55.9
	"Hot"	30	50	85.4
	"Shield"	40	150	27.8
	"Collisional"	60	150	27.8

Acknowledgments.

We are grateful to two anonymous reviewers for their thoughtful and provocative comments, which helped us improve the presentation of our results. VM was supported by grants 09-05-00258 and 09-05-91056 of Russian foundation for Basic Research. This paper is IPGP contribution XXX.

Appendix A:

More general considerations.

Let us consider a general model of the formation of regional tectonic structures under extension or compression. We assume that the lithosphere deforms by in-plane (far field) forces and/or forces applied to its base resulting from mantle dynamics. We introduce the Cartesian coordinates xOz with the axis Oz directed upwards and denote the horizontal and vertical components of the velocity field initiated within the lithosphere by tectonic forces as U and W respectively. We suppose that U slightly depends upon the vertical coordinate z, so that in a first approximation U U x t lith  ( , ) (the assumption is valid for a thin plate, see Appendix B). We then assume that within the lithosphere and sedimentary cover, the density (  ) depends upon x and z coordinates and temperature (T), neglecting density/pressure dependency:


L Ni

F Demami
email: fdemami@univ-rennes1.fr

R Rogel

A C Salaün

L Pichon

Fabrication and electrical characterization of silicon nanowires based resistors

Silicon nanowires (SiNWs) are synthesized via the Vapor-Liquid-Solid (VLS) mechanism using gold (Au) as metal catalyst and silane (SiH 4 ) as precursor gas. Au nanoparticles are employed as liquid droplets catalysis during the SiNWs growth performed in a hot wall LPCVD reactor at 480°C and 40 Pa. SiNWs local synthesis at micron scale is demonstrated using classical optical photolithography process. SiNWs grow with high density anchored at the dedicated catalyst islands. This resulting network is used to interconnect two heavily doped polysilicon interdigitated electrodes leading to the formation of electrical resistors in a coplanar structure. Current-voltage (I-V) characteristics highlight a symmetric shape. The temperature dependence of the electrical resistance is activated, with activation energy of 0.47 eV at temperatures greater than 300K.

Introduction

Owing to their physical and electrical properties, silicon nanowires (SiNWs) are currently attracting much attention as promising components for future nanoelectronic devices such as nanowire field effect transistors [1], photonic and optoelectronic devices [2]. SiNWs can be prepared by the top-down approach, using various advanced methods such as e-beam [3], AFM [4] or deep UV [5] lithography. The main disadvantage of these advanced lithographic tools with nanometer size resolution rests on the high cost generated. The bottom-up approach usually employs metal catalytic growth for preparation of SiNWs. Thus, VLS (Vapour-Liquid-Solid) growth technique uses metallic nanoparticles. The metal catalyst used is usually gold although other metals are also employed [6,7]. However, this approach suffers from the difficulty in precisely positioning the device.

In this work, we present a technological process based on localized metal deposition to activate vapour deposition synthesis of SiNWs by the VLS technique. The silicon nanowire network resulting is used to interconnect electrodes to achieve an electrical resistor in a coplanar structure.

Experiment details

The fabrication process of the SiNWs based resistors includes 2 photolithographic steps and process stages are represented on Figure 1. First, substrate (silicon wafer or glass substrate) is capped with a SiO 2 layer deposited by Atmospheric Pressure Chemical Vapor Deposition (APCVD) and acts as a buffer film. Then, a heavily phosphorous in-situ doped amorphous silicon layer is deposited by Low Pressure Chemical Vapor Deposition (LPCVD) at 550°C and 90 Pa. Subsequent solid phase crystallization is performed at 600°C under vacuum to get a highly doped polysilicon film. This film is then patterned and dry etched to define the geometry of the comb shape electrodes (interdigitated structure shown in Figure 1-a). Au thin film (≤ 5 nm) is then deposited by thermal evaporation and locally removed using a lift off technique in order to define precise location for silicon nanowires synthesis (Figure 1-b). Finally, SiNWs are grown in a hot wall LPCVD reactor using silane as precursor gas in order to link together the two electrodes. Deposition temperatures and pressure are respectively 480°C and 40 Pa. This synthesis results in a tangled growth of SiNWs insuring electrical contact between the two doped electrodes (Figure 1-c). LPCVD parameters are chosen in order to avoid silicon deposition on areas uncovered with Au thin film. Different interdigitated structures can be achieved, varying number of teeth and local deposition of gold (on each tooth or half teeth). This fabrication process can be achieved either on silicon wafer or on other substrate compatible with the maximum process temperature (600°C). I-V characteristics were performed by using a HP 4155 B semiconductor parameter analyzer. Dark current versus temperature was performed by using a Keithley 617 electrometer. Samples were placed in a cryostat under vacuum and temperature varied from 100 K to 580 K. SEM observations were also performed to show SiNWs shape and length. Figure 3 shows the I-V characteristics of the SiNWs comb shape resistors after localized SiNWs synthesis. As a reference, I-V characteristic before SiNWs growth was also plotted (curve ×) in order to insure initial electrical insulation between the two electrodes. This set of measurements highlight a good qualitative correlation between current levels measured and parameters as space between teeth (w) and the number of teeth of the comb resistor (n). Indeed, current increases as space between teeth (w) decreases and number of teeth increases. On the other hand it is difficult to make a fair performance comparison between comb shapes because the number of participating nanowires can be vastly different.

Results and discussion

-2,0 -1,5 -1,0 -0,5 0,0 0,5 1,0 1,5

Voltage (V) w=5µm,n=9 w=5µm,n=5 w=10µm,n=5 Before SiNWs synthesis Figure 4 shows the temperature dependence of the current under obscurity and the electrical resistance of the interdigitated structures. At temperatures above 300K, a well defined activation energy of 0.47 eV can be determined by fit. At temperatures below 200K, the resistance tends to saturate indicating a different conduction mechanism coming from carriers hopping between localized states. Such states can be associated to defects located either in the bulk or at the surface of the wires. Main interest of SiNWs rests on their high surface that can be sensitive to charges. Figure 5 shows results of a simple test that consists in exposing SiNWs to smoke in a cryostat. Upon exposure to smoke, current is found to dramatically increase and this phenomenon exhibits a reversibility aspect when air is re-introduced. This first demonstration serves only as a proof-of-concept for a new kind of sensor based on silicon nanowires. Indeed, SiNWs are of special interest because they are expected to play an important role as their surface can be modified with functional group.

Conclusion

A random network of silicon nanowires was synthesized via the VLS mechanism in a LPCVD reactor using Au as metal catalyst and silane as precursor gas in order to achieve SiNWs comb shape resistors between two electrodes. Temperature dependence of the electrical resistance of the structure is of activated type, with a well defined activation energy of 0.47 eV in the high temperature range. At temperatures below 100K, a hopping conduction mechanism is observed. SiNWs are expected to be used as highly sensitive sensors for a wide range of functional applications.

Figure 1 :

 1 Figure 1: Key device fabrication steps of the SiNWs based resistors. a) definition of the comb shape electrodes geometry on a SiO 2 capped substrate. b) Thin film Au local deposition after lift off technique. c) Growth of SiNWs by LPCVD.

Figure 2 -Figure 2

 22 Figure 2: a) SEM observation of SiNWs synthesized. Length can exceed 10 µm. b) SEM observation of SiNWs based resistor (Au was deposited on local areas on half teeth).

Figure 3 :

 3 Figure 3: I-V characteristics of SiNWs comb shape resistors depending on the space between teeth (w) and the number of teeth (n).

1 )Figure 4 :

 14 Figure 4: Resistance-Temperature dependence of SiNWs based devices.

Figure 5 :

 5 Figure 5: Current response of the device as a function of time upon exposure to smoke.