The Isidis Planitia impact basin on Mars is located on the north-south dichotomy boundary, bordered by Utopia Planitia and the Syrtis Major volcanic province. The basin records a long geological history of global and regional events of impact-induced, volcanic and sedimentary processes. This is evident in the presence of a high-density subsurface mass concentration, the strongest on Mars outside the major volcanic provinces. The nature of this interior structure remains poorly understood despite modelling efforts (e.g., [1-3]). Isidis Planitia’s surface also hosts the densest clustering of pitted cones [4,5]. The formation mechanism of these landforms, characterised by a conical mound with a central depression, remains debated as volcanic [6], sedimentary [4] or glacial [7].

We present an integrated approach to Isidis Planitia, showing that pitted cones are topographically constrained by surface wrinkle ridges driven by its subsurface structure. The subsurface is modelled using impact scaling laws combined with geological context to formulate a multi-layered model, which is fit to the local gravity field. Resultant structural elements are consistent with impact theory [8-10], estimated structures below Lunar basins [11,12], as well as mapped basins [13]. However, the gravity field cannot be constrained using infill, scaling laws and realistic density values. The models require mantle-like materials in the innermost parts of the basin. This element does not reconcile with expectations of impact theory nor basin infill, and is interpreted as significant post-impact plutonic intrusions.

This intrusive element is linked to a set of wrinkle ridge surface expressions with anomalous direction and dip. Two distinct formations of ridges are identified: an initial radial set of ridges and a latter concentric inward-dipping formation. This anomalous concentric set is not mirrored in Lunar basins [14,15] nor in Martian basins Utopia and Hellas [16,17]. The initial set is likely driven by regional compressive effects. The latter formation is driven by a stress field in the inner basin, which could be achieved during pluton inflation.

The pitted cones are shown to correlate with the basin topography dominated by the wrinkle ridges. The population conforms to both sets of pre-existing wrinkle ridges in distinct surface flow patterns. They are most consistent with volcanic rootless cones formed by lavas interacting with near-surface volatiles. The lava could be sourced from the intrusive magmatism, addressing the lack of other sources [6]. Overall, this study links Isidis Planitia’s subsurface structure to surface morphology. It could redefine the complex and dynamic basin, offering new insights into the active geological evolution of Mars.

The volcanic complex Tharsis Region on Mars is known for its numerous volcanoes on top of the crust, elevated topography (doming), and a long-wavelength gravity anomaly correlated with the region. Flexural modeling of the lithosphere has commonly been used to understand the relationship between observed topography, crustal structure, and gravity, but no conclusive answers have been obtained due to the ambiguity of these models. NASA's InSight mission has brought new information about the Martian lithosphere, which warrants a reanalysis of the support of the Tharsis Region. After analyzing the topography and gravity data, we found that a thin shell model of Mars matches both the observed gravity field for spherical harmonic degrees higher than 8 and the crustal thickness at Elysium determined by the InSight mission. Our thin shell flexure model uses an average crustal thickness of 55 km, crustal density of 3,050 kg/ (Formula presented.), average mantle density of 3,750 kg/ (Formula presented.), and an elastic thickness ((Formula presented.)) of 100 km. The mismatch between modeled and observed gravity field for the long-wavelengths (between (Formula presented.) degrees) correlates with the Tharsis Region, suggesting active large-scale dynamic support of the volcanic region. After modeling this dynamic support, we concluded that a substantial negative mass anomaly (hot buoyant mantle material, or depleted mantle region) in the mid mantle underneath the Tharsis Rise can explain the long-wavelength gravity residual. The remaining short-scale gravity residual gives insight to the Martian crustal density distribution and seems to correlate with geological structures of Mars. Buried mass anomalies in the subsurface of the northern polar plains seem not to be related to any geological or surface expressions, suggesting a more complex geology of the northern Martian crust than is suggested by the surface topography.

With the increasing interest in the Solar System's smaller bodies, quite a few missions have been sent to comets and asteroids, and more will be send in the near future. Due to the large distances involved, communication to command mission parameters takes a long time, which has a negative impact on operational safety. Autonomous navigation would be one of the key technologies that can make the mission more robust, safe, and cost effective. This is especially true if one considers the unknown fight environment when the spacecraft is first encountering the body. Most asteroids and comets have a very irregular shape and unknown mass distribution. Therefore, knowledge about its irregular gravity field will be directly beneficial as input to orbital corrections and manoeuvre planning. This paper addresses the estimation process of gravity-field parameters that could potentially be implemented in an autonomous navigation system. The focus is on a spherical-harmonic modelling of asteroid Eros-433, most notably outside the Brillouin sphere where the validity of the model is guaranteed. By using Kalman filtering it is shown that all degree and order coefficients up to degree 8 can be estimated with an error below 10%. This is the first step towards an autonomous navigation system that can operate in a highly-perturbed environment close to the asteroid.

Significant land uplift and horizontal motions have been recorded with Global Navigation Satellite Systems (GNSS) in areas such as Alaska, Iceland and the Northern Antarctic Peninsula (NAP) as a result of Glacial Isostatic Adjustment (GIA) due to ice melt after the Little Ice Age. Here, analysis of horizontal displacement rates can be of extra importance, as they are more sensitive to Earth properties in shallower layers than vertical displacement rates. Proper modelling of horizontal displacement rates with dedicated GIA models requires a spherical Earth with compressible rheology. However, in these small areas, the used GIA models are often incompressible using a Cartesian geometry to ease computation and in some cases allow for lateral viscosity changes or more complex rheology. We investigate the validity of modelled horizontal displacement rates using different approximations, that is using spherical or Cartesian Earth structures, and incompressible, material compressible or compressible rheology. Although the lack of self-gravity and sphericity compensate each other in the vertical, this is less the case for the horizontal. For a disc ice sheet with a radius just over 200 km and a thickness of 1000 m, differences due to sphericity are minimal and the modelled horizontal displacement rates of compressible Cartesian models differ from those simulated by a compressible spherical model by 0.63 mm a-1. Thus, compressible Cartesian GIA models can be applied for modelling horizontal displacement rates of small ice sheets like those in Alaska, Iceland and NAP. Unfortunately, the implementation of compressibility in Abaqus that we use here cannot be extended to spherical models as gravity can not be specified for a spherical body. Other modelling approaches are recommended in such cases.

The aim of this paper is to present the concept of a dedicated gravity field mission for the planet Mars, the Mars Quantum Gravity Mission (MaQuIs). The mission is targeted at improving the data on the gravitational field of Mars, enabling studies on planetary dynamics, seasonal changes, and subsurface water reservoirs. MaQuIs follows well known mission scenarios, currently deployed for Earth, and includes state-of-the-art quantum technologies to enhance the gained scientific signal.

Mars can provide unique insights into the mechanisms of planetary formation, thereby offering valuable clues to the early history of Earth and other rocky bodies. Currently, the internal structure of Mars is investigated using the instruments of the InSight lander, offering clues to its internal structure and formation. However, many questions remain unanswered, such as the existence and strength of convective mantle plumes. One major limitation of current experiments is that they provide measurements from only one point on the Martian surface. We propose the In-situ MArs Geodetic Instrument NEtwork (IMAGINE) instrument - a network of combined radiobeacon and laser retro-reflector instruments to be deployed on Martian surface using a swarm of wind-driven Mars Rovers. After being deployed on Mars, the instruments will be spread to cover significant portions of the Martian surface, such as the Tharsis region. They take advantage of already-existing ranging capabilities on orbital spacecraft used on legacy surface missions and can provide geodetic data over long periods up to several decades. The Tharsis region on Mars is uniquely suited to provide insights into the interior structure of Mars by investigating volcanic and tectonic activity. Gathering geodetic data and measuring potential ground deformation will offer vital clues on the mechanisms supporting the region. Moreover, it is possible to measure tidal deformations, providing more exact constraints on the Love number k₂ which can give further insight into the size of the planetary core and mantle properties. Next to that, the proposed network augments gravimetry of Mars through tracking orbiters. The radio beacon network also allows for the precise determination of Martian rotation, precession and nutation to gain insights into polar ice cap evolution and Martian interior structure. We have also identified numerous secondary applications for this network, namely long-term atmospheric studies using optical sensing. Using proven instruments and methods, it is possible to make measurements of the optical density and absorption characteristics of the atmosphere using the retro-reflector. The same instrument can also be used in fundamental science, validating aspects of general relativity. Lastly, the technical feasibility of the instrument is evaluated - while there is a laser retro-reflector with flight heritage fitting the requirements, creating a radio beacon and transmitter that is sufficiently light will require further development.

Several alternative gravity forward modelling methodologies and associated numerical codes with their own advantages and limitations are available for the solid Earth community. With upcoming state-of-the-art lithosphere density models and accurate global gravity field data sets, it is vital to understand the opportunities and limitations of the various approaches. In this paper, we discuss the four widely used techniques: global spherical harmonics (GSH), tesseroid integration (TESS), triangle integration (TRI), and hexahedral integration (HEX). A constant density shell benchmark shows that all four codes can produce similar precise gravitational potential fields. Two additional shell tests were conducted with more complicated density structures: laterally varying density structures and a crust-mantle interface density. The differences between the four codes were all below 1.5% of the modelled gravity signal suitable for reproducing satellite-acquired gravity data. TESS and GSH produced the most similar potential fields (<0.3 %). To examine the usability of the forward modelling codes for realistic geological structures, we use the global lithosphere model WINTERC-G that was constrained, among other data, by satellite gravity field data computed using a spectral forward modelling approach. This spectral code was benchmarked against the GSH, and it was confirmed that both approaches produce a similar gravity solution with negligible differences between them. In the comparison of the different WINTERC-G-based gravity solutions, again GSH and TESS performed best. Only short-wavelength noise is present between the spectral and tesseroid forward modelling approaches, likely related to the different way in which the spherical harmonic analysis of the varying boundaries of the mass layer is performed. The spherical harmonic basis functions produce small differences compared to the tesseroid elements, especially at sharp interfaces, which introduces mostly short-wavelength differences. Nevertheless, both approaches (GSH and TESS) result in accurate solutions of the potential field with reasonable computational resources. Differences below 0.5% are obtained, resulting in residuals of 0.076 mGal standard deviation at 250 km height. The biggest issue for TRI is the characteristic pattern in the residuals that is related to the grid layout. Increasing the resolution and filtering allow for the removal of most of this erroneous pattern, but at the expense of higher computational loads with respect to the other codes. The other spatial forward modelling scheme, HEX, has more difficulty in reproducing similar gravity field solutions compared to GSH and TESS. These particular approaches need to go to higher resolutions, resulting in enormous computation efforts. The hexahedron-based code performs less than optimal in the forward modelling of the gravity signature, especially with a laterally varying density interface. Care must be taken with any forward modelling software as the approximation of the geometry of the WINTERC-G model may deteriorate the gravity field solution.

Current seismic tomography models show a complex environment underneath the crust, corroborated by high-precision satellite gravity observations. Both data sets are used to independently explore the density structure of the upper mantle. However, combining these two data sets proves to be challenging. The gravity-data has an inherent insensitivity in the radial direction and seismic tomography has a heterogeneous data acquisition, resulting in smoothed tomography models with de-correlation between different models for the mid-to-small wavelength features. Therefore, this study aims to assess and quantify the effect of regularization on a seismic tomography model by exploiting the high lateral sensitivity of gravity data. Seismic tomography models, SL2013sv, SAVANI, SMEAN2 and S40RTS are compared to a gravity-based density model of the upper mantle. In order to obtain similar density solutions compared to the seismic-derived models, the gravity-based model needs to be smoothed with a Gaussian filter. Different smoothening characteristics are observed for the variety of seismic tomography models, relating to the regularization approach in the inversions. Various S40RTS models with similar seismic data but different regularization settings show that the smoothening effect is stronger with increasing regularization. The type of regularization has a dominant effect on the final tomography solution. To reduce the effect of regularization on the tomographymodels, an enhancement procedure is proposed. This enhancement should be performed within the spectral domain of the actual resolution of the seismic tomography model. The enhanced seismic tomography models show improved spatial correlationwith each other and with the gravity-based model. The variation of the density anomalies have similar peak-to-peak magnitudes and clear correlation to geological structures. The resolvement of the spectral misalignment between tomographic models and gravity-based solutions is the first step in the improvement of multidata inversion studies of the upper mantle and benefit from the advantages in both data sets.

Despite the ongoing advancements in low-thrust propulsion technology and the rise of all-electric satellite platforms, low-thrust spacecraft trajectory optimization remains a complex field of research. Shape-based approximations are predominant in interplanetary applications, but they are generally unsuitable for many-revolution trajectories, common in terrestrial applications. Indirect optimization methods allow for global optimization of many-revolution trajectories, but their mathematical complexity generally requires significant simplifications of the dynamical model, and they must be re-derived for any modification to the system dynamics or constraints. Conversely, direct optimization methods exhibit larger convergence radii and are flexible for application in different problems yet suffer from impractical computational times due to large design vectors. This paper presents a methodology for the optimization of low-thrust many-revolution trajectories, employing a hybrid combination of indirect and direct optimization methods. Similar hybrid approaches have been shown to be highly reliable for minimum-time trajectories. This methodology preserves similar performance while additionally enabling minimum-propellant optimization, through a mechanism that allows for coasting (non-thrusting) arcs, as well as targeting of the final geodetic-longitude. To reduce the propagation load of the methodology, we combine an orbital averaging scheme with a differential evolution algorithm, leading to a global optimization process with a practical computational effort. The analytical nature of the methodology reduces the number of optimization variables and its computational counterpart provides unmatchable flexibility for a configurable force and perturbation model as well as operational constraints fulfilment. The approach is applied to an unperturbed and a J₂-perturbed GTO-GEO transfer, revealing a 0.03% and a 0.4% error, for time- and propellant-minimization respectively, relative to the reference optimal trajectories. This proves that the method can match the performance of former hybrid approaches while additionally allowing for engine on/off switching. Moreover, the inclusion of the J₂ perturbation shows that, in contrast to indirect methods, it can accommodate modifications to the system dynamics without the need to re-derive the optimal control laws. Furthermore, a superior convergence radius of the optimization problem is demonstrated for the hybrid method, with respect to a reference indirect method, through the simultaneous optimization for minimum-propellant expenditure and final geodetic-longitude targeting. This research constitutes a significant advancement for space mission design and satellite operations, because it simultaneously harnesses the advantages of indirect and direct methods with broader flexibility than the popular indirect approaches and enhanced functionality than the former hybrid methods published in literature.

The Barents Sea is situated on a continental margin and was home to a large ice sheet at the Last Glacial Maximum. Studying the solid Earth response to the removal of this ice sheet (glacial isostatic adjustment; GIA) can give insight into the subsurface rheology of this region. However, because the region is currently covered by ocean, uplift measurements from the center of the former ice sheet are not available. The Gravity Recovery and Climate Experiment (GRACE) gravity data have been shown to be able to constrain GIA. Here we analyze GRACE data for the period 2003-2015 in the Barents Sea and use the data to constrain GIA models for the region. We study the effect of uncertainty in non-tidal ocean mass models that are used to correct GRACE data and find that it should be taken into account when studying solid Earth signals in oceanic areas from GRACE. We compare GRACE-derived gravity disturbance rates with GIA model predictions for different ice deglaciation chronologies of the last glacial cycle and find that best-fitting models have an upper mantle viscosity equal or higher than <span classCombining double low line"inline-formula">3×1020</span>&thinsp;<span classCombining double low line"inline-formula">Pa s</span>. Following a similar procedure for Fennoscandia we find that the preferred upper mantle viscosity there is a factor 2 larger than in the Barents Sea for a range of lithospheric thickness values. This factor is shown to be consistent with the ratio of viscosities derived for both regions from global seismic models. The viscosity difference can serve as constraint for geodynamic models of the area.

Next-generation planetary tracking methods, such as interplanetary laser ranging (ILR) and same-beam interferometry (SBI) promise an orders-of-magnitude increase in the accuracy of measurements of solar system dynamics. This requires a reconsideration of modelling strategies for the translational and rotational dynamics of natural bodies, to ensure that model errors are well below the measurement uncertainties. The influence of the gravitational interaction of the full mass distributions of celestial bodies, the so-called figure-figure effects, will need to be included for selected future missions. The mathematical formulation of this problem to arbitrary degree is often provided in an elegant and compact manner that is not trivially relatable to the formulation used in space geodesy and ephemeris generation. This complicates the robust implementation of such a model in operational software packages. We formulate the problem in a manner that is directly compatible with the implementation used in typical dynamical modelling codes: in terms of spherical harmonic coefficients and Legendre polynomials. An analytical formulation for the associated variational equations for both translational and rotational motion is derived. We apply our methodology to both Phobos and the KW4 binary asteroid system, to analyze the influence of figure-figure effects during estimation from next-generation tracking data. For the case of Phobos, omitting these effects during estimation results in relative errors of 0.42% and 0.065% for the C¯ ₂₀ and C¯ ₂₂ spherical harmonic gravity field coefficients, respectively. These values are below current uncertainties, but orders of magnitude larger than those obtained from past simulations for accurate tracking of a future Phobos lander, showing the need to apply the methodology outlined in this manuscript for selected future missions.