Europlanet Science Congress 2021
Virtual meeting
13 – 24 September 2021
Europlanet Science Congress 2021
Virtual meeting
13 September – 24 September 2021
Access live video chat


Planetary Dynamics: Shape, Gravity, Orbit, Tides, and Rotation from Observations and Models

Shape, gravity field, orbit, tidal deformation, and rotation state are fundamental geodetic parameters of any planetary object. Measurements of these parameters are prerequisites for e.g. spacecraft navigation and mapping from orbit, but also for modelling of the interior and evolution. This session welcomes contributions from all aspects of planetary geodesy, including the relevant theories, observations and models in application to planets, satellites, ring systems, asteroids, and comets.

Co-organized by OPS/SB
Convener: Alexander Stark | Co-conveners: Hannah Susorney, Anton Ermakov, Marie Yseboodt

Thu, 16 Sep, 10:40–11:25

Chairpersons: Hannah Susorney, Alexander Stark, Anton Ermakov

Ryota Nakano and Masatoshi Hirabayashi


We investigated the Yarkovsky-O’Keefe-Radzievskii-Paddack (YORP) effect [Rubincam 2000] acting on asteroid (162173) Ryugu. To accurately simulate the thermal condition on Ryugu, we employed a thermophysical model recently developed by Nakano and Hirabayashi [2021], which considers self-heating effect and scattering of radiation and solves the 3-dimensional heat equation via Finite Element Modeling (FEM) approach. Our result indicates that on 1 August 2018 TDB, the thermal torque (YORP torque) was acting in a direction that accelerates Ryugu’s spin rate.



Dynamics of small bodies among the near-Earth and main belt asteroids are controlled by secular effects due to thermal radiation forces and torques, known as the Yarkovsky and YORP effects [Rubincam 2000; Vokrouhlický and Čapek, 2002, Bottke et al., 2006]. Recent increase in the number of observed top-like shape objects suggests that the YORP effect is a key mechanism for their spin rates acceleration (or deceleration), which can induce mass-shedding or catastrophic disruptions [e.g., Holsapple, 2010; Hirabayashi, 2015]. Thus, modeling the Yarkovsky and YORP effects based on a proper thermophysical technique is one of the crucial steps to understand their formation and evolution mechanisms. Here, we investigate the YORP effect on asteroid (162173) Ryugu by using a thermophysical model recently developed by Nakano and Hirabayashi [2021].


Thermophysical modeling of Ryugu:

Unlike earlier thermophysical models which typically solve the 1-dimensional heat equation, the FEM approach thermophysical model solves the 3-dimensional heat equation. It considers the self-heating and scattering of radiation, similar to other models [e.g., Rozits and Green, 2011; Davidsson and Rickman, 2014]. The surface boundary condition includes: direct solar flux U, diffuse solar radiation flux W, direct self-heating u, and diffuse self-heating w [Davidsson and Rickman, 2014]. To efficiently describe the surface/subsurface thermal evolution, “layer” FEM shape model, where only the layers of tetrahedral meshes exist for the surface and subsurface, and the inside is left as cavity, is used. (See Nakano and Hirabayashi [2021] for more detail.)

We retrieved Ryugu’s ephemeris and the sun’s position vector as seen from Ryugu on 1 August 2018 TDB, using a SPICE kernel for Ryugu. The orbital elements and physical parameters used in this investigation are listed in Table 1. The layer FEM shape model was constructed based on SFM_49k_v20180804 [Watanabe et al., 2019], but the number of facets was decreased from 49k to 4.9k. Figure 1 shows the surface temperature distribution on Ryugu obtained after 25 spins (≈190 hr) at the fixed distance from the Sun. The surface temperature distribution is comparable to that observed by Hayabusa2 spacecraft [Okada et al., 2020].


The YORP effect on Ryugu:

Once the surface temperature is obtained, a recoil force df due to reflected solar and thermally radiated photons for each facet can be computed [Rozits and Green, 2012]. Assuming Lambert (isotropic) emission from the surface, it is given by [Vokrouhlický and Čapek, 2002; Rozits and Green, 2012]:

where ε is the emissivity of the surface, σ is the Stefan-Boltzmann constant, T is the surface temperature, G is reflected scattered flux, c is the speed of light, and dS is the facet area. The net torque τ acting on Ryugu can be computed by taking a cross product of the facet position r and the recoil force df for each facet and integrating it over the whole surface:

The net torque computed by the above equation can be decomposed into three meaningful YORP torque components of our interest: a spin period change component τω, a spin pole obliquity change component τε, and a precession of rotation τψ (see e.g., Bottke et al., 2006 for the equations). Figure 2 shows the time history of YORP torque components over one Ryugu’s spin. The one-spin averaged YORP torque components are found to be:



Due to its top-like asymmetric shape, Ryugu should be susceptible to the YORP effect. Our investigation with the 3-dimensional thermophysical modeling technique clearly showed that Ryugu does experience non-negligible YORP effect. While the YORP torque components change their signs during one spin (Figure 2), the one-spin averaged values were all found to be positive. In particular, the positive  implies an increase of Ryugu’s spin rate. This contradicts with the current hypothesis that Ryugu spun faster in the past and later slowed down to the current spin state [Watanabe et al., 2019]. In this investigation, however, we only explored Ryugu’s thermal condition at a particular epoch. Therefore, we cannot infer Ryugu’s secular dynamical evolution. The YORP effect may be highly sensitive to small topographic features on the surface [Statler, 2009]. The surface roughness may also affect thermal evolutions [Shimaki et al., 2020]. We plan to expand our modeling technique to investigate such effects and to further constrain the Yarkovsky and YORP effects on small bodies.


Table 1. List of parameters used in this investigation


Figure 1. Surface temperature distribution on 1 August 2018 TDB


Figure 2. YORP torque components versus time


Acknowledgments: RN and MH acknowledge support from NASA/Solar System Workings (NNH17ZDA001N/80NSSC19K0548) and Auburn University/Intramural Grant Program.



Bottke et al. (2006) Annu. Rev. Earth Planet. Sci., 34; Davidsson and Rickman (2014), Icarus, 243; Hirabayashi (2015) MNRAS, 454; Holsapple (2001) Icarus, 154; Nakano and Hirabayashi (2021) 52nd LPSC, 1292; Okada et al. (2020) Nature, 579; Rozits and Green (2011) MNRAS, 415; Rozits and Green (2012) MNRAS, 423; Rubincam (2000) Icarus, 148; Shimaki et al. (2020) Icarus, 348; Statler (2009) Icarus, 202; Vokrouhlický and Čapek (2002) Icarus, 159; Watanabe et al. (2019) Science, 364

How to cite: Nakano, R. and Hirabayashi, M.: Investigation of the YORP effect on asteroid (162173) Ryugu – An application of FEM approach thermophysical model, Europlanet Science Congress 2021, online, 13–24 Sep 2021, EPSC2021-441,, 2021.

Alfonso Caldiero et al.


We present here a method for the retrieval of the internal density distribution of a small body, via a least squares inversion of its gravity field and rotational state. Multiple solutions are averaged in order to limit the effect of the a priori density distribution. Simulations show successful density map estimations even with little to no initial information on the interior structure of the body.


1. Introduction

The gravity field of a body is a property among the most readily obtainable by a spacecraft mission. Be it by employing radio-tracking, optical, or gravimetric data, several missions have provided and will provide gravity expansions of degree 2 and higher for small bodies. In turn, the extended gravity field is an effect of the distribution of the mass inside the body. However, the inverse problem of estimating the density distribution from the gravity field is not a trivial one, and generally with non-unique solutions [1].

Here, we propose a method to obtain a unique solution for the internal mass distribution of a small body by restricting our estimation to zones of constant density, which are limited in number by the accuracy of the gravity model [2]. Where available, additional observables such as the libration amplitudes could help in the estimation of the interior structure. Averaging of the solutions obtained for different subdivisions of zones provides then the estimated density map [3].


2. Methods

The code builds from available software which allows to compute the gravity harmonics of a body approximated by a collection of cubes, each with a given density [4]. Here, the cubes are grouped in subregions of constant density, as many as there are measured coefficients. Each subregion is assigned an initial density, from which the software computes gravity coefficients matching in degree and order those observed. The difference between the coefficients thus computed and the observed ones is then minimized in a least-squares sense [5]. The solution of the regularized least-squares inversion is the density for each of the initial zones, along with its uncertainty. However, this solution is heavily dependent on the initial subdivision of the body into zones. Hence, a wide range of possible initial zones is explored, and the corresponding density solutions averaged together. 

From such an average solution, a new set of subregions is generated, with the help of a blob detection algorithm [6]. Finally, the gravity coefficients are inverted using this new zone subdivision. This step is necessary because, while providing a good distribution of the density zones, the density values of the averaged map are biased by solutions coming from inaccurate initial zones subdivisions. 

The data used in the inversion are the spherical harmonics coefficients of the gravity field expansion. The advantage of these observables is that they can be measured without the need for specific payload. Although the method can be easily adapted to process surface gravity measurements, the larger errors on the surface gravity coming from the cubes discretization make it less appealing for these type of data.


3. Discussion and outlook

Figure 1 shows the results obtained when applying the method to a fictitious density distribution for asteroid Bennu. The nominal density distribution, shown in Figure 1a, was used to produce a synthetic gravity field expansion of degree 11. An optimistic noise profile was assumed for these simulated observables, with a 1% error on each gravity coefficient. The spherical harmonics coefficients were then inverted via the iterative process described in Section 2, resulting in the density distribution shown in Figure 1b. 

With a gravity field of this resolution, the inversion method is clearly able to distinguish the 4 density anomalies inside the body. Still, the comparison between the nominal and the estimated density models reveals the current limitations of the method in accurately reconstructing the size of the anomalies or their density — which is equivalent, given the constraint on the total mass. 

Once a better agreement between the retrieved and the nominal models will be achieved with the current (high) resolution of the gravity harmonics, the synthetic gravity fields will be truncated at lower degrees, making them more realistic as outputs of spacecraft missions to small bodies. Future simulation studies will also make use of more realistic noise profiles, with the errors on the observed coefficient increasing as a function of the degree.



Figure 1: (a) density model used to simulate a degree-11 gravity field; (b) density distribution retrieved from said gravity field, without any initial assumptions about the true distribution



This work was financially supported by the French community of Belgium within the framework of a FRIA grant, and by the Belgian PRODEX program, managed by the European Space Agency in collaboration with the Belgian Federal Science Policy Office.



[1]  P. Tricarico. Global gravity inversion of bodies with arbitrary shape. Geophysical Journal International, 195(1):260–275, 2013.
[2]  Y. Takahashi and D. J. Scheeres. Morphology driven density distribution estimation for small bodies. Icarus, 233:179–193, 2014.
[3]  L.-I. Sorsa, M. Takala, P. Bambach, et al. Tomographic inversion of gravity gradient field for a synthetic Itokawa model. Icarus, 336:113425, 2020.
[4]  S. Le Maistre, A. Rivoldini, and P. Rosenblatt. Signature of Phobos’ interior structure in its gravity field and libration. Icarus, 321:272–290, 2019.
[5]  D. J. Scheeres, B. Khushalani, and R. A. Werner. Estimating asteroid density distributions from shape and gravity information. Planetary and Space Science, 48(10):965–971, 2000.
[6]  S. v. d. Walt, J. L. Schönberger, J. Nunez-Iglesias, et al. scikit-image: image processing in Python. PeerJ, 2:e453, 2014.

How to cite: Caldiero, A., Le Maistre, S., and Dehant, V.: Estimation of the interior density of a small body given its gravity field, Europlanet Science Congress 2021, online, 13–24 Sep 2021, EPSC2021-756,, 2021.

Ilona Hauser et al.

The binary asteroid 65803 Didymos-Dimorphos is the target of NASA’s DART (Cheng et al., 2018) and ESA's Hera missions (Michel et al., 2018). Hera will arrive at the asteroid system several years after the DART impact. It will carry out a detailed characterisation of the asteroids’ overall properties and will measure the outcome of the DART impact on Dimorphos. 

Didymos, the primary body, is a fast spinning asteroid with a period of only 2.26 hours, which is very close to the critical spin limit of 2.2 hours for asteroids bigger than approx. 200 m (Walsh et al., 2018). The interior structure of Didymos is crucial for our understanding of its stability. Without cohesion and with an estimated bulk density of 2.1 g/cc, Didymos is not able to keep its shape stable (Holsapple, 2001; Zhang et al., 2017) and thus, cohesive forces might be present in its structure.

The interior strength properties of asteroids determine to a large degree their collisional evolution. For small bodies, even a small level of cohesion can significantly affect the outcome of an impact  (Raducan and Jutzi, 2021). 

Previous studies have applied N-body and soft-sphere discrete element (SSDEM) codes to study the structural stability of rubble pile asteroids (Sanchez and Scheeres, 2012; Zhang et al., 2017, 2021; Ferrari and Tanga, 2020). Recently, the stability and failure modes of such objects have been studied using finite element (FEM) (Hirabayashi et al., 2020). Here we use Bern’s Smooth Particle Hydrodynamics (SPH) code (Jutzi et al., 2008; Jutzi 2015) to simulate the rotating asteroid and to investigate its physical properties. To do this, we set up Didymos according to its recent shape model (provided by the Hera Didymos Reference Model) with a bulk density of 2.1 g/cc and we spin up the body to a rotation period of 2.26 hours.

For the preliminary simulations presented here, we assume that Didymos has a homogeneous structure. We perform simulations using a range of values for the cohesion and the internal friction coefficient of the asteroid, and then for each case we evaluate Didymos’ stability.

An example of such a simulation is shown in Figure 1. For a cohesion of 100 Pa and a friction coefficient of 0.4, we find that Didymos’ shape is stable and only very little deformation takes place, as indicated by the small values in the total accumulated strain. In Figure 2, we show the case with a lower cohesion of 0.1 Pa and a coefficient of friction of 0.4. In this case, Didymos experiences large-scale deformations. To quantify the degree of deformation, we compute the cumulative strain distributions for both cases (Figure 3). For the 100 Pa cohesion case, > 90 % of the total mass experiences a total strain of < 0.03 and therefore this case is considered to be stable. On the other hand, for the case with 0.1 Pa cohesion, > 50 % of the total mass experiences a total strain of > 0.5, consistent with the large scale deformations observed in Figure 2. 

Our initial results suggest that a cohesion of ~ 10-100 Pa is required to maintain structural stability, which is consistent with recent SSDEM modeling results (Zhang et al. 2021). 

However, so far only homogeneous structures have been explored. We plan to investigate also more complex interior structures that may lead to different outcomes in terms of required cohesion to obtain an overall structural stability. 


Figure 1: Surface of Didymos after spin-up is complete, for a friction coefficient of 0.4 and a cohesion of 100 Pa. The body is shown from the same side at four different times. The body experiences only very little deformation, as indicated by the small amount of accumulated strain.

Figure 2: same as Figure 1, but this time for a friction coefficient of 0.4 and a cohesion of 0.1 Pa. The body gets deformed already before the spin-up is completed.

Figure 3: Cumulative strain distribution for the two simulations shown in Figures 1 and 2. The strain experienced by the body with a cohesion of 100 Pascal (a) is much smaller than for the body with a cohesion of 0.1 Pa (b).





This work has received funding from the European Union’s Horizon 2020 research and in- novation programme under grant agreement No. 870377



Cheng, A. F., et al. (2018) Planet. Space Sci., 157, 104-115.

Jutzi, M. et al. (2008) Icarus, 198:242–255

Jutzi, M. (2015) Planet. Space Sci., 107:3–9. 

Ferrari, F. and Tanga P. (2020) Icarus, 350.

Hirabayashi et al. (2020) Icarus, 352.

Holsapple (2001) Icarus, 154, 432-448.

Michel, P., et al. (2018) Adv. in Space Res., 62 (8), 2261-2272.

Raducan, S. D. and Jutzi, M. (2021) LPSC, (2548):1900.

Sanchez, D. P. and Scheeres, D. J. (2012) Icarus, 218, 876-894.

Sugiura et al. (2021) Icarus,

Walsh, K. J. (2018) Annu. Rev. Astron. Astrophys., 56, 593-624.

Zhang, Y., et al. (2017) Icarus, 294, 98-123.

Zhang, Y., et al. (2021) Icarus 362,

How to cite: Hauser, I., Jutzi, M., Raducan, S. D., and Ferrari, F.: Structural stability of 65803 Didymos: insights from SPH simulations, Europlanet Science Congress 2021, online, 13–24 Sep 2021, EPSC2021-483,, 2021.

Nair Trógolo et al.


(65803) Didymos is a binary asteroid that orbits the Sun having a semi-major axis of 1.64 AU and is the target of the DART (NASA) and Hera (ESA) missions. The system is made up of a 780 m diameter primary body (Didymos) and a 160 m satellite (Dimorphos), orbiting the primary with a semi-major axis of 1180 m and an orbital period of 11.9 h [1]. The primary has a rotation period of 2.26 h, very close to the limit of structural stability [2] [3]. The low density estimated for Didymos, 2170 kg/m3, shows that it is not a monolithic body but instead has a high macroporosity, typical of gravitational aggregates or rubble-piles. Numerical simulations show that they can be generated naturally in the aftermath of a catastrophic collision between asteroids, as a result of the gravitational interaction between the irregular fragments resulting from the collision [4]. The evolution under YORP spin-up may lead in the same case to oblate spheroidal shape asteroids with an equatorial bulge [5], commonly called top-shapes (eg: (162173) Ryugu [6], (101955) Bennu [7], (65803) Didymos, etc.). These bodies may rotate close to the limit of structural stability, kept together by the shear forces generated by friction between their components. The local acceleration near the equatorial regions may be directed outwards in these asteroids, allowing regolith to leave the surface [8] [9]. When this happens, particles evolve under the action of the gravitational field of the asteroid, the gravitational force of the Sun, the pressure of solar radiation, and in the case of binary asteroids, the secondary’s gravitational force intervene as well.


In this work, we study the dynamics of the particles that are ejected from the surface of Didymos when the centrifugal acceleration is large enough to overcome local gravity. The analysis is carried out from the development of a numerical code that integrates the particles’ equation of motion in a rotating frame of reference, centered on the primary asteroid. A polyhedral shape model for Didymos is considered, formed by 1000 vertices and 1996 triangular faces, in which centers particles are placed. Particle size distribution is generated by a power law n(r)=krα with an index set as α=-3.5. The environment of the asteroid is studied by computing the radial density of particles. To do this, a 3D grid is built; the surface is divided into bins of latitude and longitude and is propagated in radial bins. A process of detachment and re-entry of particles to the primary is observed, which could promote the potential formation of dust bands in the equatorial region of the asteroid. Trajectories are analyzed and the percentage of particles that re-accumulate on the secondary and the ones that completely escape the system are calculated.



[1] Benner, L. A. et al. Radar imaging and a physical model of binary asteroid 65803 Didymos. Bull. Am. Astron. Soc. 42, 1056 (2010).

[2] Pravec, P. and 30 colleagues 2008. Spin rate distribution of small asteroids. Icarus 197, 497–504.doi:10.1016/j.icarus.2008.05.012

[3] Margot, J.-L., Pravec, P., Taylor, P., Carry, B., Jacobson, S. 2015. Asteroid Systems: Binaries, Triples, and Pairs. Asteroids IV 355–374. doi:10.2458/azuuapress9780816532131-ch019

[4] Campo Bagatin, A., Alemañ, R.A., Benavidez, P.G., Richardson, D.C., 2018. Internal structure of asteroid gravitational aggregates. Icarus 302, 343–359.

[5] Cheng, B. and 7 colleagues 2021. Reconstruct-ing the formation history of top-shaped asteroids from thesurface boulder distribution. Nature Astronomy 5, 134–138.doi:10.1038/s41550-020-01226-7

[6] Watanabe, S. et al. Hayabusa2 arrives at the carbonaceous asteroid 162173 Ryugu—a spinning top-shaped rubble pile. Science 364, 268–272 (2019).

[7] Lauretta, D. S. et al. The unexpected surface of asteroid (101955) Bennu. Nature 568, 55–60 (2019).

[8] Campo Bagatin, A. 2013. Small Asteroids with “Dusty” Atmospheres?. Lunar and Planetary Science Conference.

[9] Yu, Y., Michel, P., Hirabayashi, M., Richardson,D. C. 2019. The expansion of debris flow shed from the primary of 65803 Didymos. Monthly Notices of the Royal Astronomical Society 484, 1057–1071. doi:10.1093/mnras/sty3515



How to cite: Trógolo, N., Moreno, F., Campo Bagatin, A., and Pérez Molina, M.: Analysis of the dynamical evolution of lofted particles around (65803) Didymos asteroid, Europlanet Science Congress 2021, online, 13–24 Sep 2021, EPSC2021-676,, 2021.

Maria del Pilar Caballo Perucha et al.

The Austrian contribution to Phase B2 Part 1 of HERA mission [1] was carried out in 2020 by JR, VRVis, and two science collaborators under GMV contract. It mainly consisted of the design of tools needed to determine Didymain´s shape [2], [3]. The challenges to be faced in Phase B2 Part 2 of the mission are focussed on Dimorphos 3D reconstruction and must be carried out within the timeframe January-June 2021, therefore some parts are still not finalised. These challenges consist of:

  • Parameter definition to optimize Dimorphos imaging to support HERA science return sequence
  • Rendering of Didymos images with PRo3D along COP (Core Operations Phase) trajectory
  • Dimorphos 3D reconstruction with ColMap [7] using the simulated images along the COP trajectory
  • Implementation of evaluation parameters in PRo3D [4] to assess Dimorphos´ reconstruction accuracy

Dimorphos images optimization

The main objective of this task is to obtain the best HERA AFC images of Dimorphos along the COP trajectory. At the same time, the maximal number of images for daily downlink must be considered as it is limited to 16 images per day. Due to the downlink data restriction the images must be optimised by selecting predetermined HERA acquisition times that provide the best imaging of Dimorphos. Therefore, images where Dimorphos appears completely occluded by Didymain, as well as images where the small asteroid is outside the AFC FoV, must be avoided for downlink and reconstruction. To recognise these images, the parameter angular deviation “α” between Didymain (Dd) and Dimorphos (Dm) was defined and the diagonal FoV of the HERA AFC was considered. The calculation and pseudocode to classify if an image should be rendered are displayed in Equation 1 and are based on HERA SPICE kernels.

Image rendering of Dimorphos

To generate synthetic images of Dimorphos in the most representative and faithful way, the information about the geometry of the Asteroid Frame Camera (AFC), its positions and orientations along the COP trajectory, and a celestial body of shape parameters similar to Dimorphos had to be considered. For the last, a 3D model of the Itokawa asteroid was textured and scaled along the three axes to get Dimorphos similar size. The resulting “*.opc” file can be seen on the left side of Figure 1. For Didymain, the same body approximation taken for Part 1 of the project was used (Figure 1,  right).