Search

Monitoring the evolution of freshwater water resources worldwide is challenging yet essential to prevent severe water scarcity. In that regard, the Gravity Recovery and Climate Experiment and Follow-On satellite missions (GRACE/-FO) have successfully proven to complement in-situ hydrological observations by giving access to variations in continental water storage through measurements of the Earth’s gravity field since 2002.However, the measurement noise and errors inherent to the GRACE/-FO missions design and data processing cause estimation bias and spatial leakage that limit applications to large-scale hydrological basins and complicate geophysical interpretation. Moreover, missing observations during and between the missions, and in particular the 11-month gap between missions make if difficult to monitor long-term mass variations. To overcome these issues, we build a procedure, based on a spectral analysis by Multichannel – Singular Spectrum Analysis (M-SSA), that uses spatio-temporal correlations of the GRACE/-FO time series to fill data gaps and reduce a significant portion of the distinctive noise pattern while maintaining the best possible spatial resolution. This processing reveals hydrological signals that are less well or not resolved by other processing strategies.We seek validation of the GRACE/-FO M-SSA solution by comparison with independent estimates of regional hydrological mass balance. We first develop a new forward modeling approach to derive regional hydrological balance that accounts for potential bias caused by necessary filtering of the GRACE/-FO data. In particular, we account for signal leakage of a known source while estimating regional mass balance.The Caspian Sea, due to its relative isolation, large size and significant mass variations over the past decades, offers an interesting example to validate mass balance derived from GRACE/-FO measurements. We compare space-gravity derived hydrological mass balance with independent altimetry data, accounting for thermosteric effects using in-situ data and ground deformation caused by mass redistribution, and show consistent results on the long term variations. However, we also discuss potential causes of a significant discrepancy in seasonal amplitude.

Ever since the Maule earthquake (Mw8.8, 2010), a quick vertical uplift is measured thanks to GNSS in the Andes, facing the rupture zone (~250 km to the trench). Models built for the Maule earthquake [Klein et al. 2016] have highlighted that a low-viscosity channel is required to explain the post-seismic uplift. This channel is located along the slab between 50 km and 130 km depth and has a viscosity of a few 1017 Pa.s - lower than in the asthenosphere, 1018 Pa.s. After the Illapel earthquake (Mw8.3, 2015), simple observations on GNSS time-series show that no uplift occurred in the Andes at an equivalent distance to the trench than in the case of the Maule earthquake. The subduction in the Illapel region is characterized by a flat-slab (called the Pampean flat-slab) in contrast with the normal-dipping subduction in the region of Maule.Here, we investigate what is the impact of the Pampean flat-slab on the post-seismic deformations of the Illapel earthquake. In particular, we try to understand whether the presence of the flat-slab inhibits the effect of the low-viscosity channel. For that purpose we compare GNSS vertical displacements with predictions in both regions of Maule and Illapel from a 3D spherical finite-element model that accounts for the slab geometry of the Chilean subduction zone.

Thanks to space geodesy we know with a millimetric precision how the lithosphere deforms at each stage of the seismic cycle. In particular, during the post-seismic phase, it can deform over thousands of kilometers and for decades. These deformations are partly due to viscoelastic relaxation of the asthenosphere.In a previous work, we have shown that at the temporal and spatial scale of the seismic cycle, the viscoelastic relaxation can be modeled by a linear creep law [Boulze et al. 2022]. Therefore, because of the linearity of the creep law, the superposition principle applies and the present day deformation is simply the sum of the post-seismic deformations induced by past earthquakes. Based on this result, the objective of our work is to determine what slip history is needed on the Chilean subduction interface to reproduce the current deformation of South America, which is well measured by GNSS.To investigate this challenging problem, we first develop a 3D spherical finite-element model of the Chilean subduction zone. This model covers the entire South American continent and incorporates a slab with a geometry described by Slab2.0 model [Hayes et al. 2018]. Then, we compare different ways to model the seismic cycle using the backslip theory [Savage 1983]. Finally, by comparing GPS time-series with our seismic cycle model prediction, we discuss many ingredients of the model: e.g. the viscosity of the asthenosphere (Maxwell, Burgers), the impact of a flat slab and low viscosity zones, the magnitude and extent of historical earthquakes.

In continental interiors, tectonically-driven deformation rates are low, often to the point where they are undetectable with modern geodesy. However, a range of non-tectonic surface processes, particularly relating to hydrological, cryospheric, and sedimentological mass changes, can produce strain-rates which on geologically-short timescales are substantially greater than those produced by tectonics. Here, we illustrate the problem that such transient strain rates may pose in low-strain environments by considering the impact that the growth and decay of the Fennoscandian and Laurentian ice sheets over the Holocene had on Europe and North America respectively. Induced deformation extended far beyond the periphery of the ice sheets, with the potential to impact on seismicity rates thousands of kilometres south of the maximum ice extent. We consider how the modelled non-tectonic deformation would have interacted with several known active fault systems, including the European Cenozoic Rift System and the New Madrid fault system. In low strain continental interiors, seismic hazard assessment -- crucial for the long-term planning of critical infrastructure, including nuclear waste disposal -- is often dependent on sparse information from observational and historical seismicity, and from paleoseismological studies of surface fault systems. We recommend that for a more complete seismic hazard assessment, the impact of non-tectonic transients should be considered -- both in the context of the role such transients may have played in recent seismicity, and the role they may play in seismicity to come. Whilst such consideration has previously been given to the direct impact on glacial loading in areas directly glaciated, we show that it should also be considered much more broadly.

Spatial geodesy through GNSS allows to measure with a millimetric precision the displacement of the lithosphere during the seismic cycle. The post-seismic part of this cycle can last for decades, traduced by long-lasting surface deformations. One of the main physical processes involved in the postseismic deformation is the viscoelastic relaxation in the asthenosphere. However, a long term debate persists about the involved rheology of the asthenosphere: Is the viscosity highly variable from one region to the next and is effective viscosity Newtonian (linear) or non-Newtonian (non linear)? To investigate these questions, we compare the horizontal post-seismic deformations induced by three Chilean megathrust earthquakes: Maule Mw8.8 (2010), Illapel Mw8.3 (2015) and Iquique Mw8.1 (2014). For each earthquake, we select permanent GPS stations along profiles perpendicular to the trench, extending as far as 1400 km. We calculate the ratio of the cumulative post-seismic (post) over 5 years and the coseismic (co) displacements for each station. Remarkably, at a given distance to the trench, the post/co ratios from the three earthquakes differ only slightly.What can be the interpretation of this observation in terms of rheology of the asthenosphere? First we can analyse the response of the asthenosphere in the case of homothetic earthquakes of different magnitude. The post/co ratio obeys simple analytical relationships: For a Newtonian rheology, it is simply a function of the (time/viscosity) ratio. For a non-Newtonian viscosity with a stress exponent n=3, the timescale becomes inversely proportional to M**2, where M is the moment of the earthquake. We show that these relationships are only slightly modified when the earthquakes are no longer homothetic and that the post/co ratio is a good proxy to quantify the strain-rate and stress ratio in the underlying asthenosphere. As a conclusion, the post-seismic deformation following the three Chilean earthquakes reveals very similar viscosity. In particular, a Newtonian, rather than a non-Newtonian, effective viscosity is required to explain the post-seismic deformation process.

International audience; The continuous redistribution of water involved in the hydrologic cycle leads to deformation of the solid Earth. On a global scale, this deformation is well explained by the loading imposed by hydrological mass variations and can be quantified to first order with space-based gravimetric and geodetic measurements. At the regional scale, however, aquifer systems also undergo poroelastic deformation in response to groundwater fluctuations. Disentangling these related but distinct 3D deformation fields from geodetic time series is essential to accurately invert for changes in continental water mass, to understand the mechanical response of aquifers to internal pressure changes as well as to correct time series for these known effects. Here, we demonstrate a methodology to accomplish this task by considering the example of the well-instrumented Ozark Plateaus Aquifer System (OPAS) in the central United States. We begin by characterizing the most important sources of groundwater level variations in the spatially heterogeneous piezometer dataset using an Independent Component Analysis. Then, to estimate the associated poroelastic displacements, we project geodetic time series corrected for hydrological loading effects onto the dominant groundwater temporal functions. We interpret the extracted displacements in light of analytical solutions and a 2D model relating groundwater level variations to surface displacements. In particular, the relatively low estimates of elastic moduli inferred from the poroelastic displacements and groundwater fluctuations may be indicative of aquifer layers with a high fracture density. Our findings suggest that OPAS undergoes significant poroelastic deformation, including highly heterogeneous horizontal poroelastic displacements.

The evolution of the Greenland Ice Sheet (GIS) is an important indicator of climate change and driver of sea level rise. However, providing accurate GIS ice mass balance remains a challenge today. Here, we propose to combine a unique set of geodetic measurements to improve our knowledge of the GIS spatial and temporal evolution. We attempt at reconciling satellite observations of ice volume with regional GNSS velocities estimates and time variable space gravity measurements over the 2003-2009 and 2011-2015 periods. The GIS mass variations are inferred from satellite altimetry for large ice sheets (IceSat and CryoSat-2; Sorensen et al.,2018, Simonsen et al.,2017) and digital elevation models (DEMs) generated from multiple satellite archives for peripheral glaciers (Hugonnet et al.,2020), associated with IMAU-FDM firn model (Ligtenberg et al., 2011). The spatial and temporal variations of the gravity field are given by the GRACE mission for which we use a solution where smaller wavelength signals are preserved (Prevost et al., 2019).To resolve short wavelengths load variations affecting the displacement of nearby GNSS stations, we use Green’s functions for vertical crustal displacements assuming purely elastic Earth properties (Martens et al., 2019). We first assume that the deformation is entirely due to recent ice melting and show that vertical elastic displacements predicted by our refined ice loading model, while in good agreement with observations in some regions, cannot explain observations overall. In particular, observations and model disagree in the Southeastern and the Northern parts of Greenland.We then explore potential viscoelastic deformation associated with short-term rheology of the asthenosphere induced by recent ice melting that could explain the observed GNSS displacements. We define a history of ice loading from 1900 to 2009 using both in situ and satellite altimetric measurements, compute today’s associated viscoelastic deformation for various mantle rheologies and discuss the potential contribution of ice melting since the little ice age to current observations. Remaining differences between observations and viscoelastic models may reflect a viscoelastic deformation induced by glacial isostatic adjustment. We discuss implications in terms of regional rheological constraints, and impact on estimates of present-day GIS ice mass budget.Hugonnet, R. (2020). A globally complete, spatially, and temporally resolved estimate of glacier mass change: 2000 to 2019. In EGU 2020. Ligtenberg, S. R. M., et al (2011). An improved semi-empirical model for the densification of Antarctic firn. The Cryosphere, 5, 809-819.Martens, H. R.,et al (2019). LoadDef: A Python‐based toolkit to model elastic deformation caused by surface mass loading on spherically symmetric bodies. Earth and Space Science, 6(2), 311-323.Prevost, P., et al (2019). Data-adaptive spatio-temporal filtering of GRACE data. Geophysical Journal International, 219(3), 2034-2055.Simonsen, S. B., & Sørensen, L. S. (2017). Implications of changing scattering properties on Greenland ice sheet volume change from Cryosat-2 altimetry. Remote Sensing of Environment, 190, 207-216.Sørensen, L. S., et al (2018). 25 years of elevation changes of the Greenland Ice Sheet from ERS, Envisat, and CryoSat-2 radar altimetry. Earth and Planetary Science Letters, 495, 234-241.

SUMMARYMeasurements of the spatio-temporal variations of Earth’s gravity field from the Gravity Recovery and Climate Experiment (GRACE) mission have led to new insights into large spatial mass redistribution at secular, seasonal and subseasonal timescales. GRACE solutions from various processing centres, while adopting different processing strategies, result in rather coherent estimates. However, these solutions also exhibit random as well as systematic errors, with specific spatial patterns in the latter.In order to dampen the noise and enhance the geophysical signals in the GRACE data, we propose an approach based on a data-driven spatio-temporal filter, namely the Multichannel Singular Spectrum Analysis (M-SSA). M-SSA is a data-adaptive, multivariate, and non-parametric method that simultaneously exploits the spatial and temporal correlations of geophysical fields to extract common modes of variability.We perform an M-SSA analysis on 13 yr of GRACE spherical harmonics solutions from five different processing centres in a simultaneous setup. We show that the method allows us to extract common modes of variability between solutions, while removing solution-specific spatio-temporal errors that arise from the processing strategies. In particular, the method efficiently filters out the spurious north–south stripes, which are caused in all likelihood by aliasing, due to the imperfect geophysical correction models and low-frequency noise in measurements.Comparison of the M-SSA GRACE solution with mass concentration (mascons) solutions shows that, while the former remains noisier, it does retrieve geophysical signals masked by the mascons regularization procedure.

The Earth's surface was depressed under the weight of ice during the last glaciations. Glacial Isostatic Adjustment (GIA) induces the slow recession of the trough that is left after deglaciation and is responsible for a contemporary uplift rate of more than 1 cm/yr around Hudson Bay. The present-day residual depression, an indicator of still-ongoing GIA, is difficult to identify in the observed topography, which is predominantly sensitive to crustal heterogeneities. According to the most widespread GIA models, which feature a viscosity of 2 – 3 × 10 21 Pa s on top of the lower mantle, the trough is approximately 100 m deep and cannot explain the observed gravity anomalies across North America. These large anomalies are therefore usually attributed to subcontinental density heterogeneities in the tectosphere or to slab downwelling in the deep mantle. Here, we use observed gravity gradients (GG) to show that the uncompensated GIA trough is four times larger than expected and that it is the main source of the North American static gravity signal. We search for the contribution to these GGs from mantle mass anomalies, which are deduced from seismic tomography and are mechanically coupled to the global mantle flow. This contribution is found to be small over Laurentia, and at least 82% of the GGs are caused by GIA. Such a contribution from GIA in these GG observations implies a viscosity that is greater than 10 22 Pa s in the lower mantle. Our conclusions are a plea for GIA models with a highly viscous lower mantle, which confirm inferences from mantle dynamic models. Any change in GIA modelling has important paleoclimatological and environmental implications, encouraging scientists to re-evaluate the past ice history at a global scale. These implications, in turn, affect the contribution of bedrock uplift to the contemporaneous mass balance over Antarctica and Greenland and thus the present-day ice-melting rate as deduced from the GRACE space mission. Additionally, studies of the thermo-chemical structure of the lithosphere/crust under North America that exploit gravity or geodetic data should be corrected for a GIA model, which is not the case today.

The occurrence of large earthquakes in stable continental interiors challenges the applicability of the classical steady-state ‘seismic cycle’ model to such regions. Here, we shed new light onto this issue using as a case study the cluster of large reverse faulting earthquakes that occurred in Fennoscandia at 11-9 ka, triggered by the removal of the ice load during the final phase of regional deglaciation. We show that these reverse-faulting earthquakes occurred at a time when the horizontal strain-rate field was extensional, which implies that these events did not release horizontal strain that was building up at the time, but compressional strain that had been accummulated and stored elastically in the lithosphere over timescales similar to or longer than a glacial cycle. We argue that the tectonically-stable continental lithosphere can store elastic strain on long timescales, the release of which may be triggered by rapid, local transient stress changes caused by surface mass redistribution, resulting in the occurrence of intermittent intraplate earthquakes.

On 16 September 2015, the M-w 8.3 Illapel, Chile, earthquake broke a large area of the Coquimbo region of north-central Chile. This area was well surveyed by more than 15 high-rate Global Positioning System (GPS) instruments, installed starting in 2004, and by the new national seismological network deployed in Chile. Previous studies had shown that the Coquimbo region near Illapel was coupled to about 60%. After the M-w 8.8 Maule megathrust earthquake of 27 February 2010, we observed a large-scale postseismic deformation, which resulted in a strain rate increase of about 15% in the region of Illapel. This observation agrees with our modeling of viscous relaxation after the Maule earthquake. The area where upper-plate GPS velocity increased coincides very well with the slip distribution of the Illapel earthquake inverted from GPS measurements of coseismic displacement. The mainshock started with a small-amplitude nucleation phase that lasted 20 s. Backprojection of seismograms recorded in North America confirms the extent of the rupture, determined from local observations, and indicates a strong directivity from deeper to shallower rupture areas. The coseismic displacement shows an elliptical slip distribution of about 200 km x 100 km with a localized zone where the rupture is deeper near 31.3 degrees S. This distribution is consistent with the uplift observed in some GPS sites and inferred from field observations of bleached coralline algae in the Illapel coastal area. Most of aftershocks relocated in this study were interplate events, although some of the events deeper than 50 km occurred inside the Nazca plate and had tension (slab-pull) mechanisms. The majority of the aftershocks were located outside the 5 m contour line of the inferred slip distribution of the mainshock.

International audience; Megathrust earthquakes of magnitude close to 9 are followed by large-scale (thousands of km) and long-lasting (decades), significant crustal and mantle deformation. This deformation can be observed at the surface and quantified with GPS measurements. Here we report on deformation observed during the 5 yr time span after the 2010 M w 8.8 Maule Megathrust Earthquake (2010 February 27) over the whole South American continent. With the first 2 yr of those data, we use finite element modelling (FEM) to relate this deformation to slip on the plate interface and relaxation in the mantle, using a realistic layered Earth model and Burgers rheologies. Slip alone on the interface, even up to large depths, is unable to provide a satisfactory fit simultaneously to horizontal and vertical displacements. The horizontal deformation pattern requires relaxation both in the asthenosphere and in a low-viscosity channel along the deepest part of the plate interface and no additional low-viscosity wedge is required by the data. The vertical velocity pattern (intense and quick uplift over the Cordillera) is well fitted only when the channel extends deeper than 100 km. Additionally, viscoelastic relaxation alone cannot explain the characteristics and amplitude of displacements over the first 200 km from the trench and aseismic slip on the fault plane is needed. This aseismic slip on the interface generates stresses, which induce additional relaxation in the mantle. In the final model, all three components (relaxation due to the coseismic slip, aseismic slip on the fault plane and relaxation due to aseismic slip) are taken into account. Our best-fit model uses slip at shallow depths on the subduction interface decreasing as function of time and includes (i) an asthenosphere extending down to 200 km, with a steady-state Maxwell viscosity of 4.75 × 10 18 Pa s; and (ii) a low-viscosity channel along the plate interface extending from depths of 55-135 km with viscosities below 10 18 Pa s.

International audience; We model surface displacements induced by variations in continental water, atmospheric pressure, and nontidal oceanic loading, derived from the Gravity Recovery and Climate Experiment (GRACE) for spherical harmonic degrees two and higher. As they are not observable by GRACE, we use at first the degree-1 spherical harmonic coefficients from Swenson et al. (2008, https://doi.org/10.1029/2007JB005338). We compare the predicted displacements with the position time series of 689 globally distributed continuous Global Navigation Satellite System (GNSS) stations. While GNSS vertical displacements are well explained by the model at a global scale, horizontal displacements are systematically underpredicted and out of phase with GNSS station position time series. We then reestimate the degree 1 deformation field from a comparison between our GRACE-derived model, with no a priori degree 1 loads, and the GNSS observations. We show that this approach reconciles GRACE-derived loading displacements and GNSS station position time series at a global scale, particularly in the horizontal components. Assuming that they reflect surface loading deformation only, our degree-1 estimates can be translated into geocenter motion time series. We also address and assess the impact of systematic errors in GNSS station position time series at the Global Positioning System (GPS) draconitic period and its harmonics on the comparison between GNSS and GRACE-derived annual displacements. Our results confirm that surface mass redistributions observed by GRACE, combined with an elastic spherical and layered Earth model, can be used to provide first-order corrections for loading deformation observed in both horizontal and vertical components of GNSS station position time series.

International audience; We discuss the constraints on short-term asthenospheric viscosity provided by seasonal deformation of the Earth. We use data from 195 globally distributed continuous Global Navigation Satellite System stations. Surface loading is derived from the Gravity Recovery and Climate Experiment and used as an input to predict geodetic displacements. We compute Green's functions for surface displacements for a purely elastic spherical reference Earth model and for viscoelastic Earth models. We show that a range of transient viscoelastic rheologies derived to explain the early phase of postseismic deformation may induce a detectable effect on the phase and amplitude of horizontal displacements induced by seasonal loading at long wavelengths (1,300-4,000 km). By comparing predicted and observed seasonal horizontal motion, we conclude that transient asthenospheric viscosity cannot be lower than 5 × 10 17 Pa.s, suggesting that low values of transient asthenospheric viscosities reported in some postseismic studies cannot hold for the seasonal deformation global average.

Accepted for publication in Earth and Planetary Science Letters. 2017; The September 16, 2015 Mw8.3 Illapel Earthquake occurred on a locked segment of the South American subduction in Chile. This segment ruptured during comparable size earthquakes in the past, in 1880 and 1943, suggesting a somehow regular pattern of characteristic Mw 8+ earthquakes occurring every 60 to 80 years. This recurrence is in agreement with the accumulation of elastic deformation in the upper plate due to the Nazca-South America subduction at a constant rate of 6.5 cm/year, leading to a deficit of ~4.5 meters of slip to be released every 70 years. Previous studies consistently imaged the distribution of co-seismic slip along the fault based on geodetic, seismological and far field tsunami data and all described a significant amount of shallow slip resulting in a large tsunami. In addition, some models highlighted an apparent mismatch between the modeled rake of slip and the direction of plate convergence, suggesting the buildup of large strike-slip deficit. Some of these important questions remain open. Is shallow slip really well resolved and substantiated ? Is the apparent principal direction of slip during the earthquake really required by the geodetic data ? Here, using a comprehensive analysis of continuous GPS sites (including high rate and static displacements) and new survey data from acquired over more than 50 pre-existing sites, complemented with InSAR data, we show that the 2015 rupture overlaps very well the 1943 rupture, with the absence of significant slip south of 32°S and north of 30.2°S (peninsula Lingua de Vaca). Despite the wealth of geodetic data, the shallowest part of the subduction interface remains poorly resolved. We also show that the rake of the earthquake is fully compatible with the oblique plate convergence direction (rather than perpendicular to the trench), meaning that no subsequent trench-parallel motion is required by the data. We propose that the large Low Coupling Zone (LCZ) at the latitude of La Serena revealed by present day coupling distribution is stable over at least two seismic cycles. Inside the coupled area, peak coseismic slip is located precisely offshore the highest coastal topography and elevated terraces, adding weight to a potential correlation between the seismic cycle and long term permanent deformation. Finally, we show that early post-seismic after-slip occurs mostly down-dip of co-seismic asperity(ies), extending north and south of the 2015 rupture area.

Following previous findings from ongoing GPS research in Thailand since 2004 we continue to exploit the GPS technique to monitor and model land motions induced by the Sumatra–Andaman Earthquake. Our latest results show that up to the end of 2010, Thailand has been co-seismically displaced and is subsequently undergoing a post-seismic horizontal deformation with total displacements (co-seismic plus post-seismic) ranging from 10.5 to 74.7 cm. We observed the largest horizontal displacements in the southern part of Thailand and moderate and small displacements in the central and northern parts. In addition to horizontal displacements throughout Thailand, continuous GPS measurements show that large parts of Thailand are subsiding at rates up to 1 cm/yr. It is the first time that such vertical post-seismic deformations at large distances (650–1500 km away from the Earthquake’s epicentre) have been recorded. We have investigated the physical processes leading to the observed subsidence. While after-slip on the subduction interface induces negligible or even slightly positive vertical motions, relaxation in the asthenosphere is associated with a sizable subsidence. Predictions from a 3D finite element model feature an asthenosphere with an effective viscosity of the order of 3 * 10 18 Pas, fit the horizontal post-seismic data and the observed subsidence well. This model is then used to predict the subsidence over the whole seismic cycle. The subsidence should go on with a diminishing rate through the next two decades and its final magnitude should not exceed 10 cm in the Bangkok area. The post-seismic subsidence makes it difficult to identify other geophysical signals, particularly sea level rise, when observed from tide gauge data and thus there is a need for reliable estimation of subsidence velocities. This phenomenon may also worsen coastal erosion and flooding from sea water and so cause a considerable impact on the socio-economic development of coastal and low-lying areas, especially in the greater Bangkok area.

The deformations of the overriding and subducting plates during the seismic cycle associated with large subduction earthquakes are modelled using 2D and 3D finite element techniques. A particular emphasis is put on the interseismic velocities and on the impact of the rheology of the asthenosphere. The distance over which the seismic cycle perturbs significantly the velocities depends upon the ratio of the viscosity in the asthenosphere to the period of the seismic cycle and can reach several thousand km for rheological parameters deduced from the first years of deformation after the Aceh earthquake. For a same early postseismic velocity, a Burger rheology of the asthenosphere implies a smaller duration of the postseismic phase and thus smaller interseismic velocities than a Maxwell rheology. A low viscosity wedge (LVW) modifies very significantly the predicted horizontal and vertical motions in the near and middle fields. In particular, with a LVW, the peak in vertical velocity at the end of the cycle is predicted to be no longer above the deep end of the locked section of the fault but further away, above the continentward limit of the LVW. The lateral viscosity variations linked to the presence at depth of the subducting slab affect substantially the results. The north–south interseismic compression predicted by this preliminary 2D model over more than 1500 km within the Sunda block is in good agreement with the pre-2004 velocities with respect to South-China inferred from GPS observations in Thailand, Malaysia and Indonesia. In Japan, before the Tohoku earthquake, the eastern part of northern Honshu was subsiding while the western part was uplifting. This transition from subsidence to uplift so far away from the trench is well fitted by the predictions from our models involving a LVW. Most of the results obtained here in a 2D geometry are shown to provide a good estimate of the displacements for fault segments of finite lateral extent, with a 3D spherical geometry, as long as the displacements during the seismic cycle are normalised by the coseismic displacement.

Seismic studies of the lowermost mantle suggest that the core-mantle boundary (CMB) region is strongly laterally heterogeneous over both local and global scales. These heterogeneities are likely to be associated with significant lateral viscosity variations that may influence the shape of the long-wavelength non-hydrostatic geoid. In the present paper we investigate the effect of these lateral viscosity variations on the solution of the inverse problem known as the inferences of viscosity from the geoid. We find that the presence of lateral viscosity variations in the CMB region can significantly improve the percentage fit of the predicted data with observations (from 42 to 70% in case of free-air gravity) while the basic characterisics of the mantle viscosity model, namely the viscosity increase with depth and the rate of layering, remain more or less the same as in the case of the best-fitting radially symmetric viscosity models. Assuming that viscosity is laterally dependent in the CMB region, and radially dependent elsewhere, we determine the largescale features of the viscosity structure in the lowermost mantle. The viscosity pattern found for the CMB region shows a high density of hotspots above the regions of higher-than-average viscosity. This result suggests an important role for petrological heterogeneities in the lowermost mantle, potentially associated with a post-perovskite phase transition. Another potential interpretation is that the lateral viscosity variations derived for the CMB region correspond in reality to lateral variations in the mechanical conditions at the CMB boundary or to large-scale undulations of a chemically distinct layer at the lowermost mantle.

SUMMARY Numerical simulations of 2-D Rayleigh-Benard convection are designed to study the devel- opment of convection at the base of the cooling lithosphere. A zero temperature is suddenly imposed at the mantle surface, which has initially a homogeneous temperature. The strongly temperature-dependent viscosity fluid is heated from within, in order to balance the internal temperature drift resulting from global fluid cooling. For a while, the lithosphere cools approx- imatively as a conductive half-space and the lithospheric isotherms depth grows as the square root of age. As instabilities progressively develop at the base of the lithosphere, lithospheric cooling departs from the half-space model. We propose two different parametrizations of the age of the first dripping instability, using boundary layer marginal stability or quantifying the characteristic timescale of the instability exponential growth as a function of the Rayleigh number and of the viscous temperature scale. Both parametrizations account very well for our numerical estimates of onset times, but with a slightly better adjustment of the viscous tem- perature scale dependence in the second case. The absolute value of the onset time depends on the amplitude and location of initial temperature perturbations within the box and on the initial temperature structure of the thermal boundary layer (TBL). Furthermore, thermal perturba- tions of finite amplitude located within the lithosphere (such as the ones induced by transform faults, for example) strongly reduce the age of the first dripping instability. However, the onset time parametrization derived from transient cooling experiments well adjusts the lithospheric age of the first drip instability below a lithosphere cooling perpendicularly to the ridge axis. This study emphasizes the role of the initial topography of the lithospheric isotherms on the development of instabilities.

SUMMARY We have investigated the impact of lateral viscosity variations (LVV) in the top 300 km of the mantle on the long-wavelength gravitational response of the Earth. In contrast to the previous studies we demonstrate that the LVV may play a crucial role if a model with imposed plate velocities and partial layering is considered. Assuming LVV in the top mantle associated with the continental roots and the oceanic asthenosphere, we are able to explain a significant portion of the geoid (85 per cent) and of the free-air gravity data (55 per cent) without invoking a complex radial viscosity profile or assuming a lithosphere with an unrealistically low uniform viscosity. The best fit to the data is found for models in which mass anomalies located at the boundary between upper and lower mantles reduce the mass flux across the 660 km discontinuity by a factor of 3. These models are characterized by a low average value of viscosity in the asthenosphere below the oceans (∼10 19 Pa s) and a large viscosity contrast (>100) between the young suboceanic mantle and the deep continental roots. The viscosity in the deep upper mantle is found to be around 3 × 10 20 Pa s while a value higher by at least two orders of magnitude is typical for the lower mantle. The predicted dynamic surface topography is influenced by the surface plate velocities and its amplitude amounts to a few hundred metres. Besides investigating the viscosity structures that reflect the different thicknesses of the lithosphere in continental and oceanic regions, we have searched for a large-scale pattern of lateral viscosity variation in the top 300 km that yields the best fit to the geoid. The pattern obtained by this inversion shows a good correlation with the ocean‐continent distribution and 90 per cent of the geoid and 67 per cent of the free-air gravity can be explained.

The stretching and mixing of blobs more viscous than the surrounding mantle is an important issue: The subducted lithospheric material at depth might be more viscous than the surrounding mantle, and its ability to mix affects its properties as a geochemical reservoir. In order to study the effect of the viscosity ratio λ (viscosity of the blob to the viscosity of the surrounding mantle) on the deformation we perform three types of two-dimensional numerical simulations: a pure shear simulation, a simulation of steady state convection, and a simulation of time-dependent convection. In a pure shear flow the inverse of the rate of deformation of a viscous blob is, to a first approximation, a linear function of the ratio λ/AR, where AR is the aspect ratio of the blob. For a viscous blob in a steady convection cell a periodic deformation is superimposed on the long-term deformation. The period corresponds to an overturn time along the streamline followed by the blob. The strains of this periodic deformation obey the laws found for pure shear. In a time-dependent flow, more representative of the mantle flow, we find that the average distance a* between two tracers situated within a viscous blob obeys the relationship ln(a*/a0)=At/(1+λ)2 and that it is independent of the aspect ratio of the blob. These two results are well explained by a simple model of successive strains in arbitrary directions. The viscosity contrast controls the mixing efficiency, and it has an effect much larger than predicted by earlier studies. Blobs of matter more viscous than the average mantle can therefore persist unmixed in the mantle over Earth's history.

SUMMARY The existence of continental roots and the gradual thickening of the cooling oceanic lithosphere give rise to large-scale rheological heterogeneities in the uppermost mantle. The effect of these heterogeneities on the long-wavelength geoid is investigated using a 3-D mantle flow model involving a low-viscosity asthenosphere beneath the oceanic lithosphere and the tectonically active continental regions and a thick highly viscous lithosphere beneath ancient continents. Below 400 km the mantle viscosity is laterally homogeneous with a lower mantle more viscous than the overlying layer. The mantle circulation is driven by imposed surface velocities NUVEL-1 HS2 (Gripp & Gordon 1990) and by the density anomalies inferred from the tomographic models P16B30 and S16B30 (Masters et al. 1996). The geoid heights both due to plate motion and due to internal loading differ by as much as 30 per cent between the models with and without lateral viscosity variations. In contrast to what was suggested in previous studies, spherical harmonics 2 and 3 are strongly affected by the lateral viscosity variations. These differences in the forward problem suggest that the response to the inverse problem that consists of finding the profile of viscosity as a function of depth providing the best fit to the geoid should be considerably affected by the lateral viscosity variations. The shear stresses at the base of the plates induced by the imposed surface velocities and those induced by the internal loading are sensitive to the lateral viscosity variations. This suggests that the lateral viscosity variations are very important for understanding the stress field or the forces acting on the plates. The direction of these shear stresses, which should be linked to the anisotropy direction, is very different from the plate motion direction.

SUMMARY The origin of the flattening of the seafloor bathymetry and geoid at old ages is still a matter of debate. Is this departure from the square-root-of-age law induced by dynamic uplift of the lithosphere linked to the flow in the mantle or by cooling of the lithosphere different from that predicted for an infinite half-space? The amplitude of the dynamic uplift should strongly decrease at short wavelengths and should not be correlated with the seafloor age on a local scale. Therefore, in order to eliminate the ‘dynamic signal’, the cooling of the lithosphere is studied here by means of depth–age and geoid–age relationships computed from high-pass filtered bathymetric, geoid and age fields, or from local dtopography/d and dgeoid/dage slopes. The geoid–age and depth–age relationships obtained here are in good agreement with those previously deduced from global non-filtered data. At all observation scales, a significant departure from the half-space solution is found. However, the data are not sufficiently accurate to discriminate between the Plate and CHABLIS (Constant Heat flow At the Base of the Lithospheric ISotherm) models. Our results support the existence of small-scale convection supplying heat at the base of the old oceanic lithosphere.

The transformation of olivine to spinel in a subducting slab gives rise to thermal strain, transformation strain, and buoyancy forces inside the slab. We use a finite element method to solve the mechanical equation, and we follow the evolution of the stresses inside the slab during its descent. We use two different rheologies: elastic and viscoelastic with a viscosity that depends on pressure, temperature, and grain size. The stresses due to buoyancy are found to be negligible compared to the internal stresses due to the phase transformation. When the viscoelastic rheology is used, a zone of very high differential stress coincides with the olivine-spinel phase boundary. This zone is associated with large downdip compressive stresses inside the metastable olivine wedge, and large downdip tensile stresses in the spinel region, along the wedge. If we assume that the failure criterion can be related to the differential stress, our model is in agreement with the patterns of deep seismicity. However, patterns of deep seismicity are unlikely to reflect buoyancy stresses, which in the presence of a metastable olivine wedge, are negligible compared to the internal stresses.

Most inversions of the long-wavelength geoid in conjunction with the seismic tomographic information have so far been carried out under the assumption of either purely whole mantle or perfectly layered circulation. Moreover, modeling the lithosphere as a spherical shell with a uniform low viscosity was found to yield the best fit to the observed geoid. We have tested whether a good prediction of the geoid can also be achieved by including two constraints: a semipermeable behavior of the 660-km discontinuity and surface plate velocities equal to the observed ones. The mass transfer between upper and lower mantle has been changed by imposing a surface density anomaly at a depth of 660 km, which is proportional to the mass anomaly needed to achieve perfectly layered circulation. On the top of the mantle we assume a stiff lithosphere that moves with a velocity corresponding to the observed plate motion. The viscosity only varies with depth. Considering a simple three-layer viscosity structure and changing the permeability of the 660-km interface, we have obtained a satisfactory variance reduction of the geoid data (∼75% for degrees 2–12). The best fit to the geoid is obtained if the mass transfer across the 660-km boundary is reduced to one third in comparison with the purely whole mantle model. The best fitting viscosity profile is characterized by a clearly defined asthenosphere and a viscosity increase by at least 2 orders of magnitude in the lower mantle. The amplitudes of dynamic topography predicted by our model are remarkably small (∼100 m), thus fully compatible with the observation.

The paper presents a denudational rift flank uplift model of the Western Ghats of Karnataka, India. The Cenozoic denudation of the Deccan plateau was constrained by relative dating of regional planation surface levels, in combination with preliminary apatite fission track results. The denudational history of the Western Ghats escarpment coastal foreland was constrained both by onshore fission track and offshore sediment data. Methods are briefly described. Results were used as reference data to elaborate a computer simulation model of both flexural rift flank upwarp and escarpment retreat from the K/T boundary to the present. Depending on the assumptions concerning the attributes of the thin elastic lithospheric sheet (infinite or semi-infinite), the shape and elevation of the initial Deccan topography (flat-lying, or affected by lithospheric updoming due to the Reunion mantle plume) and the position of the continental divide upon rifting, the flexural response to denudational unloading is shown to yield either a concave-up flexure of the Dharwar craton, or a convex, monoclinal downwarp of the margin. The relative merits of each model are discussed in the light of auxiliary field evidence, particularly from the neighbouring Deccan traps. However, the diverse assumptions on boundary conditions fail to entirely avoid geomorphological equifinality. © 1998 John Wiley & Sons, Ltd.

SUMMARY In order to gain some insights into the role of deep mass anomalies on surface velocity and geoid anomalies, Cartesian models involving strong lateral viscosity contrasts and where plate boundaries are mimicked by zones of very low viscosity are performed. The response of the system to mass anomalies at various depths and various distances from the plate boundary is computed. It is compared with the responses obtained using the three other methods used in global geoid models for mimicking the mobility of the lithosphere: models with layered viscosity and a uniform low-viscosity lithosphere; models computing plate velocities from horizontal force balance; and models where the surface velocities are set equal to the observed ones with plates ‘moved by the hand of God’. Contrary to what occurs in models with layered viscosity, the geoid response for the weak-zone model proposed here changes drastically from positive to negative, depending upon the position of the mass anomaly with respect to the plate boundary. This arises from the fact that only mass anomalies situated close to the plate boundary are able to move the plates efficiently. Changing the location of an upper-mantle mass anomaly by a few hundred kilometres considerably modifies the amplitude of the predicted plate velocities and geoid anomalies. Consequently, upper-mantle density anomalies deduced from global tomography analysis, or density anomalies assumed to be known a priori (subducted slabs) but truncated by a spherical harmonic decomposition should not be used in models aimed at computing plate motions and the associated geoid. When a mass anomaly is situated just below a plate boundary, the sign of the predicted geoid anomalies can be the same as for models with a layered viscosity and a low-viscosity lithosphere. However, the shape of the predicted geoid anomaly and its amplitude as a function of the depth of the mass anomaly cannot be mimicked by models with a uniform-viscosity lithosphere, whatever the lithospheric viscosity. For mass anomalies situated at a large distance from the plate boundaries, the induced plate velocities are small and the geoid is equal to that computed for a highly viscous continuous lithosphere. The ‘force-balance method’ is strongly affected by the cut-off in the harmonic decomposition. Only models with a lithospheric layer below the surface reproduce the geoid successfully for mass anomalies situated away from a plate boundary. The plate velocities for the force-balance method and for the weak-zone model proposed here are not proportional, due to differences in the structure of the plate boundary. In short, only mass anomalies very close to a plate boundary are able to influence the plate velocities, but these velocities depend strongly on the mechanical details of the plate boundary. As a consequence, we believe that it is impossible to compute plate velocities from what we know from mass anomalies in the Earth precisely enough to include them in a geoid model. On the other hand, actual plate velocities are known. We show that an accurate geoid can be computed by imposing the observed plate velocities onto the Earth's surface, even if the details of the mechanical structure of the lithosphere are not respected. It is this method that is suggested for introducing the effect of plate tectonics in global geoid models.

Geoid anomalies in the wavelength range corresponding to spherical harmonies 6 to 30 are used to bring new constraints on the deep structure of the continental lithosphere, on the nature of upper mantle seismic anomalies, and on the thermal evolution of the oceanic lithosphere. For the shallow lithospheric density anomalies considered here, the geoid anomaly is shown to follow approximately the moment law derived from isostasy. By a procedure of least squares inversion, the geoid over old cratons is found to be equal to that over relatively young (50 Myr) oceanic lithosphere. The oceanic geoid as a function of age follows the expected linear trend with flattening at old ages. These results are independent of the degree at which the spherical harmonics are truncated. If the tectosphere consisted only of a thick thermal boundary layer, marked geoid lows, with an amplitude larger than 20 m, should be present over cratons, in disagreement with the present observation. The predicted geoid agrees with the observations if the negative density anomalies associated with a depleted peridotitic layer are present to the bottom of the lithosphere. This analysis of the geoid confirms the petrological inferences that the deep roots observed seismically in the top upper mantle form a chemically differentiate reservoir. The density of these cold roots is similar to that of mantle material elsewhere, because the increase in density resulting from their lower temperature is largely balanced by the decrease in density caused by their depletion. Analysis in wavelength range corresponding to harmonics 6 to 30 in the oceans confirms that the geoid flattens at old ages, and thus supports the view that an important amount of heat is transferred at the base of the lithosphere by small-scale convection.

The most common model used for representing the evolution with age of the oceanic lithosphere is the ‘plate model’ where the temperature is set at a fixed depth, called the base of the plate. This ‘base of the plate’ has no physical meaning but this model provides a mathematical substitute for a system where small-scale convection occurs through instabilities growing at the base of the cooling lithosphere and becomes effective only below old ocean. Another possible view is that convection provides heat at the base of the lithosphere whatever the age of the overlying plate. This last process can be modeled by a Constant Heat flow Applied on the Bottom Lithospheric ISotherm (CHABLIS model). A good fit to the observables (bathymetry and geoid as function of age, and old age heat-flow) can be obtained both for plate and CHABLIS models in spite of an experimentally determined thermal expansion coefficient much larger than assumed in previous plate models. These models have important consequences for several geodynamic processes. The plate, at an age of 100 Ma is only 80 km thick for both models: melting above a hot-spot can then occur in the garnet-spinel transition field without much plate thinning. In the plate model the subsidence is stopped at an age of about 80 Ma while, according to the CHABLIS model, several hundred meters of subsidence are expected after 100 Ma. Thus the two models predict quite a different long-term pattern of subsidence in the sedimentary basins. Finally, in the CHABLIS model, the global cooling of the mantle coming from cold material eroded by secondary convection at the base of the plates is considerably larger than in plate models: it amounts to 40%, the remaining 60% being due to the subduction process.

The two earthquakes of respective magnitudes 8.6 and 8.2 that occurred off the coast of the Sumatra subduction zone on 11 April 2012 are shown to be part of a continuing boost of the intraplate deformation between India and Australia that followed the Aceh 2004 and Nias 2005 megathrust earthquakes. On 11 April 2012, two of the largest strike-slip earthquakes ever recorded — at magnitudes of 8.7 and 8.2 — occurred in the northeastern Indian Ocean, a few hundred kilometres off the coast of Sumatra. Three groups now report the analysis of seismic data from the days and months before and after these events, as well as the events themselves. Matthias Delescluse and co-authors show that these earthquakes are part of the ongoing boost of intraplate deformation between India and Australia that followed the Aceh 2004 and Nias 2005 megathrust earthquakes. They conclude that the Australian plate, driven by slab-pull forces at the Sunda trench, is gradually detaching from the Indian plate. Han Yue and colleagues show that the 11 April event involved a complex four-fault rupture lasting several minutes, followed two hours later by a magnitude-8.2 aftershock. These great ruptures on a lattice of strike-slip faults that extends through the crust and into the upper mantle represent large lithospheric deformation that may eventually create a localized boundary between the Indian and Australian plates. Fred Pollitz and colleagues show that, in the six days following 11 April, the global rate of remote earthquakes with magnitudes greater than 5.5 increased nearly fivefold, and events up to magnitude 7 seem to have been triggered. The unprecedented delayed triggering power of this earthquake may arise from its strike-slip source geometry, or because it struck at a time of an unusually low global earthquake rate and increased the number of nucleation sites that were very close to failure. Large earthquakes nucleate at tectonic plate boundaries, and their occurrence within a plate’s interior remains rare and poorly documented, especially offshore. The two large earthquakes that struck the northeastern Indian Ocean on 11 April 2012 are an exception: they are the largest strike-slip events reported in historical times1,2 and triggered large aftershocks worldwide3. Yet they occurred within an intra-oceanic setting along the fossil fabric of the extinct Wharton basin, rather than on a discrete plate boundary4,5,6,7,8. Here we show that the 11 April 2012 twin earthquakes are part of a continuing boost of the intraplate deformation between India and Australia that followed the Aceh 2004 and Nias 2005 megathrust earthquakes, subsequent to a stress transfer process recognized at other subduction zones9,10. Using Coulomb stress change calculations, we show that the coseismic slips of the Aceh and Nias earthquakes can promote oceanic left-lateral strike-slip earthquakes on pre-existing meridian-aligned fault planes. We further show that persistent viscous relaxation in the asthenospheric mantle several years after the Aceh megathrust explains the time lag between the 2004 megathrust and the 2012 intraplate events. On a short timescale, the 2012 events provide new evidence for the interplay between megathrusts at the subduction interface and intraplate deformation offshore. On a longer geological timescale, the Australian plate, driven by slab-pull forces at the Sunda trench, is detaching from the Indian plate, which is subjected to resisting forces at the Himalayan front6,8,11.

[1] The volume changes associated with a phase transition need to be accommodated by creep. Otherwise, a pressure anomaly impeding the transformation is induced. Phase changes in mantle flows are therefore delayed by this process, which occurs not only on a microscopic scale but also, as shown in the present paper, on a macroscopic (i.e., kilometric) scale. The mechanism is illustrated first with the help of a simple analytical “pipe flow” model. The broadening and average deflection of the phase transition vary with the vertical velocity of the flow V, the viscosity in the zone of phase transition η, the density jump Δρ, and the pressure interval over which the phase change occurs Δ P. They are characterized by a nondimensional number: α = (4/3)ηgV (Δρ/ΔP2). For small α, the deflection of the phase transition is equal to α times the width of the transition loop. A numerical, one-dimensional stationary model is then used to quantify how the flow is hampered by mass anomalies associated with the phase transition deflections. The response functions for harmonic loads in a self-gravitating mantle are computed. Transition deflection is shown to have a significant influence in the case of a thin and viscous discontinuity. Conversely, it may be negligible if the material within transition zones is weakened. In an estimate where we link the viscosity in the transition zone to the level of deviatoric stress, we show that sizable (≈7 km) deflections of the 670 km discontinuity can be expected as a consequence of these macroscopic volume changes associated with phase transformations.

SUMMARY The geoid and bathymetry data in the Pacific Ocean are used to shed new light on the structure of upper mantle convection. At short wavelengths (200 km), the geoid pattern is characterized by parallel lineations. Lineations however appear less clearly at larger wavelengths (400 or 600 km) where the signal due to convective thermal anomalies might be expected. The most prominent features among the short-wavelength lineations are associated with the well-known chains of seamounts and islands. Other lineations observed in the geoid coincide with computed tracks of other volcanic centres. Although these tracks are not continuous, this suggests that long-lived plumes are at the origin of the apparently short chains of seamounts. Apart from the lineations corresponding to the tracks of ‘classical’ hotspots, there are numerous similar lineations visible on the filtered geoid maps which have the same directions as hotspot tracks. These preferential directions and the gravity-to-bathymetry ratios seem to indicate that the lineations are due to magmas implaced just above or within the lithosphere. Convective plumes, fixed with respect to each other and much more numerous than the classical hotspots could be the source of these magmas. A geometrical test of this hypothesis is performed by comparing the lineation pattern with the computed tracks for a set of tentative new hotspots. In this way, numerous plumes, long-lived and fixed with respect to each other are found to provide an explanation of the observed small-wavelength geoid features.

The 2004 Aceh and 2005 Nias events are the two greatest earthquakes of the past 40 years with a total rupture of 1700 km long and a coseismic slip reaching up to 25 m. These two earthquakes have caused large stress perturbations which significantly altered seismic activity in the Sumatra-Andaman region. Using both detailed mapping of failure planes and various slip distributions, we calculate this stress change along the Sumatra-Andaman-Sagaing fault system from central Sumatra to southern Myanmar. The static Coulomb stress change Delta CFF and the observed seismic activity are in very good agreement with a Coulomb index similar to 20% greater than the one obtained for random events. Compared to previous studies, this high Coulomb Index confirms two important issues on the use of static stress change criterion: unsuited to study near-field aftershocks and only relevant for aftershocks analysis on large and mature faults at a time scale of several months. The calculated Delta CFF distribution suggests that the 2004 and 2005 earthquakes inhibit failure on the North Andaman rift and on the Sagaing fault, while failure is encouraged along the transform Andaman zone, the central Andaman rift, the West Andaman fault, the Sumatra fault system, and the offshore thrust faults west of Sumatra Island. The maximum value of similar to 15-20 bar (1.5-2 MPa) for Delta CFF is reached in the northern part of the Sumatra fault system. This high value together with the lack of major earthquake in the last 170 years result in a high seismic hazard for this region. Our results are also consistent with temporal evolution of both earthquakes' location and focal mechanism prior to and after the events. In particular, we explain the occurrence and the mechanism of seismic swarms observed in the central Andaman rift and along the west Andaman fault. Finally, our calculations reveal that the seismicity in the Andaman rift zone can only be explained if mu' > 0.5. This result leads to two end-member models: one with a constant and high fault friction and one with spatial variations, for which friction may depend on either the nature of the lithosphere (oceanic versus continental) or the fault type.

The stresses associated with large-scale tectonic deformation have three possible origins: (1) plate-boundary forces counterbalanced by viscous drag beneath the plates; (2) density heterogeneities situated within the plates (say at depths shallower than 200 km); (3) mass heterogeneities in the deep mantle. The first two are shown to be equally important for the understanding of the stress field. No topography (no vertical stress) seems to be associated with lower-mantle mass anomalies. This is most compatible with a two-layer convective mantle where the lower-mantle mass anomalies, mechanically decoupled from the lithosphere, are unable to induce tectonic stresses.

The sources of geoid and topography anomalies in a convecting upper mantle with realistic temperature- and pressure-dependent viscosity are studied. The results differ from those obtained in previous studies with layered, depth-dependent viscosity. For a Newtonian viscosity, the topography and the geoid anomalies coming from below the low-viscosity zone (i.e., deeper than 250 km) amount to 30 m and 10 cm, respectively, while they are reduced to 10 m and 10 cm for a non-Newtonian viscosity. It represents at most 6% of the complete signal and affects very little the geoid over topography ratio. The weakness of the contribution from the lower part of the cell is due to both small kernels and small lateral temperature variations at depth. Yet, a linear regression of the complete geoid signal versus topography yields an apparent compensation depth shallower than the base of the lithosphere. The damping of the vertical stress through a lithosphere with a depth-dependent viscosity implies a departure from isostasy and is the cause of the shallow apparent compensation depths found in our model.

[1] Short- and intermediate-wavelength gravity and geoid anomalies are used to provide constraints on the mechanical structure of subduction zones and on the forces involved. This study is based on 2D Cartesian dynamically self-consistent models with Newtonian or power law rheologies. We show that both strong decoupling of the two convergent plates (shear stresses of the order of 107 Pa) and weakened bending lithosphere are necessary to reproduce the observed geoid and gravity data. Good fits are found for relatively low failure stresses (≈30–50 MPa). For all models providing reasonable predictions of gravity, only a small fraction of the downgoing slab weight is transmitted to the surface plates. About 10% of the energy is dissipated in the contact zone between the two plates, 10% to 20% in the bending region, and more than 70% in the sublithospheric mantle. The basal tractions (on the order of ≈1−4 × 1012 N/m) induce a net motion of the plates, with the subducting lithosphere moving faster than predicted by the no-net-motion principle. A marked positive geoid anomaly is predicted above subduction zones at intermediate wavelength (λ = 2000–4000 km) in the case of pure whole mantle convection. Such large geoid highs are not observed. Introducing “partial layering” (i.e., mass anomalies hampering mantle flow through the transition zone) is necessary to reconcile model predictions and observations for these wavelengths.

Topography associated with large wave-length geoid is derived by two methods. First, the bathymetry corrected for sediment loading is decomposed using a least-squares method as a function of age plus terms proportional to the geoid. The topography at old ages is flatter predicted by the law. At degrees 2 and 3, very little topography is associated with the geoid. The transfer function is not even certain in sign but is never smaller than 0.6. More topography seems associated with higher degree harmonics (6 and 7). Second, the function of age alone best fitting the topography is determined. The correlation coefficients of the residual topography with the geoid are then computed. A slight correlation appears for degree 6, but not for degrees 2 and 3. These results suggest a, mechanical decoupling between the lower mantle mass heterogeneities and the lithosphere.

An infinite horizontal layer, with vertically stratified temperature and solute concentration, is considered in the case where the viscosity is exponentially dependent on temperature, and the Prandtl number is infinite. Its linear stability is investigated when the destabilizing thermal gradient acts against a stabilizing solute gradient. The analysis is performed using horizontal Fourier and vertical Chebyshev polynomial expansions. For the constant viscosity case, the laws well established in the free boundary configuration are seen to be directly suitable for the rigid one. In the variable viscosity case, characterised by a given viscosity contrast c, the scaling laws with c are settled extrapolating to the double diffusive situation the approach initiated by Stengel et al. (1982). In contrast with the constant viscosity case, the critical wave number is found to be strongly dependent on the solutal Rayleigh number in the marginal oscillatory obtained at large contrast values.

International audience; Numerical simulations of two-dimensional Rayleigh-Bénard convection are designed to study lithospheric cooling above a convective mantle. A strongly temperature-and pressure-dependent viscosity fluid is heated from below or from within. An imposed velocity at the surface of the box mimicks the plate motion between the ridge on one side and the subduction zone on the other side. As the lithosphere cools, its upper part remains rigid and therefore conductive, while its bottom part is convectively unstable. Dripping instabilities are not observed close to the ridge. Nevertheless, the material flows along the slope defined by the lower part of the lithosphere and feeds the first descending drip. Afterwards, cold downgoing instabilities develop continuously and randomly at the base of the lithosphere and are replaced by hot material from the convecting core of the box. The lithosphere continues to thicken even after the onset of the first instability. Surface heat flow, subsidence and lithospheric temperature structure obtained by the convective simulations are compared to the predictions of three conductive models: the Plate, Chablis, and modified Chablis models. These models differ by their applied bottom boundary condition which represents the lithosphere/asthenosphere convective coupling, i.e. by the presence or absence of instabilities developing at the base of the lithosphere. The conductive model which best explains the lithospheric cooling obtained by convective simulations is the modified Chablis model. In this model, a variable heat flow (depending upon the viscosity at the base of the lithosphere) is applied along an isotherm located in the lower unstable part of the lithosphere.

SUMMARY Several studies have found that the ‘subsidence rate’ of the oceanic lithosphere measured along flow lines is a linear function of the height of the ridge. This might have very important implications concerning our understanding of the cooling of the lithosphere. Here we study the relationship between the subsidence rate and the height of the ridge in the South Atlantic and East Indian Oceans. From the regression slopes of this relationship, it has been inferred that lithosphere created at shallower ridges subsides faster and to deeper levels than that generated at deeper ridges. We show however, that the topography at young ages for various seafloor segments stays positively correlated with the topography at older ages. Synthetic data sets are obtained by adding random synthetic signals with the same spectrum as the residual topography to the average depth versus age relationship. They are also characterized by the same linear relationship between ‘subsidence rates’ and the height of the ridge. This relationship is therefore not proof that lithosphere formed at shallow ridges cools faster than lithosphere formed at deeper ridges. Random signals linked with crustal thickness variations, hotspot bulges and perhaps large wavelength ‘dynamical topography’ would also yield such a relationship.

We have tested the impact of lateral viscosity variations (LVV) in the top 200 km of the mantle (tectosphere) on the interpretation of Fennoscandian postglacial uplift data in terms of 1D viscosity models. Our test consists of two steps. In the first step, we simulate the postglacial uplift for a given 2D viscosity structure involving a cratonic lithospheric root below the ice sheet and an asthenosphere below oceanic areas. In the second step, we interpret the result of this simulation in terms of a 1D viscosity model. The resultant 1D viscosity model shows a lower viscosity in the lower mantle, a higher viscosity in the shallow sublithospheric mantle and a smaller thickness of the lithosphere than the synthetic model of viscosity.

International audience; A numerical model of two‐dimensional Rayleigh‐Bénard convection is used to study the relationship between the surface heat flow (or Nusselt number) and the viscosity at the base of the lithosphere. Newtonian or non‐Newtonian, temperature‐ and pressure‐dependent rheologies are considered. In the high Rayleigh number time‐dependent regime, calculations yield Nu ∝ RaBL1/3beff−4/3 where beff is the effective dependence of viscosity with temperature at the base of the upper thermal boundary layer and RaBL is the Rayleigh number calculated with the viscosity νBL (or the effective viscosity) at the base of the upper thermal boundary layer. The heat flow is the same for Newtonian and non‐Newtonian rheologies if the activation energy in the non‐Newtonian case is twice the activation energy in the Newtonian case. In this chaotic regime the heat transfer appears to be controlled by secondary instabilities developing in thermal boundary layers. These thermals are advected along the large‐scale flow. The above relationship is not valid at low heat flow where a stationary regime prevails and for simulations forced into steady state. In these cases the Nusselt number follows a trend Nu ∝ RaBL1/5beff−1 for a Newtonian rheology, as predicted by the boundary layer theory. We argue that the equilibrium lithospheric thickness beneath old oceans or continents is controlled by the development of thermals detaching from the thermal boundary layers. Assuming this, we can estimate the viscosity at the base of the stable oceanic lithosphere. If the contribution of secondary convection to the surface heat flux amounts to 40 to 50 mW m−2, the asthenospheric viscosity is predicted to be between 1018 and 2×l019 Pa s.

[1] We investigate the effect of an endothermic phase transition on the average petrology of the mantle as function of depth. We study chemical differentiation using a numerical model in which the oceanic lithosphere, comprised of a crustal layer overlying depleted harzburgite, is injected at the ridges. Three families of tracers are used, each one with its own buoyancy, to distinguish oceanic crust, harzburgite and peridotite. The depth of the phase transition and the Clapeyron slope are assumed to be the same whatever the petrological component of the mantle. Velocities are imposed at the surface to mimic plate tectonics. The viscosity is depth and temperature dependent. The negative Clapeyron slope tends to induce a segregation of the heavy material below the phase transition and of the light material above the phase transition with a compositional jump across the transition zone.

Summary. Narrow zones of intense shear deformation, i.e. viscous slip zones, are studied analytically with a one-dimensional time-dependent model of two half-spaces of identical or contrasting rheologies and ambient temperatures in relative motion. The rheologies of the half-spaces are strongly temperature-dependent and viscous heating maintains a thin zone of high temperature, low viscosity and large strain rate. The mathematical model is used to describe the structures of slip zones at ridge and plate-boundary transform faults, major continental strike—slip faults and at the top of subducting oceanic crust. No a priori assumption about slip-zone width or shear-stress magnitude is necessary; the thermal-mechanical structure of the slip zone evolves in time and all its characteristics are self-consistently determined. Slip-zone widths and shear stresses depend on the ambient temperatures, the relative velocity, the rheology and the length of time following the onset of relative motion; for reasonable geologic times, 0.1–10 M yr for example, slip zones are generally several kilometres wide and shear stresses are several hundred bars (tens of MPa). The region of intense shear in a viscous slip zone is an order of magnitude narrower than the width of the accompanying thermal anomaly. The maximum temperature generated by viscous dissipation in a slip zone depends only on the relative motion and the creep properties of the rocks; it is independent of slip-zone age and ambient temperature. Maximum temperatures associated with frictional heating are always less than those required for partial melting. The slip zone on a descending slab is influenced most strongly by the contrast in creep behaviour between the relatively soft oceanic crustal rocks and the hard, overlying mantle rocks; as a result, the slip zone is confined entirely within the oceanic crustal layer. Oceanic crustal rocks deform so readily that frictional heating in a slip zone on a descending slab cannot by itself lead to partial melting and thermal conduction from the hotter overriding mantle must play an essential role in heating the descending crust if melting is to occur therein. Because of the increase in the mantle temperature with depth, narrow slip zones in oceanic regions probably do not exist below depths of about 100 km; they may extend to greater depths beneath continents.

Both morphologic and tectonic studies indicate that Mercury and the Moon have quite different internal histories, despite their apparently similar morphologies. The evaluation of the volcanic surfaces indicates a decreasing volcanism on Mercury at the largest impacting time, despite short and local reactivations. On the Moon, the basaltic volcanism was increasing at the same time and continued for 1 billion years. That indicates a strongly different thermal evolution for these two planetary bodies. A widespread graben pattern is present on the Moon, with a statistical dominance of radial or tangential orientation with respect to the Imbrium basin, thus suggesting a relation between this major basin and the expansion of the Moon. Azimuthal studies show that the compressive structures, observed on the stereographic covered surface of Mercury are not randomly oriented, but seem radial towards the Caloris basin, thus indicating a possible influence of this largest basin on Mercurian contraction. The qualitative and quantitative formulations of these tectonic perturbations induced by large basins will be developed in a companion paper [1].

Gravitational, topographical and seismological studies of uplift regions around hot spots on oceanic,1,2 and continental plates3–7 have suggested that considerable lithospheric thinning has taken place beneath these areas. Previous models used to study this phenomenon have been based on thermal conduction8–11 and on the dynamic interaction between large-scale plume structures12 and the lithosphere, for a purely temperature-dependent viscosity. However, their predictions of the time scale for direct erosion12 of the lithosphere by these upwellings are too long. We propose here a new dynamical mechanism which works most efficiently for a temperature and pressure-dependent rheology. Our idea is that thinning of the lithosphere is accomplished by strong secondary convection enhanced by the presence of a low-viscosity zone which grows from a plume associated with a hot spot. This process is capable of producing the rapid time scale of the swell uplift13, of the order of 10 Myr, and the rates of thinning, of order 10 km Myr−1, required by both geophysical13,14 and geochemical15 observations.

Geoid observations in the oceans suggest that lithospheric thinning is effected by small-scale convection within hot-spot material trapped in dimples a few tens of kilometers deep at the base of the lithosphere. Significant thinning can occur within 5 Ma if the viscosity of the convecting material is 10 16 Pa s. Partial melting can enhance considerably the vigour of convection and the same rate of lithospheric thinning is obtained for a viscosity about five times higher. These results are derived from convection models with a realistic temperature-dependent rheology using the mean field approximation.

The purpose of this paper is to clarify the dynamic role of lithospheric density heterogeneities, in particular with respect to mountain building and other processes of intraplate deformation. Density anomalies within or just beneath the lithosphere constitute major sources for tectonic stress fields: the product of their magnitude by their depth is shown to characterize their ability to induce deformation. This rule of the density moment directly yields the lithospheric thickening or thinning rate when applied to structures of large lateral extent. For anomalies of lateral extent that is small in comparison with their depth, the deformation is vertically inhomogeneous and has been computed with the help of simple physical models of a stratified viscous Newtonian lithosphere. The analytical treatment is based on Fourier transform. Continent-continent collision thickens not only the crust but the entire lithosphere. The cold root underlying a mountain chain induces strong regional compressive stresses able to sustain the mountain bulding process without further help from forces transmitted from far away. Thus the continental lithosphere is in a somewhat metastable mechanical state. Adiabatic, i.e. rapid, thickening (or thinning) tends to grow further once initiated. Tectonic phases of strong compression correspond to the climax of such instabilities. The response of models with cold lithospheric roots of various intensities has been computed both in two and three dimensions. They yield velocity distributions and stress fields. Instructive comparisons are made with earthquake focal mechanisms and in situ stress measurements in the Alpine and Appalachian regions. In the presence of lateral variations of the mechanical properties of the lithosphere, the tectonic style is not only function of the local topography and of the nature of its compensation. Deformations in neighbouring provinces are coupled as shown by 3-dimensional models. For example, thickening sustained by a cold lithospheric root may generate extension in peripheral zones of weakness. These last results illustrate the point that the computation of regional tectonic stresses requires the knowledge of the density anomalies within the lithosphere on the one hand, and of geometrical constraints related to lateral mechanical heterogeneities on the other.

Lithospheric readjustments after the formation of a large impact crater have been computed. The models predict tectonic perturbations over a major portion of the surface of the planet. The weight of the ejecta and the mechanical perturbations in the crater area give rise to membrane stresses. If the planetary lithosphere is in compression, the direction of maximum compression after the impact is perpendicular to a meridian (in a referential where the basin would be at the North pole). On Mercury the preferential orientation of compressive scarps radial to Caloris is explained in this manner. For a planet in extension grabens perpendicular and parallel to a meridian can appear successively as the mechanical state of the crater area is modified. On the Moon this result is consistent with the large-scale tectonics associated with Imbrium.