© The Authors 2023

1 Introduction

Understanding the physical properties of asteroids is required for modeling many processes, including space weathering, the formation of planetesimals, the entry of bolides into planetary atmospheres, granular mechanics, impact cratering, the thermal evolution of their parent bodies, activity drivers, and many others (e.g., Flynn et al. 2018). They are also key to properly modeling the long-term evolution of collisional asteroid families (Novaković et al. 2022). Furthermore, insight into the physical properties of asteroids is essential for the design of robotic, lander, and sample return spacecraft missions to small bodies (Murdoch et al. 2021).

Despite their great importance, little is known about the surface and internal properties of asteroids because most of them are difficult to constrain from remote observations. For instance, surface properties, such as cohesion and porosity, can be deduced by estimating thermal inertia. This, however, requires infrared observations (Alí-Lagoa et al. 2020), which are generally difficult to obtain for small asteroids. Consequently, a reliable estimation of thermal inertia is available only for a limited number of objects (see, e.g., Delbo’ et al. 2007, 2015; Harris & Drube 2016; Marciniak et al. 2019). Though this situation has begun to change (MacLennan & Emery 2021; Hung et al. 2022), new data on asteroid thermal properties and alternative methods for their determination are still of great importance.

Recently, Fenucci et al. (2021) proposed a statistical method for estimating the surface thermal conductivity of near-Earth asteroids (NEAs). The procedure is based on a comparison between the model-predicted and the measured values of the Yarkovsky effect, and, as such, it relies mainly on ground-based observations. We recall here that a generally similar idea was proposed by Rozitis & Green (2014, see also Rozitis & Green 2011), who estimated the thermal properties of asteroids with a given Yarkovsky drift using a thermophysical model (TPM). However, though the approaches share some conceptual similarities, the TPM-based method requires data such as shape models and thermal-infrared observations. Therefore, it can provide accurate estimates of the thermal properties of individual asteroids, but for a limited number of objects for which the necessary information is available. Our Monte Carlo (MC) model presented here, on the other hand, generally has limited accuracy for individual objects if only population-based parameters are used.

In addition to that, our MC model also works well in certain individual cases, such as extremely fast rotators. The key point here is that the magnitude of the Yarkovsky effect depends on a temperature gradient across the surface. For fast-rotating objects, such a gradient can be present only in the case of low surface thermal inertia, which provides an additional constraint on the model and thus allows the reliable estimation of the thermal inertia. This opportunity was already taken advantage of by Fenucci et al. (2021) to estimate the thermal properties of the small super-fast-rotating NEA (499998) 2011 PT. Here, we follow the same logic and estimate the thermal inertia of another super-fast-rotating object, namely 2016 GE1.

The low thermal inertia of a small super-fast-rotating asteroid, (499998) 2011 PT, found by Fenucci et al. (2021) was generally unexpected. Such findings may point to either a ruble-pile internal structure or the presence of a dust layer at the surface. Though generally possible, both scenarios are relatively unexpected for super-fast-rotating bodies. Asteroids with a rotation period shorter than 2.2 h are known, but they are typically small and, therefore, thought to be rocky, monolithic asteroids. This is because strength-less objects, such as rubble piles, should start disintegrating once their rotation period approaches the rotational disruption limit of 2.2 h (Pravec & Harris 2000). Still, this theory has known exceptions, suggesting that even ruble-pile asteroids are not completely strength-less. For instance, Zhang et al. (2021) studied the asteroid (65803) Didymos¹ and showed that it should have a bulk cohesion on the order of at least 10 Pa to maintain its structural stability. Regarding the presence of a dust layer on the surface of fast-spinning asteroids, theoretical works suggest that even such objects could potentially maintain small dust particles and gravel on their surface if there is relatively weak cohesion (Sánchez & Scheeres 2020). However, we do not have any direct evidence for this yet. In this respect, the extended Hayabusa 2 mission is planned to rendezvous with asteroid 1998 KY26, a small super-fast-rotating asteroid, in 2031 (Hirabayashi et al. 2021), and it will provide a better understanding of these intriguing objects.

Building upon previous works and results, here we study in detail asteroid 2016 GE1, another super-fast-rotating NEA that shares some similarities with asteroid 2011 PT, studied in Fenucci et al. (2021). Both objects are small and rotate extremely fast: 2011 PT in about 10 min and 2016 GE1 in only 34 s. Interestingly, we constrain the thermal inertia of 2016 GE1 to extremely low values with a high probability, similar to the case of 2011 PT. As mentioned above, these findings point to a subpopulation of super-fast-rotating asteroids that are either dust-covered or of high micro-porosity. However, we also discuss possible alternative explanations.

The paper is organized as follows. In Sect. 2, we report the Yarkovsky effect measurement obtained via orbit determination. In Sect. 3, we describe the method used for the estimation of the thermal properties, while in Sect. 4, a basic testing and verification of the model is presented. The results of the method applied to 2016 GE1 are reported in Sect. 5. In Sect. 6, we discuss objects similar to 2016 GE1 and the implications of our results. Finally, we summarize our conclusions in Sect. 7.

2 Preliminary consideration

2.1 Orbit determination and Yarkovsky effect detection

The orbital parameters of 2016 GE1 provided by the Jet Propulsion Laboratory Small-Body Database² (JPL SBDB) are reported in Table 1, together with their uncertainties. We also independently performed the orbit determination by using the free OrbFit software³, version 5.0.8. The Minor Planet Center (MPC)⁴ reports a total number of 127 observations for 2016 GE1, obtained on the nights of 2016 April 2 and 2019 April 5, coinciding to the close approaches with the Earth happened at a distance of 0.00356 au and 0.00743 au, respectively. All the available observations were used for the computation of the orbit. The dynamical model used for the orbital fit is similar to the one described in Del Vigna et al. (2018), and includes the gravitational forces of the Sun, the eight planets, the Moon, the 16 most massive main-belt asteroids, and Pluto. The masses and the positions of these bodies are all computed by using the JPL ephemerides DE431 (Folkner et al. 2014). In addition, we added the relativistic effects of the Sun, the planets, and the Moon expressed as a first-order post-Newtonian expansion. The Yarkovsky effect is modeled as in Farnocchia et al. (2013), as an acceleration along the transverse direction of motion aa46160-23taa46160-23 of the form $$a_t=\frac{A_2}{r^2}\widehat{t},$$(1)

where r is the distance from the Sun. The parameter A₂ is determined together with the orbital elements by fitting the model to the observations through a least-square procedure (see, e.g., Milani & Gronchi 2009). The orbit determination algorithm is also endowed with an automatic outlier rejection procedure, described in Carpino et al. (2003).

The orbital elements and the value of A₂ obtained with our orbital fit are reported in Table 1. Of the total 127 observations provided by the MPC, only 1 was rejected as an outlier. The weighted root mean square (RMS) of the astrometric residuals resulted to be 0.522 arcsec, which is only slightly smaller than the RMS of 0.567 arcsec obtained by fitting the orbit of 2016 GE1 without the Yarkovsky effect. The orbital parameters that we determined are in good agreement with those provided by the JPL SBDB, in the sense that they are the same within 2σ. It should be noted also that OrbFit provides slightly smaller uncertainties. The values of A₂ also agree within 2σ uncertainty, with OrbFit giving the smallest nominal semimajor axis drift. The value of the signal-to-noise ratio is 3.3 for the JPL solution, and 3.4 for the OrbFit solution, suggesting that the detection of the Yarkovsky effect is positive. The semimajor axis drift associated with the value of A2 obtained from astrometry is also reported in Table 1.

2.2 Preliminary constraints on 2016 GE1’s thermal inertia

In the previous subsection, we verify the positive detection of the Yarkovsky effect and the corresponding induced drift in the semimajor axis da/dt, which was found to be of considerable magnitude. This situation is somewhat unusual for a super-fast rotator such as 2016 GE1, as the Yarkovsky effect mechanism requires a temperature gradient across the surface in order to be effective. Therefore, only specific surface thermal properties may be able to produce the measured semimajor axis drift.

For this reason, we performed rough preliminary constraints on surface thermal inertia. To this purpose, we use the analytical implementation of the Yarkovsky effect in our model (see Fenucci et al. 2021). Instead of using the entire distribution of the input parameters, we defined some extreme values and provided each parameter as a single constant value. This allows the whole range of possible thermal inertia to be constrained.

In order to define the range of possible values of the surface thermal inertia of 2016 GE1, we tried different combinations of the minimum and maximum values of the parameters given in Table 2. By doing so, we considered that some parameters are correlated with changes in thermal inertia, while others are anticorrelated. For instance, larger values of the Yarkovsky drift are compatible with lower values of the thermal inertia and vice versa. As the drift is inversely proportional to the mass of the object, increasing the size or density of the asteroid, while keeping the other parameters fixed, results in lower thermal inertia.

By estimating the thermal inertia from analytical Yarkovsky formulation, we found that the thermal inertia solution is not possible for many of the combinations of the extreme values of the parameters. Nevertheless, we found a number of successful estimations. We highlight that in all of these solutions, we obtained Γ < 50 m⁻² K⁻¹ s^−1/2, which is a strong indication that the surface thermal inertia of 2016 GE1 is very low. Building on this interesting indication, we used the more complex statistical model to better constrain the thermal inertia of GE1. This model also relies on the semi-analytic implementation of the Yarkovsky effect that takes into account the orbital eccentricity (see Sect. 3), which is very important for the considered object (see Table 1).

3 Monte Carlo model and input parameters

We used the MC method developed by Fenucci et al. (2021) to estimate the thermal properties of 2016 GE1. The method is based on the comparison between the Yarkovsky drift measured from astrometry (see, e.g., Farnocchia et al. 2013; Del Vigna et al. 2018; Greenberg et al. 2020), and the model-predicted value (see, e.g., Vokrouhlický 1999; Bottke et al. 2006; Vokrouhlický et al. 2017).

In practice, the MC model by Fenucci et al. (2021) searches for input parameters so that a theoretically predicted value of the Yarkovsky effects best matches a measured value. In this respect, we recall that the Yarkovsky effect depends on several orbital and physical parameters: the orbital semimajor axis a, the orbital eccentricity e, the diameter D, the density ρ, the thermal conductivity K, the heat capacity C, the obliquity γ, the rotation period P, the absorption coefficient α, and the emissivity ε. Suppose all but one parameter are fed to the model as the inputs. In that case, the remaining parameter could be determined, provided that a measurement (da/dt)_m of the Yarkovsky effect is available.

Among all the parameters, the thermal conductivity K is the most uncertain one, because it strongly depends on the type of materials present at the surface of the asteroid, and it can vary by several orders of magnitude (Delbo’ et al. 2015). Therefore, it is the one to be determined by the model. The model versus observed Yarkovsky drift equation, $$\left(\frac{\text{d}a}{\text{d}t}\right)\left(a,e,D,\rho ,K,C,\gamma ,P,\alpha ,\epsilon \right)={\left(\frac{\text{d}a}{\text{d}t}\right)}_{\text{m}},$$(2)

is solved for K on a set of parameters randomly sampled from the input distributions, and a probability density function (PDF) is reconstructed from the output sample of K.

3.1 Semi-analytical Yarkovsky model

In Fenucci et al. (2021), the left-hand side of Eq. (2) was computed using the analytical model by Vokrouhlický (1999), which assumes a spherical shape of the asteroid, a circular orbit, and a linearization of the surface boundary condition. Despite the circular model being appropriate for (499998) 2011 PT, this may not be the case for many objects for which the Yarkovsky effect has been determined through astrometry, because NEAs generally reside at moderately to high eccentricity orbits.

For this reason, here we implemented a Yarkovsky model that takes into account the effect of the eccentricity in the orbit of the asteroid. The instantaneous osculating semimajor axis drift caused by the Yarkovsky effect is given by $${\left(\frac{\text{d}a}{\text{d}t}\right)}_{\text{i}}=\frac{2}{n^{2}a}{f}_{\text{Y}}\cdot v,$$(3)

where a is the semimajor axis of the asteroid orbit, n is the mean motion, v is the heliocentric orbital velocity, and f_Y is the instantaneous value of the Yarkovsky acceleration. The term f_Y is computed via an analytical model described in Vokrouhlický et al. (2017), which assumes a spherical shape of the asteroid and a linearization of the surface boundary condition, and it is given by $$f_{\text{Y}}=f_{\text{Y},}{}_{\text{d}}+f_{\text{Y},\text{\hspace{0.17em}s}},$$(4)

where f_{Y, d} and f_{Y, s} are the diurnal and the seasonal component, respectively. The diurnal component is expressed as $$f_{\text{Y},\text{\hspace{0.17em}d}}=\kappa \left[\left(n\cdot s\right)s+{\gamma }_1\left(n\times s\right)+{\gamma }_{2}s\times \left(n\times s\right)\right].$$(5)

In Eq. (5), n = r/r is the heliocentric unit position vector, and s is the unit vector of the asteroid spin axis. In addition, $$\kappa =\frac{4\alpha }{9}\frac{SF}{mc},$$(6)

where S = 4πR² is the cross section of the asteroid, R is the radius, F is the solar radiation flux at a heliocentric distance r, m is the asteroid mass, and c is the speed of light. The coefficients γ₁ and γ₂ are expressed as $$\begin{array}{l}{\gamma }_1=-\frac{k_1\text{Θ}}{1+2k_2\text{Θ}+k_3{\text{Θ}}^2},\hfill \\ {\gamma }_2=-\frac{1+k_2\text{Θ}}{1+2k_2\text{Θ}+k_3{\text{Θ}}^2},\hfill \end{array}$$(7)

where $\text{Θ}=\sqrt{\rho KC{\omega }_{\text{rot}}}/\left(\epsilon \sigma T_{\star }^3\right)$ is the thermal parameter, ω_rot is the rotation frequency, σ is the Stefan-Boltzmann constant, and T_⋆ is the subsolar temperature, which is given by $4\epsilon \sigma T_{\star }^4=\alpha F$. The coefficients k₁, k₂, and k₃ are positive analytic functions of the rescaled radius $R_{\text{d}}^{\prime }=R/l_{\text{d}}$, where $l_{\text{d}}=\sqrt{K/\left(\rho C{\omega }_{\text{rot}}\right)}$ is the penetration depth of the diurnal thermal waves. The analytic expressions of the k_i, i = 1, 2, 3 coefficients can be found in Vokrouhlický (1998, 1999).

The seasonal component is given by $$f_{\text{Y,s}}=\kappa \left[{\overline{\gamma }}_1\left(n\cdot s\right)+{\overline{\gamma }}_2\left(N\times n\right)\cdot s\right]s,$$(8)

where N is the unit vector normal to the orbital plane, and ${\overline{\gamma }}_1,{\overline{\gamma }}_2$ have the same expressions as Eq. (7), but evaluated with thermal parameter $\overline{\text{Θ}}=\sqrt{\rho KCn}/\left(\epsilon \sigma T_{\star }^3\right)$, and with rescaled radius $R_{\text{s}}^{\prime }=R/l_{\text{s}}$ where $l_{\text{s}}=\sqrt{K/\left(\rho Cn\right)}$ is the penetration depth of the seasonal thermal waves.

The average Yarkovsky drift, da/dt, is then obtained by averaging the instantaneous Yarkovsky drift of Eq. (3) over an orbital period: $$\frac{\text{d}a}{\text{d}t}=\frac{1}{2\pi }{\displaystyle {\int }_0^{2\pi }{\left(\frac{\text{d}a}{\text{d}t}\right)}_{\text{i}}\text{\hspace{0.17em}d}\ell },$$(9)

where ℓ is the mean anomaly. The integral at the right-hand side of Eq. (9) is numerically computed with the trapezoid rule (see, e.g., Stoer & Bulirsch 2002). To this purpose, we used a fixed step in the eccentric anomaly u, which is then translated into a step in ℓ. This is done to secure a proper sampling of the orbit around the perihelion.

3.2 Defining the input parameters

The next step toward modeling the Yarkovsky effect is to generate the most likely probability distribution of the input parameters on which the effect depends. As the availability of these parameters and their uncertainties could be very different, we divided the input parameters into three categories.

The first group includes the parameters known with high accuracy. In such cases, instead of providing the probability distribution of a parameter to the model, we used only a single value, which is the nominal value of the parameter. The second group contains parameters that are determined for the individual asteroids, but their uncertainties are not negligible. Such parameters are modeled assuming a Gaussian distribution with a standard deviation corresponding to the estimated errors of the parameters. Finally, the third group includes the parameters that are, in most cases, not available for individual objects. In these cases, their probability distribution is derived from a population-based distribution.

3.2.1 Fixed parameters

Changes within 3σ in the semimajor axis a and eccentricity e produce negligible fluctuations in the left-hand side of Eq. (2), and therefore they are kept at their nominal values. On the other hand, the heat capacity C, the emissivity ε, and the absorption coefficient α are all unknown. However, the plausible range for each of these parameters is narrow compared to the uncertainty of other relevant quantities. Therefore, instead of providing a full distribution of these parameters, we fixed their values in each simulation and provided results with a few different values. Typical heat capacity values assumed for asteroids are in the range 600−1200Jkg⁻¹ K⁻¹ (Farinella et al. 1998; Delbo’ et al. 2015; Piqueux et al. 2021). Therefore, we performed simulations for a few values from this interval.

For the emissivity, ε, we adopted 0.984 as a nominal value corresponding to the mean value of measurements performed on meteorites (Ostrowski & Bryson 2019). Additionally, we tested how much the results change if a value of ε = 0.9 is assumed.

The absorption coefficient is defined as α = 1 − A, where A is the Bond albedo of the object. As we found that in the case of 2016 GE1, it does not affect the main conclusions (see Sect. 5.1), the absorption coefficient α was set to 1 in all our simulations. A test has been performed to verify that assuming a value of 0.9 does not change the results significantly (see Sect. 5.1).

3.2.2 Modeling 2016 GE1-based parameters

These parameters are modeled according to their values determined for the asteroid 2016 GE1. All these are assumed to be Gaussian distributed, with a mean value equal to the nominal estimated value and standard deviation equal to the 1σ uncertainty of the measurement. The value of the measured semimajor axis drift (da/dt)_m has been estimated through orbit determination, and it is reported in Table 1. The absolute magnitude H, although not explicitly present in Eq. (3), is needed for the model by Fenucci et al. (2021) in order to construct a population-based distribution of the diameter D and of the density ρ (see Sect. 3.2.3). The measured absolute magnitude of 2016 GE1 is H = 26.7, and we used the uncertainty of 0.5 given by OrbFit after the convergence of orbital fit.

The rotation period reported in the Asteroid Light Curve Database (LCDB; Warner et al. 2009) is P = 0.009438 h, which corresponds to about 34 s. The light curve was obtained with an exposure time of 10 s (Warner 2016), and the quality code U reported in LCDB corresponds to 2, implying an uncertainty of about 30%. Therefore, we used a σ = 0.3 × P for this parameter. In fact, the U = 2 code flag might also imply that the period could be wrong by an integer multiple. However, the short period mentioned above was recently confirmed also by Ghosal et al. (2022), which gives us some confidence that it is accurate. Therefore, despite this limitation, 2016 GE1 is very likely the NEA with the shortest rotation period for which the Yarkovsky effect has been determined so far. Nevertheless, in Sect. 5 we analyzed how the estimated thermal inertia changes with rotation period, and the limitations this might impose on the results.

3.2.3 Modeling population-based parameters

The density ρ and the diameter D of 2016 GE1 are both unknown. Therefore, we use the population-based distribution model by Fenucci et al. (2021). The model combines the NEA orbital distribution by Granvik et al. (2018) and the NEA albedo distribution by Morbidelli et al. (2020), and provides a bi-variate distribution of the couple (ρ, D). Given the orbital elements a, e, i and the absolute magnitude, H, of an object, the population model by Granvik et al. (2018) provides the probability of an NEA originating from each main-belt source region, which are: the ν₆ secular resonance region, the 3:1, 5:2, and 2:1 mean-motion resonances regions with Jupiter, the Hungaria region, the Phocaea region, and the Jupiter-family comets (JFCs) region. We note that the absolute magnitude value of 2016 GE1 is beyond the limit of H = 25 of validity of the model by Granvik et al. (2018); therefore, the source-route probabilities (reported in Table 3) are extracted by a linear interpolation. These probabilities are then combined with the NEA albedo distribution by Morbidelli et al. (2020), and a PDF p_pv for the albedo p_V is determined first (see Fenucci et al. 2021, for details).

To obtain a distribution of (ρ, D), we sampled the albedo according to its PDF p_pv. For each point of the sample, we produced a value of diameter D by using the conversion formula (see, e.g., Bowell et al. 1989; Pravec & Harris 2007) $$D=\frac{1329\text{\hspace{0.17em}km}}{\sqrt{p\text{v}}}{10}^{-H/5}.$$(10)

The same albedo value is used to generate a value of the density ρ. To this end, we divided the albedo into three categories and associate an asteroid complex with each of them: p_V ≤ 0.1 is associated with the C-complex, 0.1 < p_V ≤ 0.3 with the S-complex, and p_V > 0.3 with the X-complex.⁵ A value of the density ρ is then generated according to the class in which the selected albedo value falls in. The density of each group is assumed to be log-normal distributed, with the average and the standard deviation listed in Table 4.

Figure 1 shows the joined distribution of (ρ, D) obtained for 2016 GE1, where the correlation given by the albedo can be seen by the fact that smaller size is associated with a larger density (moderate and large albedo cases), while the larger size is associated with a smaller density (low albedo case). The marginal PDFs of ρ and D (i.e., the distribution of the set of all the possible ρ and D values alone) are also shown in Fig. 1. The most likely value and the median of the density are almost the same, and they are about 2490 kg m⁻³. The most likely value of the diameter is 12 m, while the median value is 14m.

The obliquity γ is also unknown. Therefore, we assume it to be distributed according to the NEA obliquity distribution determined by Tardioli et al. (2017). This distribution has a 2:1 ratio between retrograde and prograde rotators. However, the model always rejects values of γ that are not compatible with the sign of the measured Yarkovsky drift because, in this case, solutions to Eq. (2) cannot be found.

Fig. 1 Input density, ρ, versus the diameter, D, distribution for 2016 GE1. The blue histograms at the top and right show the marginal distributions of ρ and D, respectively.

4 Basic model testing and verification

An essential step in the application of any new method is its testing and verification. To this purpose, we used asteroid (101955) Bennu, the target of NASA’s OSIRIS-REx sample return mission, for which all relevant input and output parameters are well constrained, allowing us to test how different assumptions or unknown parameters affect the results.

The values of orbital and physical parameters of Bennu used in our model are given in Table 5. We note that there are different estimations of Bennu’s thermal inertia. Based on the Spitzer Space Telescope measurements of the Bennu’s thermal emission, Emery et al. (2014) derived thermal inertia of 310 ± 70 Jm⁻² K⁻¹ s^−1/2. A global-average thermal inertia of 350 ± 20 Jm⁻² K⁻¹ s^−1/2 was estimated by Dellagiustina et al. (2019), using the OSIRIS-REx approach-phase thermal emission light curves and the encounter-based shape model. More recently, Rozitis et al. (2020) analyzed the data from the OSIRIS-REx Thermal Emission Spectrometer (OTES) and the OSIRIS-REx Visible and InfraRed Spectrometer (OVIRS), and derived mean thermal inertia values of 300 ± 30 and 320 ± 30Jm⁻² K⁻¹ s^−1/2 from OTES and OVIRS, respectively. The authors also found spatial variations in thermal inertia, with larger values closer to the equator.

In principle, Bennu’s thermal inertia values from the literature are very similar and generally consistent. We primarily compare our test results with values from Rozitis et al. (2020), but conclusions remain broadly the same if the other values are used for comparison.

To test our model, we adopted the following strategy. Since for a typical asteroid we have much less data than for Bennu, we started with just basic properties. Then, in each subsequent step, we added some new parameters with better constraints. In this way, we expect to show that the model provides meaningful results even with only basic knowledge about the object, but that the results become more accurate when additional knowledge is available. Therefore, we performed four tests. In the first test, we assumed that only basic information about the object are available: the orbit (semimajor axis and eccentricity), absolute magnitude and rotation period. In the second test, albedo and diameter are included as well. The third test also includes Bennu’s density, while in the final fourth test, we added the obliquity. The information about the test are summarized in Table 6, while the results are also shown in Fig. 2.

Before we discuss the results obtained for Bennu, we recall that the thermal inertia solution is typically bi-modal, resulting in two prominent peaks in the obtained distribution. While both solutions are possible in principle, based on the empirical understanding, we favor the right peaks. Therefore, in what follows, we refer to the right peaks in the distribution of thermal inertia as our nominal results. In some cases, depending on the exact goal, we may speak about maximum (or minimum) values of the thermal inertia, limiting the results from one side rather than providing the most likely value.

The result of the first test (Γ = 592 ± 257 Jm⁻² K⁻¹ s^−1/2), which is significantly based on NEA population models, is skewed toward higher values and the corresponding uncertainty is large. This shows that without any prior knowledge of the physical parameters, the results obtained for a single asteroid are unreliable in a general case. However, despite the limitations of the large uncertainty, the result is still statistically compatible with the values found in the literature, except the one derived by Rozitis et al. (2020) from the OTES instrument (see Tables 5 and 6). Adding knowledge about the albedo and diameter in the second test resulted in a somewhat improved result and reduced uncertainty (Γ = 505 ± 189 Jm⁻² K⁻¹ s^−1/2). Including knowledge of density in the third test changed the situation significantly. The obtained thermal inertia of Γ = 346 ± 64 Jm⁻² K⁻¹ s^−1/2 is now fully in line with the high-accuracy measurments. Interestingly, when the information about the obliquity is added in the fourth test, the result is moved further from the referent value, though still plausible.

Based on the tests presented above, we conclude that our MC-based model for thermal inertia determination could be useful even when only basic information about an object are known. In this case, however, the uncertainty of the result could be large, and the model may not be fully appropriate for individual objects. Noticeable exceptions from this are rapidly rotating objects (see Sect. 5.1).

On the other hand, if most of the input parameters are sufficiently well-known, the model provides accurate and reliable results even in the case of an individual object.

Fig. 2 MC model tests on asteroid Bennu. The panels show the distributions of thermal inertia solutions for different input parameters assumed to be known. The number of known parameters increases clockwise, starting from the bottom-left panel. For additional details on the test, see Table 6. The gray area in each panel marks the interval of the thermal inertia of Bennu as estimated by Rozitis et al. (2020).

5 Estimated thermal characteristics of 2016 GE1

We performed the MC estimation of the thermal conductivity K for four fixed heat capacity values, namely C = 600,800, 1000, and 1200 Jkg⁻¹ K⁻¹. A random sample of one million combinations of the input parameters was used. Thermal conductivity solutions were searched for in the range between 10⁻⁸ Wm⁻¹ K⁻¹ and 500 Wm⁻¹ K⁻¹, which we believe to be more than appropriate considering the known variety of materials composing asteroids. For each K solution, we also computed the corresponding thermal inertia as $$\text{Γ=}\sqrt{\rho KC}.$$(11)

Moreover, we ran the estimation of K with the semimajor axis drift solution given by the JPL SBDB first, and then ran it once again with the solution we obtained with OrbFit.

Figure 3 shows the distributions of K and Γ obtained for the different values of heat capacity C. Blue histograms refer to the results obtained with the JPL SBDB orbital solution, while red histograms to those obtained with our orbital fit.

All distributions have the same properties. For a thermal conductivity smaller than ~10⁻³ Wm⁻¹ K⁻¹ (thermal inertia smaller than ~100Jm⁻² K⁻¹ s^−1/2), two peaks always occur with a high probability. Values for K between ~0.01 and ~1 Wm⁻¹ K⁻¹ and for Γ between ~100 and ~1000 Jm⁻² K⁻¹ s^−1/2 are extremely unlikely. Finally, for values K ≳ 1 Wm⁻¹ K⁻¹ (Γ ≳ 1000 Jm⁻² K⁻¹ s^−1/2), there is a low probability tail in which two more peaks occur.

The distributions are almost independent of the heat capacity C. Therefore, we always refer to the results for C = 600 J kg⁻¹ K⁻¹ in the following, unless explicitly stated otherwise. We performed the Kolmogorov–Smirnov test to check

whether the distributions obtained for the two solutions for the semimajor axis drift are identical. The null hypothesis was rejected with a significance level of 5%. Nevertheless, the two distributions are very similar in each of the cases shown in Fig. 3, and the results show only minor differences.

The two high probability peaks at low thermal conductivity appear at $$\sim 3-6\times {10}^{-6}\text{W\hspace{0.17em}}{\text{m}}^{-1}{\text{K}}^{-1}$$

and at $$\sim 1-2\times {10}^{-4}\text{W\hspace{0.17em}}{\text{m}}^{-1}{\text{K}}^{-1},$$

while the corresponding peaks in thermal inertia are at $$3\text{\hspace{0.17em}J\hspace{0.17em}}{\text{m}}^{-2}{\text{K}}^{-1}{\text{s}}^{-1/2}$$

and at $$18\text{\hspace{0.17em}J\hspace{0.17em}}{\text{m}}^{-2}{\text{K}}^{-1}{\text{s}}^{-1/2}.$$

We note that these peaks are located at extremely low K, and they are almost an order of magnitude smaller than those obtained for 2011 PT (Fenucci et al. 2021, see also Fig. 9). On the other hand, the two low probability peaks at high thermal conductivity are at ~10 and ~250 Wm⁻¹ K⁻¹, corresponding to the thermal inertia of ~4 450 and ~18 000Jm⁻² K⁻¹ s^−1/2, respectively. The presence of these two low probability peaks is due to the fact that the measured versus predicted Yarkovsky drift equation in Eq. (2) has either three or four different thermal conductivity solutions for certain combinations of input parameters. To show an example of this behavior, we computed the Yarkovsky drift on a grid in density ρ and thermal conductivity K, by fixing the other parameters to D = 9.8 m, C = 800 Jkg⁻¹ K⁻¹, γ = 133.77°, and P = 0.009438 h. Figure 4 shows the contour plot of the computed semimajor axis drift. The levels corresponding to the da/dt solution computed with OrbFit, together with those at 1σ uncertainty, are highlighted in red. From this figure, it can be appreciated that the measured versus predicted Yarkovsky drift equation has either three or four solutions in the interval K ∈ [10⁻⁸, 500] W m⁻¹ K⁻¹, for certain values of the density, ρ.

To give quantitative constraints of the thermal parameters, we fit the thermal inertia distributions by using the kernel density estimation, and then we computed the probability $$\begin{array}{c}{P}_1=P\left(\Gamma <100\text{\hspace{0.17em}\hspace{0.17em}J\hspace{0.17em}}{\text{m}}^{-2}{\text{K}}^{-1}{\text{s}}^{1/2}\right),\\ P_2=P\left(\Gamma >1000\text{\hspace{0.17em}\hspace{0.17em}J\hspace{0.17em}}{\text{m}}^{-2}{\text{K}}^{-1}{\text{s}}^{1/2}\right).\end{array}$$(12)

The probability P₁ was always ~0.92, with negligible differences between the two solutions of the semimajor axis drift. On the other hand, the probability P₂ was always ~0.08, implying that solutions with high thermal inertia are unlikely. Therefore, the results show that the thermal inertia of 2016 GE1 is probably very low, which is unexpected for the extremely fast-rotating asteroid.

Fig. 3 Distributions of the thermal parameters of 2016 GE1 for different values of heat capacity, C. The first column shows the distributions of the thermal conductivity, K, and the second column the distributions of the thermal inertia, Γ. Blue histograms are the results obtained using the orbital solution provided by the JPL SBDB, while red histograms are the results obtained using our solution obtained with OrbFit.

Fig. 4 Estimated semimajor axis drift for 2016 GE1, obtained for D = 9.8m, C = 800Jkg−1 K−1, γ = 133.77°, and P = 0.009438h. The red level curves represent the Yarkovsky drift measured from astrometry with OrbFit and the corresponding 1σ uncertainty.

5.1 Robustness of the results

The results presented suggest that 2016 GE1 has very low thermal inertia. In this subsection we discuss how reliable such a conclusion is. The presented MC-based model for asteroid thermal inertia estimations generally depends on a set of input parameters. As demonstrated in Sect. 4, when these parameters are well known (or at least most of them), the model could provide good results for individual objects.

On the other hand, when only basic information such as orbit, absolute magnitude, and Yarkovsky drift are known, the model relies on population-based models of input parameters. In this case, the model could still provide useful results to model the thermal inertia of a population of asteroids. However, the results for an individual object are uncertain and generally unreliable. The situation with 2016 GE1 is very similar to that case, except that the rotation period is known, in addition to the orbit and absolute magnitude. This leads us to that question of why then the obtained results should be considered reliable.

The asteroid 2016 GE1 is a rapid rotator with a rotation period of only about 34 s. At the same time, a significant A₂ acceleration associated with the Yarkovsky effect has been measured. A temperature gradient across the surface must be present for this effect to work. However, in the case of such a rapid rotation, the temperature gradient could exist only in the case of low thermal inertia. As a result, the rotation period strongly constrains the range of acceptable thermal inertia. In the case of 2016 GE1, only about 16% of input parameter combinations are accepted in our simulations as possible. This is largely due to solid constraints from the rotation period.

As explained in Sect. 3.2.2, despite some uncertainties, we believe the short rotation period of 2016 GE1 is determined reliably enough. Nevertheless, we tested what would happen if the rotation period is 2, 5 and 10 times longer than the measured one. In Fig. 5 we show how the resulting thermal inertia of 2016 GE1 depends on the rotation period. Increasing the rotation period shifts the distribution of thermal inertial to the right, that is to say, toward larger values. Still, even for a 10× longer rotation period, more than 90% of the solution suggests thermal inertia below 300 Jm⁻² K⁻¹ s^−1/2, which can be considered low.

There are indeed small peaks in the TI distribution associated with values of Γ > \000 Jm⁻² K⁻¹ s^−1/2. These values are incredibly high to the point of physical implausibility, and other studies have found very few objects with such high thermal inertia estimates (Hung et al. 2022). We, therefore, discarded them as highly improbable.

We also tested how the assumed values of the emissivity ε and the absorption coefficient α may affect the results. We found that for a reasonable range of these parameters, the changes in results are small and cannot affect general conclusions about the low thermal inertia of 2016 GE1 (see Fig. 6). Additionally, we investigated how the assumption of entirely random input distribution of the density ρ in the range 1000–3500 kg m⁻³ would change the result. Again, we found that the resulting thermal inertia values are not much different. Moreover, they are even shifted toward lower values (bottom panel of Fig. 6).

To conclude this part, we show how the results change when the eccentric and circular Yarkovsky models are used. As already mentioned in Sect. 3.1, the Yarkovsky model based on the assumption of a circular orbit may not be suitable due to the large eccentricity of the orbit. In Fig. 7, we show how the results differ for two Yarkovsky models when all other parameters are equal. The thermal inertia solutions in the case of the circular model favor even smaller values. Two peaks at low Γ are slightly shifted to lower values, while two smaller peaks found in the eccentric model at high values of thermal inertia disappear in the circular model. This demonstrates that the circular model is not fully suitable for NEAs that are in moderately to highly eccentric orbits. However, in the case of 2016 GE1, the results seem to be mainly driven by its rapid rotation. They, therefore, do not vary much for different input parameters or the Yarkovsky model.

With this in mind, we concluded that our estimate of the 2016 GE1 low thermal inertia is robust, even though we cannot very accurately estimate its nominal value. Nevertheless, we showed that thermal inertia of GE1 cannot exceed Γ = 300 Jm⁻² K⁻¹ s^−1/2. The real thermal inertia is likely even smaller, with a probability of > 90% to be below 100 J m⁻² K⁻¹ s^−1/2.

5.2 Posterior distribution of the input parameters

In addition to calculating the thermal inertia, we kept track of all combinations of input parameters for which at least a single solution of Eq. (2) was found. If no solutions are found, the measured semimajor axis drift cannot be achieved for the chosen combination, and therefore these determined values of the physical parameters are not representative for 2016 GE1.

Figure 8 shows the 2D distribution of (ρ, D) and their marginal distributions obtained for the simulations using the semimajor axis solution provided by the JPL SBDB and C = 600 J kg⁻¹ K⁻¹. The results in all other cases are similar. The correlations low density-large diameter and high density-small diameter are still present in the output distributions, as shown in the top panel of Fig. 8. The median value and the most likely diameter are 9.8 and 9.2 meters, respectively, which are both smaller than the values obtained for the input distribution in Fig. 1. This means that the solutions of Eq. (2) are more likely to be found for small diameters. From this distribution, we determined that 2016 GE1 has a diameter ranging from 5 m to 15 m with a probability of ~0.87 and is smaller than 20 m with a probability of ~0.94.

The density distribution peaks at about 2020 kg m⁻³, with a median value of ~2100kgm⁻³, both lower than the corresponding values of the input distribution. The probability of 2016 GE1 having a density of less than 1200 kg m⁻³ or more than 3000kgm⁻³ is small, ~0.06 in each case. We note, however, that the density distribution is two-peaked, as expected, with each peak generally corresponding to one of two main taxonomic complexes of C- and S-type asteroids. Therefore, if asteroid 2016 GE1 is found to belong to the C-type, the larger values of the density should be rejected. Additionally, if 2016 GE1 is an X-complex low-albedo but a high-density object, the assumed correlation between diameter and density would be broken, and the analysis of the posterior parameter distribution would be meaningless. We consider, however, this scenario unlikely.

Fig. 5 Dependence of 2016 GE1’s thermal inertia estimation on the rotation period. The assumed period increases clockwise from the bottom-left panel. The results are shown for the nominal period solution of 34 s, as well as for 2, 5, and 10 times longer periods, as indicated in the plots.

6 Discussion

6.1 Possible explanations for low thermal inertia

The values of Γ < 100 J m⁻² K⁻¹ s^−1/2 are obtained with a probability of ~92 percent. Moreover, we have tested the stability of the result with respect to various input parameters and verified their reliability. Therefore, the results presented above strongly suggest a very low thermal inertia of the rapidly rotating asteroid 2016 GE1.

These thermal inertia values are not new for asteroids (see, e.g., Delbo’ et al. 2015). For example, the largest asteroids, such as (1) Ceres or (4) Vesta, have very low thermal inertia values of Γ < 50 Jm⁻² K⁻¹ s^−1/2, as measured by Leyrat et al. (2012) and Rognini et al. (2020), and this is due to the fine regolith present at their surfaces, which in turn is made possible by the relatively high gravity of these objects.

Another cause for low thermal inertia has been recently discovered during the OSIRIS-REx and Hayabusa 2 missions, and it lies in the high micro-porosity of boulders. The ground-based estimated thermal inertia of Γ = 3\0 ± 70Jm⁻² K⁻¹ s^−1/2 of Bennu (Emery et al. 2014) suggested a fine regolith-covered surface, which was not found when OSIRIS-REx arrived at the asteroid (Lauretta et al. 2019). Rozitis et al. (2020) explained the measurements with the high micro-porosity of Bennu’s surface boulders, and suggested that this could be a characteristic property of C-type asteroids. Furthermore, Cambioni et al. (2021) found that the thermal inertia of Bennu’s rocks is positively correlated with the local surface abundance of fine regolith. Similar results were obtained for Ryugu, where the global thermal inertia was estimated to be Γ = 225 ± 45 J m⁻² K⁻¹ s^−1/2 (Shimaki et al. 2020), but no fine regolith-covered surface was found (Watanabe et al. 2019). In situ images taken by the Mobile Asteroid Surface Scout (MASCOT; Ho et al. 2017) lander did not detect fine regolith on the boulders. Models of the temperature variations did not fit the measurements when the regolith layer was taken into account, and the low thermal inertia was again explained by the high porosity of the boulders (Grott et al. 2019).

The spatially unresolved thermal inertia of comets is also very low, typically below 50 Jm⁻² K⁻¹ s^−1/2 (e.g., Marshall et al. 2018; Groussin et al. 2019), although spatially resolved data indicate surface variations. In the case of comet 67P/Churyumov-Gerasimenko, for example, most values are in the range 10–170 J m⁻² K⁻¹ s^−1/2. The thermal inertia of smooth terrains covered with deposits is lower (typically lower than 30 Jm⁻² K⁻¹ s^−1/2) than those of exposed consolidated terrains (typically larger than 110 Jm⁻² K⁻¹ s^−1/2) (Groussin et al. 2019).

This then raises the questions of what can explain 2016 GE1 and what the most likely scenario is. Low thermal inertia is unexpected for super-fast-rotating asteroids. On the one hand, this is because the fast rotation should clear out the surface of loose regolith. On the other hand, a large porosity often points toward a rubble-pile structure, which is also unexpected because the internal strength of such bodies might be too weak to maintain the integrity of the body under high rotational acceleration. Despite their improbability, however, none of the scenarios can be discarded, as recent models suggest that they are still possible (Sánchez & Scheeres 2020; Hu et al. 2021). The integrity of the fast rotators could be maintained by a low yield stress on the order of 25 Pa. This level of yield stress can only be explained by assuming that the particles in asteroids are small enough (on the order of a few micrometers) to form a cement matrix (glue) between the larger particles (fragments) (Persson & Biele 2021).

Still, the most likely values for the thermal inertia found for 2016 GE1 by using the semi-analytical Yarkovsky model are significantly below 100 Jm⁻² K⁻¹ s^−1/2, smaller even than the thermal inertia of the boulders on Bennu and Ryugu. For the single boulder on Ryugu investigated by the MASCOT lander, Hamm et al. (2022) estimated the thermal inertia of ${256}_{-3}^{+4}\text{J\hspace{0.17em}}{\text{m}}^{-2}{\text{K}}^{-1}{\text{s}}^{-1/2}$, corresponding to an expected porosity of ${46.7}_{-0.4}^{+0.3}\%$. This means that, if 2016 GE1 was a boulder, than it would be more porous than those on Bennu or Ryugu. To get a general idea of what porosity would explain the obtained values of thermal conductivity and inertia, we applied an empirical relationship between conductivity and porosity given by Grott et al. (2019) (see also Henke et al. 2016) for meteorites. The employed empirical relation has its own limitations and, in particular, it is unreliable for large porosity. Nevertheless, it can give a first indication of the degree of porosity in 2016 GE1 needed to explain its thermal properties. We found that values of thermal inertia Γ ≤ 20 J m⁻² K⁻¹ s^−1/2, compatible with two the most prominent peaks (see, e.g., Fig. 3), are only feasible for porosity ≥ 70% and density ρ ≤ 1000 kgm⁻³. Therefore, a highly porous, Ryugu-boulder-like object is generally consistent with our result. If so, it also means that most of the density solutions should be discarded, and only those generally compatible with C-type asteroids can be considered. This opens the possibility of small and super-fast rotators having a common origin as anomalously low-thermal-inertia boulders similar to those found on Ryugu, which may be worth exploring in the future since they are the most similar to the primordial planetesimals (Sakatani et al. 2021).

On the other hand, Sánchez & Scheeres (2020) developed a model to study under what rotational conditions an asteroid can keep thin regolith on the surface, assuming that the asteroid has a monolithic internal structure. The authors found that regolith can survive even at very small rotation periods, especially in regions at high latitudes. Given the values of low thermal inertia we found at 2016 GE1, which are the most consistent with the thermal inertia of dust-covered asteroids and the results of Sánchez & Scheeres (2020), a plausible explanation for the low thermal inertia we find at 2016 GE1 is a dust layer.

In addition to the two explanations discussed above, another cause of low thermal inertia has emerged recently, namely the cracked surface. Ishizaki et al. (2023) analyzed samples from asteroid Ryugu and found that the thermal inertia of the samples is about 3.5 times larger than the observed thermal inertia of the asteroid Ryugu’s surface. The authors suggested that this difference in thermal inertia between millimeter- and centimeter-sized returned samples and boulders could be due to the presence of large-scale cracks caused by meteor impacts (e.g., Ballouz et al. 2020) and thermal stresses (Molaro et al. 2020; Delbo et al. 2022) on a scale larger than several hundreds of micrometers in rocks and boulders on Ryugu. Therefore, the low thermal inertia of 2016 GE1 could also be due to the cracked surface.

Fig. 6 Dependence of the resulting thermal inertia estimation on the input parameters: comparison with the nominal results. As a guide, in each plot, the results of our nominal estimation are shown as a dotted red histogram. The top panel shows the obtained thermal inertia distribution when the absorption coefficient, α, is set to 0.9 (blue histogram). The middle panel shows the results obtained assuming an emissivity, ε, of 0.9 (black histogram). The bottom panel shows how the thermal inertia distribution changes when a uniform input distribution of density, ρ, is assumed (green histogram).

Fig. 7 Resulting thermal inertial distribution, Γ, with an eccentric (blue histogram) and a circular (red histogram) Yarkovsky model. In both cases, the heat capacity was set to C = 600 J kg−1 K−1. See the main text for more details.

Fig. 8 Posterior (resulting) distribution of the parameters of GE1. The main plot shows the density, ρ, against the diameter distribution, D. The upper and right histograms show the marginal distributions of ρ and D, respectively. The results were obtained using the eccentric Yarkovsky model and a heat capacity of C = 600 J kg−1 K−1. The dashed and dash-dotted lines mark the modes and medians of the distributions, respectively.

6.2 A population of low-thermal-inertia, super-fast-rotating asteroids?

2016 GE1 is not the only super-fast-rotating asteroid with low thermal inertia. As already mentioned, the (499998) 2011 PT shares similar thermal inertia properties, and it has a diameter of about 35 m and rotates in only 11 min. The comparative results for the two objects are shown in Fig. 9. 2016 GE1 is an extreme case of a small and super-fast rotator with extremely low thermal inertia, even lower than that found for 2011 PT. Nevertheless, the values of both objects are low and consistent with a dust-covered surface or a boulder with high micro-porosity. We note that the surface reflectivity properties of 2016 GE1 and (499998) 2011 PT are still unknown. Unfortunately, 2016 GE1 will not pass close to the Earth again in the next future, and therefore, it will not be observable due to its small size. On the contrary, 2011 PT will reach visual magnitudes smaller than 24 in June 2023, 2026, and 2029, and therefore, although still challenging, it may be observed with large-diameter telescopes. Any information in this regard would help us better understand the nature of these unusual objects and unravel their mystery.

In addition, preliminary results obtained by Petković et al. (2021) on the NEA 1998 KY26, the target of the extended Hayabusa 2 mission, showed that it might also have a relatively low thermal inertia similar to that of 2011 PT. Further, by searching for objects with a determined A₂ in the JPL SBDB, many NEAs with |da/dt| > 0.007 au My⁻¹ (a value comparable with the semimajor axis drift of 2011 PT) can be found, and most of them have an absolute magnitude H > 24. Even though information about their rotational state is not always known, the possibility that such small objects are fast rotators is still high, and therefore such a fast semimajor axis drift would still be explained by low thermal inertia. It is also important to note that these A₂ measurements obtained by orbit determination are compatible with the Yarkovsky effect, except for a very small number of cases (Farnocchia et al. 2023). These facts open up the possibility for the existence of a new class of small super-fast-rotating NEAs with low thermal inertia.

Fig. 9 Comparative distributions of thermal inertia for the NEAs 2016 GE1 (dotted blue line) and 2011 PT (dashed red line). In both cases, the results were obtained with the eccentric Yarkovsky model and a heat capacity of C = 600 J kg−1 K−1.

7 Summary and conclusions

In this work we first performed a basic verification of our recently developed statistical MC method for determining asteroid thermal properties, showing that, under the particular cases of superfast rotators with large measured da/dt due to the Yarkovsky effect, we can constrain the surface thermal inertia.

We then used the model to constrain the thermal inertia of the small super-fast rotator 2016 GE1. This NEA has a diameter of less than 20 m, and its rotation period is estimated to be 34 s. We show that the thermal inertia of GE1 cannot exceed Γ = 300 J m⁻² K⁻¹ s^−1/2. The real thermal inertia is likely even smaller, with a probability of >90% of being below 100 J m⁻² K⁻¹ s^−1/2. The extensive testing of different input parameters confirmed the robustness of the result. Therefore, the thermal inertia was constrained to low values with a high probability.

We propose three possible interpretations for the extremely low thermal inertia of 2016 GE1: (i) a high micro-porosity; (ii) the presence of a layer of fine regolith on the surface; (iii) a cracked surface material. Therefore, this work, together with the work of Fenucci et al. (2021), not only demonstrates the usefulness of this alternative method for constraining the thermal properties of asteroid surfaces, but also opens up the possibility of the existence of a potentially new class of NEAs with super-fast rotation and low thermal inertia. Future characterizations and in situ explorations are needed to better understand the physical properties of such objects. In this context, the extended Hayabusa 2 mission will visit the small super-fast-rotating asteroid 1998 KY26 (Hirabayashi et al. 2021) and is expected to provide new insights into these very small asteroids.

Acknowledgements

We appreciate the support from the Planetary Society STEP Grant, made possible by the generosity of The Planetary Society’ members. M.F. and B.N. also acknowledge the MSCA ETN Stardust-R, Grant Agreement no. 813644 under the European Union H2020 research and innovation program.

References

All Tables

All Figures

	Fig. 1 Input density, ρ, versus the diameter, D, distribution for 2016 GE1. The blue histograms at the top and right show the marginal distributions of ρ and D, respectively.
In the text

Fig. 2 MC model tests on asteroid Bennu. The panels show the distributions of thermal inertia solutions for different input parameters assumed to be known. The number of known parameters increases clockwise, starting from the bottom-left panel. For additional details on the test, see Table 6. The gray area in each panel marks the interval of the thermal inertia of Bennu as estimated by Rozitis et al. (2020).

In the text

Fig. 3 Distributions of the thermal parameters of 2016 GE1 for different values of heat capacity, C. The first column shows the distributions of the thermal conductivity, K, and the second column the distributions of the thermal inertia, Γ. Blue histograms are the results obtained using the orbital solution provided by the JPL SBDB, while red histograms are the results obtained using our solution obtained with OrbFit.

In the text

	Fig. 4 Estimated semimajor axis drift for 2016 GE1, obtained for D = 9.8m, C = 800Jkg−1 K−1, γ = 133.77°, and P = 0.009438h. The red level curves represent the Yarkovsky drift measured from astrometry with OrbFit and the corresponding 1σ uncertainty.
In the text

	Fig. 5 Dependence of 2016 GE1’s thermal inertia estimation on the rotation period. The assumed period increases clockwise from the bottom-left panel. The results are shown for the nominal period solution of 34 s, as well as for 2, 5, and 10 times longer periods, as indicated in the plots.
In the text

Fig. 6 Dependence of the resulting thermal inertia estimation on the input parameters: comparison with the nominal results. As a guide, in each plot, the results of our nominal estimation are shown as a dotted red histogram. The top panel shows the obtained thermal inertia distribution when the absorption coefficient, α, is set to 0.9 (blue histogram). The middle panel shows the results obtained assuming an emissivity, ε, of 0.9 (black histogram). The bottom panel shows how the thermal inertia distribution changes when a uniform input distribution of density, ρ, is assumed (green histogram).

In the text

	Fig. 7 Resulting thermal inertial distribution, Γ, with an eccentric (blue histogram) and a circular (red histogram) Yarkovsky model. In both cases, the heat capacity was set to C = 600 J kg−1 K−1. See the main text for more details.
In the text

Fig. 8 Posterior (resulting) distribution of the parameters of GE1. The main plot shows the density, ρ, against the diameter distribution, D. The upper and right histograms show the marginal distributions of ρ and D, respectively. The results were obtained using the eccentric Yarkovsky model and a heat capacity of C = 600 J kg−1 K−1. The dashed and dash-dotted lines mark the modes and medians of the distributions, respectively.

In the text

	Fig. 9 Comparative distributions of thermal inertia for the NEAs 2016 GE1 (dotted blue line) and 2011 PT (dashed red line). In both cases, the results were obtained with the eccentric Yarkovsky model and a heat capacity of C = 600 J kg−1 K−1.
In the text