Benchmarking dynamic topography across geodynamical codes

During the ESA funded 4D Dynamic Earth project, different sensitivity studies are performed to understand the applicability of current ground and satellite datasets available to study the dynamical behavior of the solid Earth, in particular the complete mantle. This project is a joint effort between ESA and many European universities and is lead by Delft University of Technology (https://4ddynamicearth.tudelft.nl/). 

The project consists of ten work packages, many of them relying on some form of forward geodynamical modelling. Given the diversity of participants multiple codes are used in the project: a 2D axisymmetric Python code developed by C.T. at the Utrecht University, a 3D Matlab code developed by O.O-G. and the 3D massively parallel C++ community code ASPECT.

One recurring quantity that is of paramount importance for some work packages is dynamic topography, i.e. the outer surface expression to dynamic mantle flow. We have therefore designed a simple isothermal experiment of an anomalous sphere present in the mantle of a planet (the core is ignored as is customary in whole-Earth geodynamic modelling). The sphere itself can be positively or negatively buoyant, and the mantle can be isoviscous or characterized by a radial viscosity profile. Boundary conditions at the core-mantle boundary and at the surface are either no-slip or free-slip. 

Dynamic topography calculations involve the radial stress which is derived from the primitive variables velocity (actually, its gradient) and pressure which are found to be sensitive to mesh size in both radial and lateral directions. We therefore report on the root mean square velocity, the surface strain rate, stress and dynamic topography and the gravity anomaly for a range of experiments. Our objective is two-fold: characterize the accuracy of our codes and provide the community with a benchmark. 

All three codes are Finite Element codes and all rely on the Taylor-Hood element but they are also quite different with respect to meshing and solver architecture. Nevertheless we find that all measured quantities converge within approx. 1% for radial resolutions of at least 30km.