This summarises the findings from the 4D Deep Dynamic Earth Science Meeting on 28-29th of September 2021 (https://www.3dearth.uni-kiel.de/en/4d-deep-earth-science-meeting). To understand the dynamics of the Earth a complete picture of the whole mantle is necessary. The 3D Earth project in ESA’s Support to Science Element showed successful possibilities of a joint study across multiple disciplines toward the construction of a thermochemical model of the upper 400 km of the Earth. For this, satellite data (GOCE and Swarm) significantly improved models from seismology. A next phase is envisioned, in which whole 4D Earth models are developed to understand the link between the deep Earth and processes at its surface.

Most of Earth’s upper mantle structure has already been properly characterised by the 3D Earth project. However, both the lower mantle and upper mantle contribute to the longest wavelengths of the gravity field, which is also affected by the dynamic deformation of density discontinuities (Earth’s surface and Core Mantle Boundary-CMB). Therefore, separating the gravity contributions of from the upper and lower mantle, as commonly done, as well as explaining the topography of the Earth’s density boundaries requires a better understanding of the lower mantle structure and dynamics and its interaction of the upper mantle model. New studies using fundamental and higher mode surface waves are needed to constrain the large-scale and fine-scale structure of the lithosphere, asthenosphere, mantle transition zone and shallow lower mantle. Integrated petrological forward and inverse modelling can be an effective means of quantitatively combining the satellite data with the different types of seismic data.

The lower mantle contains huge structures or provinces featuring low seismic velocities (LLSVP), possibly related to core-mantle interactions, but also seem to play a role in global mantle flow. Can we improve the characterisation of the LLSVP and core-mantle boundary (CMB) structures, noting the placement of paleo-slabs in the lowermost mantle, by combining different datasets (gravity and various seismic constraints)? For example, normal modes (and perhaps body wave sampling of the CMB) can give us the large-scale structure of the deep lower mantle and the core. Normal mode splitting and body wave reflection techniques can provide constraints on anisotropy (important for understanding mantle flow patterns). Furthermore, can we constrain the lower mantle electrical conductivity of these structures by combining satellite magnetic and gravity data with global seismology plus mineral/high pressure physics? Large uncertainty is present about these LLSVPs, for example are the African and Pacific LLSVPs similar in terms of their composition and dynamics?

To really understand the deeper Earth, we need to study the time-dependent (4D) nature of the structures that we detect with seismology, gravity, and magnetics. Understanding the dynamics of the system will enable us to better characterise these structures and will also allow us to utilise time-dependent observations from the geological record, placing constraints on the mantle viscosity structure. 3D viscosity characterisation should be a major output of a future 4D dynamic Earth study. Currently, we cannot consistently link dynamic topography models (which include whole mantle density and viscosity, integrated in mantle flow modelling) with residual isostatic topography estimates (inferred from lithospheric structure based on satellite and terrestrial data, e.g. 3D EARTH). It is still unclear how observations of the present-day Earth are useful for understanding 4D Deep Earth processes that are manifested at Earth’s surface. For instance, subduction-driven flow is a first order feature that must be included in a whole Earth model. In addition, other dynamic solid Earth processes, like sea-level change, glacial isostatic adjustment, polar wander, and surface deformations in general are controlled by the viscosity distribution of the Earth. A new multi-disciplinary effort is needed to consistently couple a comprehensive whole Earth model with present-day dynamic surface processes.

Along these lines a new initiative was suggested as outcome of a Science Meeting held in September 2021. A future project extending 3D Earth outcomes should be split into two main phases. I) A first explorative phase focused on the assessment of the sensitivity of the different datasets: surface wave seismology, normal modes, satellite gravity, and satellite magnetic field observations for probing the solid Earth; ii) A second phase focused on generating a complete and consistent mantle model, in parallel with the geodynamic model, CMB bottom-up probing, and surface processes studies. A feedback loop is advised between the three application studies and the overall structure study. These studies could run in parallel but should have close cooperation and timely interactions like those implemented within ESA's 3D Earth project. The end goal of the initiative would be to assimilate global terrestrial and satellite geophysical and geodetical data to construct a whole-Earth model that is consistent and able to model the major dynamical processes of Earth.

In this contribution, we present a global study of the crustal structure with emphasis to cratons. In an inverse scheme, satellite gravity gradient data of the GOCE mission are inverted for the Moho depth, exploiting laterally variable density contrasts based on seismic tomography. Our results are constrained by an active source seismic data base, as well as a tectonic regionalization map, derived from seismic tomography. For the global analysis, we implement a moving window approach to perform the gravity inversion, followed by interpolating the estimated density contrasts of common tectonic units with a flood-fill algorithm.
The estimated Moho depth and density contrasts are especially interesting for the cratons of the Earth. Our results reveal a surprising variability of patterns with average Moho depth between 32-42 km, reflecting an individual tectonic history of each craton. Statistical patterns of Moho depth and density contrasts are discussed for the individual cratons and linked to their stabilization age. For example, Australia shows the lowest average Moho depth (32.7 km), indicating early stabilization in the Archean and removal of a dense lower protocrust. This observation matches well with receiver function studies. The globally inverted Moho depth is validated by gridded seismic Moho depth information, which shows that for many cratons the inverted Moho is within expected uncertainties of the seismic Moho. Especially in remote areas, which pertains to many of the African and South American Cratons, the Moho depth can now be treated with a higher confidence than previously assumed. Moreover, the formerly connected cratons of South America and Africa are analyzed and discussed in a Gondwana reconstruction. Here, the once-connected West African and Amazonian Cratons have a shallow Moho depth, indicating that only little tectonic activity occurred during the Phanerozoic. The tectonically-linked Congo and Sao Francisco Cratons have intermediate Moho depths, with the Congo Craton having a slightly shallower Moho depth. This could reflect dynamic support of the upper mantle on the African side.

Uncertainties in GIA estimates from forward models contribute significantly to the uncertainties in ice-sheet mass balance and sea level budget studies. These models rely on a past ice loading history and a lower mantle viscosity structure, both of which are, and will remain in the future, poorly constrained. In addition, models that can address lateral variations in mantle viscosity (i.e. 3-D) are in their infancy. An alternative approach is to seek a data driven solution to GIA, with the aid of vertical land movement (VLM) estimates from GPS and gravity data from GRACE. The challenges of obtaining a data-driven GIA estimate include solving a geophysical inversion that does not have a unique solution and often requires a-prior information and several approximations to achieve. For example, VLM is not only due to GIA and present day mass loading (PDML) but also to local processes including tectonics and subsidence, while GRACE measures some, generally, unknown combination of signals due to GIA and PDML.

In this work, a novel geophysical framework is developed to utilise GPS and GRACE data to solve for GIA. Our method relies on geophysical relations between geopotential and vertical land movement (VLM) caused by GIA and PDML. For example, the elastic response of the solid Earth to a positive PDML results in a negative VLM, while a positive GIA mass signal in GRACE leads to a positive VLM. We use these relations to express GPS observed VLM and GRACE observed gravity field anomalies, expressed as VLM, in terms of the sum of GIA and PDML. The method is first tested in a closed-loop synthetic experiment where we impose a uniform synthetic GPS coverage globally. It is then applied to the Nevada Geophysical Lab database of GPS VLM and GRACE spherical harmonic coefficients provided by ITG Graz. Here, we define GIA as any viscous response of the mantle to past loading, irrespective of when that took place. The resulting GIA estimate differs significantly from commonly used GIA models (such as ICE-6G) over, in particular, Alaska and central Greenland. Since the forward models do not consider viscous deformation after the end of the Little Ice Age, they will underestimate this term over Alaska and potentially in Greenland where there is evidence for a monotonic decline in volume since 1900 (Kjeldsen, Korsgaard et al. 2015) although the mantle response is expected to be slower than for Alaska. We also discuss the uncertainties, caveats and limitations of our method and its implications. The GIA product is publicly available at one degree grid resolution.



Kjeldsen, K. K., N. J. Korsgaard, A. A. Bjork, S. A. Khan, J. E. Box, S. Funder, N. K. Larsen, J. L. Bamber, W. Colgan, M. van den Broeke, M. L. Siggaard-Andersen, C. Nuth, A. Schomacker, C. S. Andresen, E. Willerslev and K. H. Kjaer (2015). "Spatial and temporal distribution of mass loss from the Greenland Ice Sheet since AD 1900." Nature 528(7582): 396.

The Earth's magnetic field displays variations on a broad range of time scales from years to hundreds of millions of years. The last two decades of global and continuous satellite geomagnetic field monitoring have considerably enriched the knowledge on the rapid physical processes taking place in the Earth’s outer core. Identification of axisymmetric torsional Alfvén waves with subdecadal periods from observatory and satellite data has given access to an averaged intensity of the magnetic field in the Earth's core interior. A significant part of the rapid signal, however, resides in non-axisymmetric motions. Their origin has remained elusive as previous studies of magnetohydrodynamic waves in the Earth's core mainly focused on their possible signature on centennial time scales. Here, we identify non-axisymmetric wavelike patterns in the equatorial region of the core surface from the observed geomagnetic variations. These wavelike features have large spatial scales, interannual periods in the vicinity of 7 years, amplitudes reaching 3 km/yr and coherent westward drift at phase speeds of about 1500 km/yr. We interpret and model these flows as the signature of Magneto-Coriolis (MC) modes.
Their identification offers a way to probe the magnetic field structure inside Earth’s core. It follows from our work that there is no need for a stratified layer at the top of the core to account for the rapid geomagnetic field changes.

Understanding the mechanisms leading to giant earthquakes and identifying when such events become imminent remains a challenge. Over the long term, seismic hazard assessment is based on the study of the recurrence intervals of earthquakes along active fault zones. Over the short term and thanks to satellite data, various deformation transients have been detected before large subduction events. Among them, we have evidenced in a previous study an anomalous gravity gradient signal during the months before the March 2011 Tohoku earthquake. We showed that it could reflect an episode of extension of the subducted slab prior to the event, generating the giant earthquake as the deformation propagated from depth to surface. These results stress the importance of the slab pull force in driving plate motions.

Here, we conduct a systematic retrospective analysis all over the globe of time series of GRACE-reconstructed gravity gradients truncated in February 2011, in order to test whether the gravity gradient variations before the March 2011 Tohoku earthquake can be identified in real-time as singular and potentially originating from the solid Earth. For that, we increase the angular resolution in our analysis of the gravity gradients and extract fast temporal signals of large amplitude, closely aligned with the orientation of the Northwestern Pacific subduction. Among all the sources of gravity variations, solid Earth signals can be recognized from a distinctive spatial pattern at the considered timescales, different from the fingerprints of water cycle sources and preferably common to different GRACE gravity solutions. For an automatic detection, we set-up a method to assess the coherency between two sets of gravity gradient signals, observed or predicted from geophysical models. It allows us to select robust signals in a chosen GRACE solution based on their coherency with the signals of another gravity solution, and test their sensitivity to modelled hydrological, atmospheric or oceanic sources. We present and discuss the results of these analyses applied to two sets of GRACE gravity models : the GRGS03 and the CSR06 solutions, as well as their respective ocean dealiasing models. Beyond the case of the Tohoku earhquake, this approach can be systematically applied to the monitoring of major plate boundaries, to identify in real-time gravity variations potentially linked with sudden changes in the deeper slab motions.

Thanks to the increasing availability of data, coming from various satellite missions, a great variety of high-resolution high-accuracy global gravity models have been developed within the last years thus paving the way to new studies ranging from the global to the local scale. It is known that gravity field data can be proficiently used to identify the shape of geological structures characterized by different densities and one of the principal applications regards the estimate of the boundaries of different crustal geological provinces, i.e. regions of the crystalline crust which can be considered homogenous from the mass density point of view (e.g. continental, oceanic or transitional crust, etc.). The classical approach for such kind of applications mainly consists in qualitative interpretation techniques that require an expert operator to supervise, guide and interpret the results. The common practice consists in fact in applying some filters to Bouguer gravity anomalies in order to remove from the data the effects of the shallowest sedimentary structures and of the mantle, respectively, thus isolating the signal due to density variations within the crust. At this point, on the reduced signal, a general visual matching with respect to a certain geological model, defined from the knowledge on the area is applied.
The principal issue of these techniques is the strong dependency of the results from the operator who performs the study, which can choose different filters that can highlight different geological features, also influencing the visual matching and final classification.
In this work we present a first application of an automatic interpretation algorithm, based on a Bayesian probabilistic approach, aimed at retrieving a map of geological crustal provinces in the Mediterranean Sea Area from gravity data. In particular, the proposed classification method requires the definition of a set of a-priori information, mainly derived from freely available geophysical data, that are used to reduce the initial observations and to prepare a rough model of the crustal provinces. Differently from the classical technique, within the proposed Bayesian classification algorithm we strip from the observed gravity the signals of main known geological horizons. Moreover, we replace gravity anomalies with second radial derivative of the gravitational potential, thus making the solution more robust with respect to deep density variations and finally, we exploit an automatic Bayesian probabilistic approach to classify the different geological crustal provinces. The result of this application is an enhanced map of the principal crustal provinces completed by an estimate of the predicted accuracy, as well as an estimate of the average density associable to each crustal province.
The proposed algorithm is applied in the Mediterranean Sea Area exploiting the global gravity field model XGM2019e, thus obtained a map of the main crustal domains. The obtained refined map presents adjustments in the boundaries with respect to the a-priori model of the order of 150-200 km, which can be an important source of information for the scientific community, allowing for a better definition of the different crustal domains. Moreover, the proposed approach, being based on freely available datasets and being a fully automatic procedure, can be in principle enlarged to study crustal provinces at continental and even global scale.