Asteroids discovered in the Baldone Observatory between 2017 and 2022: The orbits of asteroid 428694 Saule and 330836 Orius

3.1 Starting orbital elements

To compute starting orbital elements of the asteroid Saule, we used the OrbFit software with the JPL DE431 ephemerides, weighting and selecting observations according to the NEODyS-2^[1], 17 additional massive asteroids, and the Yarkovsky effects. We searched for the non-gravitational parameter A 2 , computed directly from observations. A 2 is a non-gravitational transverse acceleration parameter. We used the publicly available OrbFit v.5.0.5 and 5.0.7 software^[2]. These versions can compute orbits using dynamical parameters connected to the non-gravitational perturbations. We calculated the non-gravitational parameter A 2 . We used all 593 observations published by the IAU MPC: https://minorplanetcenter.net/db_search, from 2008 07 25.36573 to 2022 08 19.308671.

According to the OrbFit v.5.0.7, we used the error model vfcc17 of Veres et al., (2017). Table 3 presents the starting nominal Keplerian orbital elements of the 428694 Saule (2008 OS9) asteroid computed using the mentioned OrbFit v.5.0.7 software.

It is worth noting that two other interesting recent methods of computations of the non-gravitational parameter A 2 are in Perez-Hernandez and Benet (2022) and Dziadura et al., (2022).

3.2 Starting physical parameters

Next, we used the mercury software of Fenucci and Novakovic (2022). As stated in https://github.com/Fenu24/mercury, this is a modified version of the N-body code mercury by Chambers (1999) (https://github.com/Fenu24/mercury/tree/master/doc), which includes the Yarkovsky and YORP effects for the dynamics of small solar system objects.

The Yarkovsky effect was discovered by the Polish-Russian civil engineer Ivan Osipovich Yarkovsky: https://en.wikipedia.org/wiki/Yarkovsky_effect.

We have a diurnal effect and a seasonal effect. Both effects produce a shift in the semimajor axis of the asteroid. The influence of the Yarkovsky effect on the motion of an object is evident over a long time, for example, in the motion of members of asteroid family e.g., Novakovic et al., (2022), Wlodarczyk and Leliwa-Kopystynski (2018), studying retrograde orbits in Wlodarczyk (2022b) or dangerous asteroids (Wlodarczyk 2020). To study the motion of asteroids and comets with the influence of the Yarkovsky effect, we use the non-gravitational parameters A 2 and d a ∕ d t . Both values of A 2 and d a ∕ d t can be determined during orbit determination as a seventh parameter together with six other orbital parameters, i.e., Keplerian elements. In Wlodarczyk (2018), several non-gravitational parameters A2 and d a ∕ d t were computed. The advantage of the method used in the present work is that it simultaneously shows the evolution of P , gamma, and d a ∕ d t for the 2008 OS9 asteroid clones. The authors of the mentioned mercury software have published the first papers: Novakovic and Fenucci (2022) and Fenucci and Novakovic (2022).

As physical parameters in starting file yarkovsky.in, we used typical values of this type of asteroid, as was partially presented in Fenucci and Novakovic (2022):

• the density ρ (kg/m³) = 2,120;

• the thermal conductivity K (W/m/K) = 0.001;

• the heat capacity C h (J/kg/K) = 800;

• the diameter D (m) = 600;

• the obliquity γ ( ° ) = for the first clone = 0, and of the following 100 clones differ from the previous one by 180 ∕ 10 1 ° = 1.78 2 ° ;

• the rotation period P (h) = 8.43; taken from our previous work Cernis et al., (2010)

• the absorption coefficient α ; (usually set to 1);

• the emissivity ε (usually set to 1);

• yorp_flag = 1

• stoc_yorp_flag = 1

• step_user = 1.0

• dt_out = 1.0

• one year.

An assessment of the diameter of the asteroid Saule is in our previous article by Cernis et al., (2010), where we used absolute magnitude H = 19.42 and different values of geometric albedo pV. The diameter of pV = 0.04 (C-type asteroid) is 867 m; for pV = 0.20 (S-type asteroid), it is 388 m. The true diameter of 2008 OS9 is between these two values.

In computations of the Yarkovsky and YORP effects, we accept the asteroid has a mean diameter of 600 m, and its mean albedo is 0.15. Hence, we take the mean diameter of the asteroid Saule in our computation.

We can check NEOWISE observations of NEOs for some infrared data on Saule. Sometimes, this can help in the determination of the geometric albedo. In our case, we did not check NEOWISE data. We assume the albedo pV in the interval 0.05 ÷ 0.50 , as is shown, for example, in Gustafsson et al., (2019). Hence, we can compute the possible interval of the effective diameter of Saule: having both H and pV estimated as: H = 19.39 ± 0.5 mag, pV = 0.27 ± 0.22 (both are systematic uncertainties, not standard deviations), one can proceed and estimate the effective diameter, D , which is 0.2 < D < 1 km .

From a discussion above, it follows that the effective diameter of Saule can be determined only approximately. For this reason, we will do computations for a set of values presented in Table 4.

In turn, Table 4 indicates the calculation of the diameter of the asteroid Saule based on the range of observed H for 19.34, 19.39, and 19.44 mag and possible albedo, pV of 0.49, 0.27, and 0.05. For the diameters obtained in this way, we calculated the rotation periods of the asteroid for different bulk density densities of 1,000, 1,500, 2,000, and 2,500 kg / m 3 and for two mass cohesive trenches 10 and 100 Pa, respectively. We calculated the diameters of asteroid Saule according to the formula from Bowell et al., (1978).

We then calculate the critical rotation period given by Hu et al., (2021):

where ρ is the bulk density, C is the bulk cohesion, and D is the diameter. According to Fenucci et al., (2021), k = 0.9114 is a constant parameter, see Supplementary Material Pravec et al., (2010).

Results are presented in Table 4. It turned out that the minimum periods of critical rotation are just over 2 h, for a density of 1,000 kg/ m 3 and bulk cohesion 100 Pa, and maximum rotation periods of over 30 h occur for bulk density 2,500 kg/ m 3 and bulk cohesion 10 Pa.

3.3 Long-time orbital and physical evolution

First, we computed orbital elements and non-gravitational elements of 101 clones of the asteroid 428694 Saule (2008 OS9) for epoch JD2459600.5 and the planets from Mercury to Neptune. We used the parameter σ LOV = 3 according to the multiple solution methods (Milani 2006, Milani et al., 2005a, Milani et al., 2005b), i.e., the starting values of d a ∕ d t were obtained from the A 2 parameter and its uncertainties.

Next, the starting orbital elements of each clone are treated as inputs in one small.in the file, and all of them, one yarkovsky.in a file is included, the same for all clones. These are the inputs to the main mercury integrator program. Of course, you can associate each clone’s orbital elements with the same yarkovsky.in file. This way, we will have 101 input files with 101 data, yielding 101 × 101 calculation results. This would be a method close to a random method like the Monte Carlo method. However, this requires long, laborious calculations, and the line of variation (LOV) method of Milani avoids these difficulties. It does not require so many calculations.

Our additional computational tests showed that the computational results of one run of 101 clones yield similar to 101 × 101 computational results. So, the initial simplification of the number of integrations does not affect the result of calculations of the evolution of rotation, spin, and the non-gravitational parameter d a ∕ d t . Next, we integrated the equation of motions of asteroids and planets 0.5 Myr forward.

Figure 1 presents starting non-gravitational parameter d a ∕ d t of 101 clones of the asteroid 428694 Saule (2008 OS9) x computed from the OrbFit software. They are from ( − 0.004 to +0.006) au/Myr.

Figure 1

Starting non-gravitational parameter d a ∕ d t of 101 clones of the asteroid 428694 Saule (2008 OS9).

Figure 2 presents the physical parameters of asteroid 428694 Saule (2008 OS9) after 0.5 Myr forward integration. They are the rotation period P (h), the obliquity γ ( ° ), and the value of d a ∕ d t connected with the non-gravitational parameter A 2 . Only those clones with P ≤ 10 h are marked in the left panel. There are 51 of them, so just over half. The other 50 have a rotation period, P ≥ 1,000 h, and hence were assumed to have P = 1,000 h. According to the classification of the program, they slow down.

Figure 2

Physical parameters of asteroid 428694 Saule (2008 OS9) after 0.5 Myr forward integration. Parameter d a ∕ d t is multiplied by 1 0 4 .

Hence, this state has a minor significance for the Yarkovsky effect, the semimajor axis drift d a ∕ d t goes to 0. We can see symmetry and possible correlation between spin axis gamma and non-gravitational parameter d a ∕ d t after 0.5 Myr forward integration. Those with gamma = 0 ° may give d a ∕ d t = 0 or d a ∕ d t > 0 . In comparison, those with gamma = 18 0 ° may give d a ∕ d t negative or equal to 0. Hence, having the value of the gamma parameter, we can see something about the non-gravitational parameter d a ∕ d t .

Comparing d a ∕ d t in Figure 2 (bottom panel), with d a ∕ d t from Figure 1, we can see that the amplitude of changes of the d a ∕ d t of the clones of the asteroid 428694 Saule (2008 OS9) after 0.5 Myr is about 20 times smaller than the starting one, that is, from a range of about ± 0.004 au/Myr to about ± 0.0008 au/Myr.

It is interesting that after a period of 0.5 Myr, the clones’ d a ∕ d t values are around the borderline values − 0.0008 au/Myr and +0.0008 au/Myr, and also around d a ∕ d t = 0 .

Hence, the dynamics of changes in the d a ∕ d t value decreases, and the clones seem to slow down, stay close to a particular achieved place, and move slower in both directions of their semimajor axis. The second part of the clones is around the value d a ∕ d t = 0 au/Myr, meaning that they practically do not change their positions on the semimajor axis related to the Yarkovsky effect. Such stratification of clones on the d a ∕ d t axis begins to appear after about 100,000 years and is visible after 0.5 Myr.

Interesting is the distribution of the clonal rotation period. As shown in the left panel of Figure 2, the rotation period of the clones P is generally between 1 and 2 h, from the original value P , of 8.43 h. So, the asteroid clones sped up their rotation period by almost nine times.

Given the diameter of the asteroid Saule (D = 600 m) and its assumed density of 2,120 kg/ m 3 , the asteroid is most likely a rubble pile held by gravity forces. The shortest rotation period for such objects is P = 2.1  h, at which they experience resurfacing (e.g., creation of a typical top-hat shape), mass shading, or even rotational fission into several smaller objects. This fact has been considered by Fenucci and Novakovic (2022), who set a limit at some critical period P crit = 2.34  h for his asteroid with a diameter of 2 km.

They write: “When the rotation period becomes smaller than P crit we assume that a fission event takes place and we re-initialize the spin state. In this process, we do not simulate production of a binary system, and we assume that the mass lost is small enough to not significantly change the equivalent diameter of the object. During the fission event, we assume the obliquity g a m m a to be unchanged. On the other hand, the spin rate is decreased according to the momentum carried away by the ejected mass.” As a result, the rotation period is never smaller than P crit . Note that during test runs of the Mercury N-body code, periods were never smaller than this critical value, see Figure 4 in Fenucci and Novakovic (2022).

Figure 3

Left panel: Rotation period P vs spin axis of all clones of asteroid 428694 Saule (2008 OS9); 50 of them reached rotation period P = 1,000 h. Right panel: dependence of rotational period shorter than 80 h of 51 clones vs spin axis. Both panels present results of 0.5 Myr forward integration.

Similarly interesting is the distribution in the inclination of the rotation axis of the clones to the ecliptic plane. Initially, we assumed that the neighboring clones differ by 180 ∕ 10 1 ° in the gamma angle, i.e., in the inclination of the rotation axis. There is a visible polarization to 0, 90, and 18 0 ° .

Hence, our computations showed that the Yarkovsky effect plays a significant role after 0.5 Myr, and even signs of it are already visible after 0.1 Myr.

Figure 3 presents the rotation period P vs. spin axis of all clones of asteroid 428694 Saule (2008 OS9), 50 reached rotation period P = 10,000 h (left panel) and dependence of rotational period shorter than 80 h of 51 clones vs. spin axis (right panel). Both panels present results of 0.5 Myr forward integration. It is visible that clone groups near 0, 90, and 18 0 ° . In particular, this is visible in the case of clones that have reached their period’s rotation equal to 1,000 h.

Figure 4

Known physical parameters of 31079 asteroids according to the MPC.

Table 5 shows the number of the first group of clones for the three spin-axis groups. It turned out that they take only three spin-axis values: 0, 90, and 18 0 ° (Table 6).

3.4 Known rotation periods of asteroids

According to the recent reliable database of the asteroid rotation periods: https://minplanobs.org/mpinfo/datazips/LCLIST_PUB_2023JUN.zipthere are 32778 asteroids with defined physical parameters.

According to the database, they introduce the U code, which assesses the quality of the period solution, not necessarily of the data per se, as stated there.

If we selected asteroids with criterion U ≥ − 2 , we got 31,079 asteroids. Figure 4 presents their rotational period vs. diameter. For visibility, only every fifth asteroid is presented.

Additionally, if we count objects which fulfill the following criteria:

H : < 20 mag,

P : < 2  h we got only a few selected asteroids. Hence, the asteroid 2008 OS9, NEA, with starting q = 0.5 au, may have during his evolution the theoretical rotation period around 1 h as we computed using the mercury software.

4.1 Starting orbital elements

KC and IE discovered Asteroid 330836 Orius (2009 HW77) at Baldone on 2009-04-25. An orbit type of asteroid is a Distant object. To compute starting orbital elements of the asteroid Orius, we used the OrbFit software with the JPL DE431 ephemerides, weighting and selecting observations according to the (NEODyS-2)^[3], 17 additional massive asteroids. We did not search for the non-gravitational parameter A 2 because of the long distance of the asteroid to the Sun and the high absolute magnitude of the asteroid. We used all 85 observations published by the IAU MPC from 2002 02 12.50951 to 2012 05 19.44276.

According to the OrbFit v5.0.7, we used the error model vfcc17 (Veres et al., 2017). Table 5 presents the initial nominal Keplerian orbital elements of the 330836 Orius (2009 HW77) asteroid.

4.2 Long-time orbital evolution

The dynamical lifetime is the time between the start of integration and ejection from the integration. The studied asteroid 330836 Orius (2009 HW77), like other Centaurs, are planet-crossing object and has relatively short dynamical lifetimes of around 1 Myr (Horner et al., 2004) and short Lyapunov time, e.g., of about 4,260 year for Centaur asteroid (463,368) 2012 VU12 (Wlodarczyk et al., 2017).

Due to the lack of sufficiently precise physical parameters of the asteroid 2009 HW77, to study its orbital evolution we used SWIFT software, initially developed by H. Levison (https://www.boulder.swri.edu/∼hal/swift.html).

We included perturbations of all planets from Mercury to Neptune. We computed clones using the multiple solution methods for the 1 σ uncertainty around the nominal orbit of the asteroid. We computed the time evolution of all 1,001 starting clones during 100 Myr of forward integration. The integration is similar to that of Wlodarczyk et al., (2011). Clones that reached the ejection distance (1,000 au) or that impacted the Sun or a planet were removed from our integration. The computed half-time, i.e., when half of the clones are rejected, is 4.45 Myr, as shown in Figure 5. Wlodarczyk et al., (2011) computed a half-time of 5.25 Myr in our previous paper. The differences result from a more accurate orbit based on a longer observation arc. Also, previously starting orbital elements of the asteroid are computed only with added Ceres, Pallas, and Vesta, now 17 additional massive asteroids.

Figure 5

The Centaur 2009 HW77. The behavior of the lifetime of the remaining clones has an exponential shape. Hence we can compute the half-time of 1,001 clones in the forward 100 Myr RA15 integration. The computed half-time is 4.45 Myr.