Physics Research International
Volume 2014 (2014), Article ID 561680, 13 pages
http://dx.doi.org/10.1155/2014/561680

Research Article

The Potentially Dangerous Asteroid (101955) Bennu

I. Włodarczyk

Chorzów Astronomical Observatory, WPKiW, 41-500 Chorzów, Poland

Received 6 June 2014; Accepted 17 November 2014; Published 15 December 2014

Academic Editor: Sergey B. Mirov

Copyright © 2014 I. Włodarczyk. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

We computed impact solutions of the potentially dangerous asteroid (101955) Bennu based on 569 optical observations from September 11.40624 UTC, 1999 to January 20.11189 UTC, 2013, and 29 radar observations from September 21, 1999, through September 29, 2011. Using the freely available OrbFit software package, we can follow its orbit forward in the future searching for close approaches with the earth, which can lead to possible impacts up to 2200. With the A2 nongravitational parameter in the motion of the asteroid (101955) Bennu we computed possible impact solutions using different JPL planetary and lunar ephemerides and different number of additional massive perturbed asteroids. The possible impact path of risk for 2175 is presented. Additionally, we computed possible impact solutions using the normal places method of the selection of Bennu’s astrometric observations. Moreover, we computed time evolution of the mean orbital elements and the orbital nodes of Bennu 5 kyr in the backwards and 1 kyr in the future using the Yarkovsky effects. We computed the mean motion and secular orbital resonances of the Bennu. We also computed the influence of the JPL planetary and lunar ephemerides DE403, DE405, DE406, DE414, and DE423 on the close approaches of the asteroid (101955) Bennu with the earth.

1. Introduction

Discovered on September 11, 1999, by the Lincoln Laboratory Near-Earth Asteroid Research Team at Socorro, Bennu was an Egyptian mythological figure associated with Osiris, Atum, and Ra. This minor planet is the target of the OSIRIS-REx sample return mission (http://www.nasa.gov/osiris-rex). OSIRIS-Rex’s Touch-and-Go Sampler evokes Bennu’s image as a heron. The spacecraft is scheduled to launch in 2016, reach 101955 Bennu in 2019, and return samples to earth in 2023.

It is a potential Earth impactor and now, as of March 17, 2014, is listed on the Sentry Risk Table of the JPL NASA with the third highest rating on the Palermo Technical Impact Hazard Scale (PS = −1.71 (cumulative) and −2.32 (maximum)), after the asteroid (29075) 1950 DA (PS = −0.83 (cumulative and maximum)) and 2007 (−1.56 and −1.57), respectively (http://neo.jpl.nasa.gov/risk/). Moreover, the asteroid (101955) Bennu is listed on the second place on the CLOMMON2 site ( (maximum)) provided by the NEODyS (http://newton.dm.unipi.it/neodys/index.php?pc=4.1).

Up to now, there were many methods of computing possible impact solutions of the asteroid (101955) Bennu. In Milani et al. [1] there are computed eight potential earth impacts of the asteroid (101955) Bennu between 2169 and 2199. They found the best fit value of  au/My which give the minimum value of the rms for all observations of the asteroid (101955) Bennu. Moreover, they found possible impacts between 2169 and 2199 using different values of in the range  au/My.

Vokrouhlický et al. [2] investigate the possibility of detecting the Yarkovsky effect via precise orbit determination of the near-earth asteroids on example of the asteroids 6489 Golevka, 1620 Geographos, 1566 Icarus, and 1998 with different values of the parameter . Their figures present 3 uncertainty ellipsoids with and without the Yarkovsky effects projected onto the range () versus range rate () plane for the close approach of different asteroids.

Farnocchia et al. [3] computed for asteroid (101955) 1999 Bennu  au/My and  au/My. Tables 2 and 3 list some physical parameters and nongravitational parameters, and , for 37 asteroids.

Vokrouhlický et al.[4] present explicit, analytic formulas for the solar radiation dynamical effects on the orbits of two near-earth asteroids: (1566) Icarus and (6489) Golevka. They plot the 3 confidence ellipses of the Icarus and Golevka orbit uncertainty onto the space of radar observables: the geocentric distance (in km) and the rate of change of the geocentric distance (in km/day) computed for the nominal model, that is, without the radiation effects, and for the extended model, that is, including the radiation effects. These simulated orbit displacements will be useful for the search of the nongravitational parameters in the motion of asteroids.

Farnocchia and Chesley [5] collect equations for two nongravitational parameters, , that is, nongravitational radial acceleration parameter, and , that is, nongravitational transverse acceleration parameter. They also computed value of for the asteroid (29075) 1950 DA from theory and using astrometric observations.

Farnocchia et al. [3] concluded that the Yarkovsky effect for the asteroid (99942) Apophis with two adopted values of :  au/ and  au/ cannot be detectable with the radar apparition in 2013. They also list 25 main belt asteroid perturbers and associated values which we used in our orbital computations with the OrbFit software. denotes mass expressed as a product of the mass () and gravitational constant () in units km³/s².

Chesley et al. [6] contain a detailed analysis of the relevance of the different components of the dynamical model, a determination of Bennu’s density and mass, and a rigorous statistical analysis of the potential impacts with associated impact probabilities. We used their nongravitational parameter  au/d² in our computations of possible impact solution for the asteroid (101955) Bennu.

Moreover, to compute possible impact solution of the asteroid (101955) Bennu with the earth we take into account the OrbFit software with different JPL planetary and lunar ephemerides, perturbations of different number of additional massive asteroids, weighting and selection of observations according to the NEODyS, the error model based on Chesley et al. [7], the Yarkovsky effects, and the normal places method.

There are no differences in observational material between us and Chesley et al. [6]. But in some case, when the observational material is dense, they used new weights. Similarly, in our work we applied additionally the normal places method; see Section 3.1.

Chesley et al. [6] present computations based on two JPL ephemerides: DE405 and DE424.

We computed possible impacts with many different JPL ephemerides: JPL DE403, DE405, DE406, DE414, and DE423. Moreover, we added 25 perturbing asteroids when Chesley et al. [6] used 16 asteroids. Also our computations of the earliest possible impacts were made for 0, 4, 16, and 25 massive perturbing asteroids—see Table 2.

We also draw the path of risk of the asteroid Bennu in 2175.

We show that errors of propagated orbital elements of the asteroid Bennu are great, mainly in 2175, after CA with the earth in 2135. We computed these errors and present in Table 6. Hence impact hazard after 2135 can only be explored through statistical means, as Chesley et al. [6] state.

Our computations of CAs of the asteroid Bennu with the earth till 2135 give almost the same results as Chesley et al. [6].

Moreover, we extend orbital evolution of the asteroid Bennu from 2000 to 2140 as was explored in Chesley et al. [6] to time-span 5000 yr backwards and 1000 yr forward from the osculating epoch 2014 May 23.0 = JDT 2456800.5. Position of the asteroid Bennu in its orbit is chaotic but the orbit in the space is stable over several thousands of years as was depicted in Section 6 and in Figure 2 where the stable mean values of orbital elements are presented.

2. The Initial Orbital Elements of Asteroid (101955) Bennu

Currently, March 17, 2014, there are 10708 near-earth asteroids ( au): 832 Atens with orbits similar to that of 2062 Aten ( au;  au), 5345 Apollos with orbits crossing the earth’s orbit similar to that of 1862 Apollo ( au;  au), and 4531 Amors with orbits similar to that of 1221 Amor (1.017 au  au). They are listed at the Minor Planet Center (MPC) (http://www.minorplanetcenter.net/iau/lists/Unusual.html) and at the JPL NASA (http://neo.jpl.nasa.gov/stats/).

The MPC and the JPL NASA classified asteroid (101955) Bennu as an Apollo-class object and additionally as a potentially hazardous asteroid. Bennu belongs to one of 1457 potentially hazardous asteroids (PHA) (http://www.minorplanetcenter.net/iau/lists/Dangerous.html). According to the Minor Planet Center (MPC), the PHA are objects with brighter than and an earth MOID less than 0.05 au. MOID, Minimum Orbit Intersection Distance, is the minimum distance between the orbit of the earth and the minor planet. According to the JPL NASA (http://neo.jpl.nasa.gov/risk/a101955.html), its absolute magnitude is 20.9, with a diameter of about 0.560 km and a rotation period equal to 4.288 h. According to Chesley et al. [6] the mass of Bennu is equal to  kg and the bulk density is  kg/m³.

Table 1 lists the computed orbit of the asteroid (101955) Bennu published by different sources: the MPC in M.P.O. 251037 (http://www.minorplanetcenter.net/iau/ECS/MPCArchive/2013/MPO_20130127.pdf), the NEODyS-2 (http://newton.dm.unipi.it/neodys/index.php?pc=1.1.0&n=101955), the JPL Small-Body Database (http://ssd.jpl.nasa.gov/sbdb.cgi?sstr=101955;orb=1), and the author.

Table 1: The Keplerian orbital elements of asteroid (101955) Bennu from different current sources.

(au)

(deg)
(deg)
(deg)
(deg)

MPC, rms = 0.44′′

Epoch 2013 Nov. 04.0 = JDT 2456600.5

1.1260425
0.2036916
6.03472
2.03947
66.28127
237.78960

NEODyS-2, rms = 0.4496′′

Epoch 2014 May 23.0 = JDT 2456800.5

1.1260011
0.2037022
6.03486
2.03688
66.29364
42.75424

JPL 85, rms = 0.4482′′, = −4.664 × 10 au/d

Epoch 2011 Jan. 01.0 = JDT 2455562.5

1.1263910
0.2037451
6.03494
2.06087
66.22308
101.703944

IW, rms = 0.4470′′, = −4.618 × 10 au/d

Epoch 2011 Jan. 01.0 = JDT 2455562.5

1.1263910
0.2037451
6.03493
2.06087
66.22308
101.703942

Table 2: The asteroid (101955) Bennu. The earliest possible impacts computed with different number of perturbing asteroids. The JPL DE406 Solar system model. = −4.618  ×   au/d².

Date, UTC
Sigma
Impact
Number of

LOV
probability
asteroids

2175/09/25.167
−0.642

0

2175/09/25.153
−0.317

0

2175/09/25.154
−2.674

4

2175/09/25.171
−2.349

4

2175/09/25.164
−2.344

16

2175/09/25.159
−2.654

25

2175/09/25.162
−2.329

25

2175/09/25.150
1.340

NEODyS

2175/09/25.160
1.050

NEODyS

2176/09/24.389
0.949

0

2176/09/24.403
−1.083

4

2176/09/24.389
−1.079

16

2176/09/24.389
−1.063

25

2176/09/24.390
2.460

NEODyS

2180/09/24.392
−1.236

0

2180/09/24.358
−0.820

0

2180/09/24.355
−2.848

16

2180/09/24.351
−2.832

25

2180/09/24.360
0.890

NEODyS

2180/09/24.390
0.530

NEODyS

Table 1 presents orbital elements: : semimajor axis, : eccentricity, : inclination, : longitude of the ascending node, : argument of perihelion, and : mean anomaly. The orbits are computed from 569 astrometric positions from which 8 observations were rejected as outliers, and also on 29 radar observations with 4 observations rejected as outliers. The orbits are based on observations from 1999 Sept. 11.40624 UTC to 2013 Jan. 20.11189 UTC.

We can see in Table 1 that orbits of the asteroid (101955) Bennu computed by the JPL and the author are almost the same. Only differences are in value of the inclination, of about 0.00001 deg, that is, only 0.04′′, and in the mean anomaly, 0.00002 deg and 0.07′′, respectively. Differences are connected with the different Solar system model, DE431 by the JPL, and DE406 by the author, the number of the additional perturbing massive asteroids, 16 by the JPL, 25 by the author, and different values of the nongravitational transverse acceleration parameter,  au/d² by the JPL, and  au/d² by the author.

Our nongravitational parameter and masses of the 25. additional perturbing massive asteroids are taken from Chesley et al. [6]. Their orbit of the asteroid (101955) Bennu is close to the JPL orbit solution 85 (Nonlinear Yarkovsky Model) and to the JPL 87 (Transverse Yarkovsky Model) as is shown in Table 3 and is similar to our orbit. Hence, we used our orbit from Table 1, as a starting orbit in suitable orbital computations. Moreover, we can compute different orbital behavior of the asteroid (101955) using our nongravitational model and the freely available OrbFit software v.4.2 (http://adams.dm.unipi.it/~orbmaint/orbfit/). This new version includes the new error model based on Chesley et al. [7]. In all our computations, we follow the same method of the weighting and selection of observations that is being used by the NEODyS site [8].

To compute and propagate the orbit, the OrbFit software used internal control of propagation methods (multistep, Runge-Kutta-Gauss, and Everhart). Usually we used automatic control (multistep for main belt, Everhart for high eccentricity, and/or planet crossing).

3. Possible Impact Solutions

To compute possible impact solutions we used different number of perturbing asteroids, different JPL planetary and lunar ephemerides, and additionally the normal places method for selecting of observational material of the asteroid Bennu. We computed possible impact solutions of Bennu up to 2199.

First, we computed impact solutions of the asteroid (101955) Bennu using different number of perturbing asteroids. In each case we computed 2001 VAs of the Bennu using 3  uncertainty and multiple solution method of [9, 10]. Also, we used the JPL DE406 and the nongravitational parameter  au/d². Starting position and velocities of the perturbing asteroids were computed from the ASTDyS-2 (http://hamilton.dm.unipi.it/astdys/), and their masses were taken from Chesley et al. [6]. Also nongravitational parameter was taken from Chesley et al. [6].

Table 2 lists the earliest possible impacts of the asteroid (101955) Bennu using different number of perturbing asteroids. 0 denotes no asteroids, 4 as perturbers are asteroids: (1) Ceres, (2) Pallas, (4) Vesta, and (10) Hygiea, additions for 16 asteroids: (3) Juno, (6) Hebe, (7) Iris, (15) Eunomia, (16) Psyche, (29) Amphitrite, (52) Europa, (65) Cybele, (88) Thisbe, (511) Davida, and (704) Interamnia, and additions for 25 asteroids: (11) Parthenope, (14) Irene, (56) Melete, (63) Ausonia, (135) Hertha, (259) Aletheia, (324) Bamberga, (419) Aurelia, and (532) Herculina (see Table 6 in Chesley et al. [6]). NEODyS denotes published possible impact solution by the NEODyS (http://newton.dm.unipi.it/neodys/index.php?pc=1.1.2&n=101955) where the NEODyS’s Risk Table has not been computed with the standard CLOMON2 software but by processing the output from a Montecarlo run in the 7-dimensional space of the orbital elements and the secular perturbation on semimajor axis.

Table 2 lists computed possible impact solution: the calendar date (UTC) of the potential impact, Sigma LOV, that is, the coordinate along the Line Of Variations (LOV)—the further from zero, that is, for nominal orbit, the less likely the impact, impact probability, and the number of perturbing asteroids.

It is visible that the computed impact solutions without perturbing asteroids differ from solutions with perturbed asteroids. Also, there are small differences in possible impact solutions using different number of asteroids. Our results with 25 perturbing asteroids are in good agreement with the results of the NEODyS.

Next, we computed impact solutions of the asteroid (101955) Bennu using different JPL ephemerides. We selected five JPL planetary and lunar ephemerides from (ftp://ssd.jpl.nasa.gov/pub/eph/planets/README.txt).

They are according to the description as follows.

DE403: it was created in May 1993, includes nutations and librations, covers JED 2305200.5 (1599 Apr. 29) to JED 2524400.5 (2199 Jun. 22), and is fit to planetary and lunar laser ranging data [11].

DE405: it was created in May 1997, includes both nutations and librations, Referred to the International Celestial Reference Frame, and covers JED 2305424.50 (1599 Dec. 09) to JED 2525008.50 (2201 Feb. 20).

DE406: it was created in May 1997, includes neither nutations nor librations, and spans JED 0625360.5 (−3000 Feb. 23) to 2816912.50 (+3000 May 06) This is the same integration as DE405, with the accuracy of the interpolating polynomials, has been lessened to reduce file size for the longer time-span covered by the file.

DE414: it was created in May 2005, includes nutations and librations, covers JED 2414992.5 (1899 Dec. 04) to JED 2469872.5 (2050 Mar. 07), and is fit to ranging data from MGS and Odyssey through 2003 [12].

DE423: it was created in February 2010, includes nutations and librations, covers JED 2378480.5 (1799 Dec. 16) to JED 2524624.5 (2200 Feb. 02), and is intended for the MESSENGER mission to Mercury.

The last released ephemerides are not studied in our work.

DE430: it was created in April 2013, includes librations and 1980 nutation, Referred to the International Celestial Reference Frame version 2.0, and covers JED 2287184.5 (1549 Dec. 21) to JED 2688976.5 (2650 Jan. 25).

DE431: it was created in April 2013, includes librations and 1980 nutation, covers JED −0.3100015.5 (−13200 Aug. 15) to JED 8000016.5, (17191 Mar. 15).

The above given ephemerides were downloaded from (ftp://ssd.jpl.nasa.gov/pub/eph/planets/Linux/).

Table 3 lists the earliest possible impacts of the asteroid (101955) Bennu using different JPL ephemerides (DE406, DE403, DE405, DE414, and DE423). It is visible that using different JPL planetary and lunar ephemerides we computed almost the same possible impact solutions. Moreover, they are similar to the NEODyS results.

Table 3: The asteroid (101955) Bennu. The earliest possible impacts computed with different JPL ephemerides. 25 perturbing asteroids were added. = −4.618  ×   au/d².

Date, UTC
Sigma
Impact
The JPL

LOV
probability
ephemerides

2175/09/25.159
−2.654

DE406

2175/09/25.162
−2.329

DE406

2175/09/25.168
−2.667

DE403

2175/09/25.167
−2.342

DE403

2175/09/25.153
−2.683

DE405

2175/09/25.171
−2.357

DE405

2175/09/25.153
−2.643

DE414

2175/09/25.159
−2.317

DE414

2175/09/25.154
−2.527

DE423

2175/09/25.153
−2.201

DE423

2175/09/25.150
1.340

NEODyS

2175/09/25.160
1.050

NEODyS

2176/09/24.389
−1.063

DE406

2176/09/24.396
−1.076

DE403

2176/09/24.389
−1.092

DE405

2176/09/24.400
−1.052

DE414

2176/09/24.389
−0.936

DE423

2176/09/24.390
2.460

NEODyS

2180/09/24.351
−2.832

DE406

2180/09/24.352
−2.845

DE403

2180/09/24.358
−2.861

DE405

2180/09/24.351
−2.821

DE414

2180/09/24.360
0.890

NEODyS

2180/09/24.390
0.530

NEODyS

3.1. The Normal Places Method

Next we computed possible impact solutions of the asteroid (101955) Bennu using for the observational material the normal places method.

Looking for the observational material of the asteroid (10195) Bennu we find many close packed astrometric observations from one observatory. Hence they can have considerable influence on selection and weighting of the observational material and on accumulation of the observational errors from a given observatory. To prevent these effects we used the normal places method.

As was explained in Chesley et al. [13], in the case of the asteroid (101955) 1999 RQ36 Bennu, there are batches containing an excess of observations from a single observatory in a single night. To reduce the effect of these batches to a preferred contribution of 5 observations per night, they relaxed the weight by a factor , where is the number of observations contained in the batch. Similar method is used in Chesley et al. [6].

According to the Minor Planet Center’s (MPC) Guide to Minor Body Astrometry (http://www.minorplanetcenter.net/iau/info/Astrometry.html#numobs) making more than three observations per objects per night—it is a waste of observer’s time and rarely helps the orbit solution.

We used the method of computing normal places based on the online service provided by the Minor Planet Center (MPC) (http://www.minorplanetcenter.net/iau/VideoObs/VideoNormalPlaces.html). This service is a tool to compute normal places for video observations which must be reported as “normal” positions as the MPC states. A normal position is a position that represents an average of several positions, referenced to a specified time.

First we gathered observations for a given observatory. For each observatory we collect observations made during one to two days. Hence, we got 93 groups of observations with more than two observations in each group. Next we computed normal places for these groups. We computed normal places for minimum 3 observations in each group. Hence we got 93 normal places, 14 observations without computing normal places because there were only two observations in each group, and 29 radar observations remaining without normal places method. Then we computed orbit of the asteroid (101955) Bennu using its observations selected and computed with the normal places method and got rms = 0.4249′′. Without normal places method we have rms = 0.4470′′. Both rms-es values are computed for the orbit using JPL DE406, 25 perturbed asteroids and nongravitational parameter  au/d² taken from Chesley et al. [6].

Using the normal places method we computed many possible impacts between 2161 and 2199 with probabilities of these impacts, between in 2192 and . The earliest impact solutions with the greatest probabilities are listed in Table 4. It is interesting that impact solutions computed with the selection and weighting method used by the NEODyS give similar results as in the case of the normal placed method for years 2175–2180.

Table 4: The asteroid (101955) Bennu. The earliest possible impacts computed with the normal places method. The JPL DE406, 25 perturbing asteroids, and = −4.618  ×   au/d² are used.

Date, UTC
Sigma
Impact

LOV
probability

2167/09/25.186
0.849

2168/09/24.436
1.288

2172/09/24.407
−1.346

2175/09/25.162
−1.068

2175/09/25.159
0.248

2176/09/24.404
−1.147

2176/09/24.389
0.321

2179/09/25.152
−1.463

2179/09/25.151
0.573

2179/09/25.135
0.600

2180/09/24.397
−1.650

2180/09/24.435
−1.611

2180/09/24.365
−0.498

2180/09/24.314
−0.486

2180/09/24.310
−0.301

2180/09/24.346
−0.289

2180/09/24.393
0.196

2180/09/24.353
0.220

2180/09/24.382
0.758

4. Impact Orbits

From Table 2 we take two earliest possible impacts of the asteroid (101955) Bennu in 2175 computed using the OrbFit software with the error model, 25 additional perturbing massive asteroids, the multiple solution method with and 2001 VAs, the JPL DE406 Solar system model, and nongravitational parameter  au/d². The weighing and selection method of the NEODyS was used.

The first possible impact, called a, is on 2175 Sept. 25.159 UTC with probability of about ; the second one, b, occurs on 2175 Sept. 25.162 UTC with probability of about .

We computed impact orbit for the above given moments in 2175.

Impact orbits are the orbits computed for the initial epoch and for 7 days before day of impact. They hit the earth.

To compute possible impact solutions of the asteroid (101955) Bennu in Table 2 we used 2001 VAs using uncertainty and the multiple solution method. First we take the orbital elements of the VA which has the closest approach to the earth around two dates of the possible impact, a and b.

They are so-called prognostic impact orbits. Orbital elements with 1  uncertainties of these orbits are listed in Table 5. We can see that these orbits differ only by the small values; for example, mean anomaly differs by 0.00234′′.

Table 5: The asteroid (101955) Bennu. The initial prognostic impact orbits for two possible impacts in 2175.

(au)

(deg)
(deg)
(deg)
(deg)

Epoch 2014 May 23.0 = JDT 2456800.5

Impact a—2175 Sept. 25.159 UTC

1.126001049342
0.203702111
6.0348667
2.0368635
66.2936367
42.7542947

Impact b—2175 Sept. 25.162 UTC

1.126001049356
0.203702117
6.0348660
2.0368642
66.2936374
42.7542940

Table 6: The asteroid (101955) Bennu. The initial and 7 days before impact orbit for 2175: cases a and b. For comparison orbit before 2135 CA with the earth—case b.

(au)

(deg)
(deg)
(deg)
(deg)
Case

Epoch 2014 May 23.0 = JDT 2456800.5

1.126001049342
0.203702111
6.0348667
2.0368635
66.2936367
42.7542947
a

1.126001049356
0.203702117
6.0348660
2.0368642
66.2936374
42.7542940
b

Epoch 2175 Sep. 18.0 = JDT 2515722.5

1.136
0.2038
6.33
359.7
65.3
308.6
a

1.137
0.2043
6.30
359.7
65.2
308.7
b

Epoch 2135 Jan. 1.0 = JDT 2500852.5

0.19726899
6.2436831
0.33356
70.70832
83.980
b

Next, we compute with the multiple solution method 201 VAs, that is, 100 orbits on both sides of the LOV around these orbits with from 0.01 to 0.0000001, propagate them 20 days after the moment of the possible impact in 2175, and search for VAs which hit the earth. Hence, we got a sample of orbits which hit the earth. They are so-called impact orbits. Table 6 lists orbital elements of one of these orbits for the initial epoch and 7 days before impact in 2175, that is, for MJD 56600 and MJD 115722, for the case a and b. Additionally, orbital elements with their errors are presented computed for the epoch before CA with the earth in 2135.

We can see in Table 6 that errors of propagated orbital elements are great, mainly in 2175, after CA with the earth in 2135. We observe strong scattering of nearby orbits and so the subsequent impact hazard after 2135 can only be explored through statistical means, as Chesley et al. [6] state.

Table 7 shows CAs with the earth of impact orbits in Table 6 before the possible impact in 2175.

Table 7: The asteroid (101955) Bennu. CAs of impact orbits with the earth before possible impact in 2175.

Date of CA
Distance to the earth

(au)

Case a

2037/02/11.56138
0.09871340

2043/02/09.76131
0.09662668

2054/09/30.04165
0.03929920

2060/09/23.02529
0.00500806

2068/02/15.13979
0.07047181

2080/09/22.02693
0.01552646

2087/10/02.02561
0.05252972

2102/02/17.42098
0.08023926

2135/09/25.45803
0.00241115

2152/09/20.61574
0.03825595

2159/02/05.93591
0.08727171

Case b

2037/02/11.56136
0.09871338

2043/02/09.76131
0.09662668

2054/09/30.04165
0.03929922

2060/09/23.02529
0.00500808

2068/02/15.13981
0.07047169

2080/09/22.02779
0.01551729

2087/10/01.99090
0.05237801

2102/02/17.42176
0.08026213

2135/09/25.47238
0.00253219

2142/02/07.21875
0.09924421

It is difficult to compute these impact orbits in 2175 because of their CAs with the earth in 2060 at 0.0050 au and in 2135 at 0.0025 au. These deep CAs can be called “key holes.”

Table 8 in Chesley et al. [6] lists Bennu deterministic earth approaches closer than 0.05 AU according to JPL solution 76. One of them, CA in 2135 Sep 25.40942, occurs between 0.000819 au and 0.003549 au, with nominal value 0.002009 au. It is the nominally sublunar distance encounter in 2135. This deep close approach leads to chaotic motion of the asteroid Bennu after 2135, as Chesley et al. [6] state. Our results are inside this region of CA in 2135. Our computations show propagation of error of orbital elements after 2135 as is presented in Table 6.

5. Path of Risk

Next, using previously computed impact orbits for 2175, we can draw the paths of risk for that year.

The impact orbits can be used to compute the path of risk which is a locus of possible positions for an impact event on the earth’s surface [14]. This can be done using various existing software. We used the mercury Integrator Package v.6.1 software by Chambers [15] and the DE405/WAW Solar system ephemerides [16]. We used the Mercury package with RADAU algorithm and with accuracy parameter = .

As mentioned in the previous section, impact orbits from Table 5 (cases a and b) are cloned using the multiple solution method of Milani by giving in our case a of 0.02 and using 100 orbits on both sides of the LOV. Using the Mercury integrator and as starting orbital elements all 201 cloned orbits (VAs), for the epoch 7 days before impact in 2175, that is, for the epoch 2175 Sep. 18.0, we computed the minimum distances between the asteroid (101955) Bennu and the Earth around possible day of impact in 2175.

The results for this swarm of clones hitting the earth are almost the same as those computed directly from the OrbFit software and the JPL DE406 Solar system model. We then draw a narrow corridor containing the possible impacts.

The path of risk of the possible impacts in the year 2175 for cases a and b is shown in Figure 1.

Figure 1: The paths of risk where the asteroid (101955) Bennu could impact in the year 2175. The central path is connected with the impact of smaller probability, that is, case a, and represents only the central part of this path. Two bottom paths are connected with the possible impact of greater probability, that is, case b, without the connecting circumpolar part.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

Figure 2: Time evolution of the mean orbital elements and the mean position of the ascending and descending nodes of the 21 VAs of the asteroid (101955) Bennu 5000 yr backwards and 1000 yr forward from the starting epoch 2014-May-23.0. Positions of the asteroid (101955) Bennu are marked every 2000 days.

There is a cloud of orbit solutions, most of which do not impact the earth. Only a small part of that cloud gives impact solutions. Most of the clones spread out along track.

Note that path of risk is computed for 3 uncertainty. Both ends of the path of risk have lower impact probability than the central places.

Impact a goes through the Gulf of Mexico and is connected with the impact of smaller probability and represents only the central part of this path.

Two bottom path risks are connected with the possible impact of greater probability—case b. Path of risk goes close to the Geographic South Pole. We draw only both ends of this path of risk without the connecting circumpolar part.

6. Time Evolution

Figure 2 presents time evolution of the mean orbital elements and the mean position of the ascending and descending nodes orbital nodes of the 21 VAs of the asteroid (101955) Bennu 5000 yr backwards and 1000 yr forward from the osculating epoch 2014 May 23.0 = JDT 2456800.5. The starting orbital elements are computed with the OrbFit software using the multiple solution method [9, 10] for 1  uncertainty region.

We used the JPL DE406 Solar system model and nongravitational parameter  au/d². The weighing and selection method of the NEODyS was used. We also used the new error model “cbm10” based upon Chesley et al. [7] and 25 additional perturbing massive asteroids according to Chesley et al. [6].

Positions of the asteroid (101955) Bennu are computed every 2000 days. Next, we computed mean value of a given orbital element computed from all 21 orbital elements of each VA for the same date. Hence, we got 1091 mean values of each orbital element for each 2000 days’ step in the time-span (−5000, +1000) yr from the starting epoch, t0 = 2014-May-23.0.

6.1. Time Evolution of Semimajor Axis, Eccentricity, and Inclination

Figure 2 shows that the value of mean semimajor axis changes from about 1.146 au 5000 yr ago to about 1.123 au 1000 yr in the future. Hence mean orbital period of the asteroid (101955) changes from about 448 days to about 434 days, respectively, that is, mean orbital period decreases of about 14 days, that is, of about 3.1%. The value of mean eccentricity changes its value from about 0.226 5000 yr in the past and stays around present value of about 0.202 to the end of 1000 yr integration in the future. The value of mean orbital inclination is in the range 4.55° 5000 yr ago to 5.88° 1000 yr in the future. Generally, the orbit of the asteroid Bennu is stable in the time-span (−5000, +1000) yr from the starting epoch, t0 = 2014-May-23.0.

6.2. Time Evolution of the Orbital Nodes

In Figure 2 we can see that the mean longitude of the ascending node, , changes slowly its value from about 77° 5000 yr ago to 349° 1000 yr in the future, that is, changing its value with velocity of about −52.8′′/yr. The NEODyS for the asteroid (101955) Bennu (http://newton.dm.unipi.it/neodys/index.php?pc=1.1.6&n=101955) gives value , the frequency of time-variation of the longitude of the ascending node, equal to −50.674′′/yr. We computed this value from the mean longitude of the ascending node of all 21 VAs. It is interesting that our computed value of the time-variation of the longitude of the ascending node is close to the value of computed by the NEODyS Team.

Mean argument of perihelion, , changes in this time range its value from about 356° through 0° 4750 yr ago to 79° 1000 yr in the future, that is, changing its value between 4750 yr backwards and 1000 yr in the future with velocity of about +49.46′′/yr.

Mean value of the longitude of perihelion () changes its value from about 73° to 68°. The mean frequency of circulation of in the range 5000 backwards and 1000 yr in the future is about −3.0′′/year, where the NEODyS value of ′′/year.

Both parameters and are connected with so-called proper elements and proper frequencies. They are changing very slowly and are almost constant over long intervals of time. The detailed discussion is presented in Section 9.

6.3. Time Evolution of the Positions of the Nodes

Figure 2 shows that mean position of the ascending node of the asteroid (101955) Bennu changes its value from about 0.89 au 5000 yr ago to 1.4 au 1000 yr in the future. According to the NEODyS the ascending node-earth separation is now of about −0.00561 au.

Position of the descending node of the asteroid (101955) Bennu starts from about 1.4 au 5000 yr in the past to about 1.1 au in the next 1000 yr. According to the NEODyS the descending node-Earth separation is now of about +0.17756 au. In both positions of the nodes sign + denotes position of node outside the orbit of the earth − inside the orbit of the earth.

From both values of the position of the nodes it is clear that possible impacts of the asteroid (101955) Bennu with the earth can occur close after present time, according to the run of the position of the ascending node, and maybe after 1000 yr in the future according to the position of the descending node. These prognostic dates of possible impact are based on the mean value of both positions of nodes not for individual possible impactors (VAs).

6.4. Mean Motion Resonances

Figure 3 presents positions of the 1201 VAs of the asteroid (101955) Bennu 5000 yr backwards and 1000 yr forward.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

Figure 3: Orbital elements of the 1201 VAs of the asteroid (101955) Bennu 5000 yr backwards and 1000 yr forward from the osculating epoch 2014-May-23.0.

In Figure 2 we can see that the mean value of the semimajor axis is almost the same for a long time in two ranges. Between 500 and 5000 yr ago semimajor axis of Bennu has a value of about 1.146 au and between 500 and 1000 yr in the future has a value of about 1.123 au. At these intervals of time, the asteroid (101955) Bennu is locked in two main 5 : 6 mean motion resonance (MMR) with the earth and 1 : 2 with Venus, respectively.

7. Time Evolution 5000 yr Backwards and 1000 yr Forward

Additional computation of the time evolution of the orbit of the asteroid (101955) Bennu was made using the OrbFit software and the JPL DE406, 25 perturbing massive asteroids, biased error model based on Chesley et al. [7] and follows the same method in the weighting and selection of observations that is being used by the NEODyS site [8]. The computations were made with the nongravitational parameter  au/d² according to Chesley et al. [6].

Masses of 25 perturbing asteroids were taken from Farnocchia et al. [3]. We computed starting positions of these asteroids and their perturbations using the ASTDyS base of the initial orbital elements of asteroids (http://hamilton.dm.unipi.it/astdys/) and the OrbFit software.

We propagated 1201 VAs computed with the multiple solution method of Milani et al. [9, 10] for the 3 uncertainty around the nominal orbit of the asteroid (101955) Bennu 5000 yr backwards and 1000 yr forward.

Figure 3 presents value of orbital elements and position of the orbital nodes, of 1201 VAs of the asteroid (101955) Bennu after 5000 yr backwards evolution (triangles) and 1000 yr in the future (crosses). Horizontal lines denote their initial values.

Semimajor axis: 5000 yr ago semimajor axes of all 1201 VAs are slightly above their initial value. Instead, after 1000 yr in the future, they are mainly below the initial value of semimajor axis with many values of semimajor axis which are greater or smaller the initial value. This is almost the same results as in the case of 5000 yr backwards and 1000 yr forwards integration as presented in Figure 2 for time evolution of mean values of semimajor axes.

Generally, values of all orbital elements and positions of the ascending and descending nodes of all VAs of Bennu have greater dispersion 1000 yr in the future than in backwards because of the CAs with the earth in that time—see Section 8. It is interesting that all VAs are inside in our Solar system model during 6000 yr integration.

It is worth noting that Figures 2 and 3 present similar behavior of the orbital elements and the position of the nodes of the asteroid (101955) Bennu despite using different value of the uncertainty, 1 and 3, respectively, and mean and temporary values of the orbital elements of VAs.

8. Close Approaches to the Planets

Next we searched for the CAs computed earlier in the previous section during studying time evolution of the orbital elements of the Bennu.

Computed CAs of 1201 VAs of Bennu with the earth during 1000 yr forward integration are presented in Figure 4. In the left panel we present distances of Bennu to the earth in the next 200 yr. It is visible that till about 2140 there are only several CAs. They are connected with every VA. After 2140 there are many chaotically placed CAs in the range from about 0 to 0.1 au. In the right panel we detected that the bottom bound of CAs goes higher and distances between the earth and Bennu are greater; that is, possible impacts with the earth are almost excluded. The most possible impacts may occur in the time range of about 2140–2250.

(a)

(b)

Figure 4: CAs with the earth of the 1201 VAs of the asteroid (101955) Bennu 200 yr and 1000 yr forward from the osculating epoch 2014-May-23.0.

9. Secular Orbital Resonances of the Asteroid (101955) Bennu

From the NEODyS data set (http://newton.dm.unipi.it/neodys/index.php?pc=5), we can take the computed values of proper elements for the NEOs, that is, the precession rate of perihelion (the frequency of perihelion), that is equal to the sum of changes in the argument of perihelion, , and in the longitude of the ascending node, , of the orbit of the asteroid, expressed in arc seconds per year and the precession rate of the ascending node (the frequency of the ascending node), that is equal to the change in the longitude of the ascending node, , of the orbit of the asteroid, expressed in arc seconds per year.

The NEODyS (http://newton.dm.unipi.it/neodys/index.php?pc=1.1.6&n=101955) provides the computed values of and for the orbit of asteroid (101955) Bennu (′′/yr, ′′/yr).

A particular secular resonance occurs when the frequency of time-variation of the longitude of perihelion, , or the longitude of the ascending node, , of an asteroid becomes close to the orbital frequency of a planet. According to Milani and KnežEvić [17] a resonance is a resonant combination of the form , where is integer, and , , and denote the frequency of the longitude of perihelion and the frequency of the longitude of the ascending node of the asteroid and of Saturn (′′/yr, ′′/yr).

An asteroid in secular resonance with a planet will precess at the same rate as the planet. Over long times, that is, a million years, or so, a secular resonance can change the eccentricity and inclination of the asteroid [18].

Figure 5 presents the plane of proper frequency of the longitude of perihelion versus proper frequency of longitude of the ascending node for 10 214 near-earth asteroids computed by the NEODyS as for January 12, 2014. Lines indicate the position of the 1, 2, and 3 secular resonances. The frequencies and shown in Figure 5 are calculated according to the synthetic theory developed in Knežević and Milani [19] and are based on Nobili et al. [18] and Knežević and Milani [20].

Figure 5: The plane of proper frequency of the longitude of perihelion versus proper frequency of longitude of the ascending node for 10 214 near-earth asteroids. Lines indicate the positions of the secular resonances. Position of the asteroid (101955) Bennu is marked by a star.

Figure 5 shows that almost all the NEOs are inside the 1, 2, and 3 resonances where lines labeled with 1, 2, and 3 are positions of exact resonances.

Resonance 3 bounds the population of almost all known NEOs. Two proper elements of the asteroid (101955) Bennu, and , are placed near the center of all known NEAs. We can see that there are many NEAs which have greater values of , between (−200, +200)′′/yr and , in the range (0, −1000)′′/yr.

It is interesting that our orbital evolution of mean eccentricity and mean argument of perihelion of the asteroid (101955) Bennu in Figure 2 is similar to its evolution figure of two secular elements, eccentricity and argument of perihelion, as presented in the NEODyS site (http://newton.dm.unipi.it/neodys/index.php?pc=1.1.6.1&n=101955). Eccentricity changes its value between 0.20 and 0.22, and argument of perihelion is in the range 50 to 200°.

10. The Influence of the Different JPL Planetary and Lunar Ephemerides on the CAs of the Asteroid (101955) Bennu with the Earth

Table 8 lists date (TDB) and values of distances of the asteroid (101955) Bennu to the earth using different JPL planetary and lunar ephemerides.

Table 8: CAs of the asteroid (101955) Bennu with the earth using different JPL ephemerides.

Date (TDB)
Distance (au)
JPL ephemerides

2037/02/11.56138
0.09871325
DE403

2037/02/11.56137
0.09871324
DE405

2037/02/11.56139
0.09871324
DE414

2037/02/11.56139
0.09871324
DE423

2043/02/09.76129
0.09662666
DE403

2043/02/09.76129
0.09662666
DE405

2043/02/09.76129
0.09662666
DE414

2043/02/09.76129
0.09662666
DE423

2054/09/30.04168
0.03929934
DE403

2054/09/30.04168
0.03929934
DE405

2054/09/30.04168
0.03929933
DE414

2054/09/30.04168
0.03929933
DE423

2060/09/23.02527
0.00500821
DE403

2060/09/23.02527
0.00500821
DE405

2060/09/23.02527
0.00500820
DE414

2060/09/23.02528
0.00500820
DE423

2068/02/15.13995
0.07047081
DE403

2068/02/15.13995
0.07047080
DE405

2068/02/15.13995
0.07047082
DE414

2068/02/15.13992
0.07047086
DE423

2080/09/22.03401
0.01545136
DE403

2080/09/22.03405
0.01545092
DE405

2080/09/22.03394
0.01545204
DE414

2080/09/22.03363
0.01545533
DE423

2080/09/(22.03132 ÷ 22.03663)  (0.01542356 ÷ 0.01547990) au

2087/10/01.74349
0.05128111
DE403

2087/10/01.74185
0.05127374
DE405

2087/10/01.74603
0.05129252
DE414

2087/10/01.75828
0.05134743
DE423

2087/10/(01.64017 ÷ 01.85015)  (0.05081501 ÷ 0.05175733) au

2102/02/17.42477
0.08038730
DE403

2102/02/17.42478
0.08038790
DE405

2102/02/17.42477
0.0803863
DE414

2102/02/17.42471
0.08038179
DE423

2102/02/(17.42405 ÷ 17.42481)  (0.08034178 ÷ 0.08041831) au

2135/09/25.58891
0.00353782
DE403

2135/09/25.58977
0.00354536
DE405

2135/09/25.58759
0.00352619
DE414

2135/09/25.58122
0.00347058
DE423

2135/09/(25.53554 ÷ 25.64539)  (0.00307325 ÷ 0.00403380) au

2141/09/26.48864
0.01531402
DE403

2141/09/26.23228
0.01343206
DE405

2141/09/26.89761
0.01823168
DE414

2141/09/29.12483
0.03233813
DE423

2141/09/(16.49276 ÷ 53.82944)  (0.00492126 ÷ 0.09701912) au

All computations were made using the OrbFit software for the nominal orbit of the asteroid (101955) Bennu with the nongravitational parameter  au/d² and 25 additional perturbing massive asteroids according to Chesley et al. [6]. The weighing and selection method of the NEODyS was used. We also used the new error model “cbm10” based upon Chesley et al. [7].

Only from 2080 there are significant differences in the values of the distance of the asteroid (101955) with the earth caused with the used different JPL ephemerides.

In 2080 differences in distances to the earth computed from different JPL ephemerides are below 0.1 earth radii, in 2087, 1.7; 2102, 0.15; 2135, 1.3; and 2141, 440 earth radii. Appropriate differences in time of CAs are 2080, 36 s; 2087, 25 m; 2102, 6 s; 2135, 12 m; and 2146, 2.9 d.

The differences are connected mainly with the deep CA in 2080. Table 8 shows that the results of computations of CAs of the nominal orbit depend on the used JPL ephemerides.

Next we computed 10 VAs of the asteroid Bennu on each side of the Line of Variation (LOV) using the multiple solution method and propagated them to 2200. We used the JPL DE406 because it is the main planetary and lunar ephemerides used in our work in computing impact solutions of the asteroid (101955) Bennu. The results are presented in Table 8 starting from CAs in 2080. Each suitable line lists 1 uncertainty in time and distance of CAs computed for all 21 VAs.

It is visible that all differences in computed distances and moments of CAs of the asteroid Bennu and the earth using different JPL ephemerides are inside 1 uncertainty of the multiple solution of the asteroid Bennu; that is, error of propagation of the orbital elements of the asteroid Bennu has a main role in the motion of asteroid.

11. Summary

Using the freely available OrbFit Software Package and the nongravitational parameter in the motion of the asteroid (101955) Bennu taken from Chesley et al. [6] we computed possible impact solutions using different JPL planetary and lunar ephemerides and different number of additional massive perturbed asteroids. The possible impact paths of risk for 2175 are presented.

We computed that the most possible impacts may occur in the time range of about 2140–2250.

Additionally, we computed possible impact solutions using the normal places method of the selection of Bennu’s astrometric observations.

Moreover, we computed time evolution of the mean orbital elements and the position of the orbital nodes of Bennu 5000 yr in the backwards and 1000 yr in the future using the Yarkovsky effects. It is interesting that the asteroid (101955) Bennu is temporarily locked in the 5 : 6 mean motion resonance (MMR) with the earth and 1 : 2 with Venus. Also the asteroid Bennu is located in the phase-space of the 1, 2, and 3 secular resonances.

We also computed the influence of the JPL planetary and lunar ephemerides DE403, DE405, DE406, DE414, and DE423 on the close approaches of the asteroid (101955) Bennu with the earth.

We also find that all differences in computed distances of the asteroid Bennu to the earth using different JPL planetary and lunar ephemerides are inside 1  uncertainty of the multiple solution of the asteroid Bennu. Hence error of propagation of the orbital solution is greater than the differences connected with using different JPL planetary and lunar ephemerides.

Additional observations, optical and radar are needed to determine values of the non-gravitational parameters connected with the Yarkovsky/YORP effects.

Conflict of Interests

The author declares that there is no conflict of interests regarding the publication of this paper.

Acknowledgments

The author would like to thank Grzegorz Sitarski from the Space Research Center, the Polish Academy of Sciences in Warsaw, for many helpful discussions. The author would also thank Davide Farnocchia for his useful comments. The author’s thanks go to the researchers from the OrbFit Consortium for their free software and source code. The author also thank an anonymous reviewer for important suggestions that helped to significantly improve the final version of the paper.

References

A. Milani, S. R. Chesley, M. E. Sansaturio, F. Bernardi, G. B. Valsecchi, and O. Arratia, “Long term impact risk for (101955) 1999 RQ3 6,” Icarus, vol. 203, no. 2, pp. 460–471, 2009.

D. Vokrouhlický, A. Milani, and S. R. Chesley, “Yarkovsky effect on small near-earth asteroids: mathematical formulation and examples,” Icarus, vol. 148, no. 1, pp. 118–138, 2000.

D. Farnocchia, S. R. Chesley, P. W. Chodas et al., “Yarkovsky-driven impact risk analysis for asteroid (99942) Apophis,” Icarus, vol. 224, no. 1, pp. 192–200, 2013.

D. Vokrouhlický, “Direct solar radiation pressure on the orbits of small near-earth asteroids: observable effects?” Astronomy and Astrophysics, vol. 362, no. 2, pp. 746–755, 2000.

D. Farnocchia and S. R. Chesley, “Assessment of the 2880 impact threat from Asteroid (29075) 1950 DA,” Icarus, vol. 229, pp. 321–327, 2014.

S. R. Chesley, D. Farnocchia, M. C. Nolan et al., “Orbit and bulk density of the OSIRIS-REx target Asteroid (101955) Bennu,” Icarus, vol. 235, pp. 5–22, 2014.

S. R. Chesley, J. Baer, and D. G. Monet, “Treatment of star catalog biases in asteroid astrometric observations,” Icarus, vol. 210, no. 1, pp. 158–181, 2010.

I. Wlodarczyk, K. Cernis, R. P. Boyle, and V. Laugalys, “Discovery and dynamical characterization of the Amor-class asteroid 2012 XH₁₆,” Monthly Notices of the Royal Astronomical Society, vol. 438, no. 3, Article ID stt2382, pp. 2621–2633, 2014.

A. Milani, S. R. Chesley, M. E. Sansaturio, G. Tommei, and G. B. Valsecchi, “Nonlinear impact monitoring: line of variation searches for impactors,” Icarus, vol. 173, no. 2, pp. 362–384, 2005.

A. Milani, M. E. Sansaturio, G. Tommei, O. Arratia, and S. R. Chesley, “Multiple solutions for asteroid orbits: computational procedure and applications,” Astronomy and Astrophysics, vol. 431, no. 2, pp. 729–746, 2005.

W. M. Folkner, P. Charlot, M. H. Finger et al., “Determination of the extragalactic-planetary frame tie from joint analysis of radio interferometric and lunar laser ranging measurements,” A&A, pp. 287–279, 1994.

A. S. Konopliv, C. F. Yoder, E. M. Standish, D.-N. Yuan, and W. L. Sjogren, “A global solution for the Mars static and seasonal gravity, Mars orientation, Phobos and Deimos masses, and Mars ephemeris,” Icarus, vol. 182, no. 1, pp. 23–50, 2006.

S. R. Chesley, M. C. Nolan, D. Farnocchia et al., “The trajectory dynamics of near-earth asteroid 101955 (1999 RQ36),” in Proceedings of the Conference on Asteroids, Comets, Meteors (ACM '12), LPI Contribution 1667, 6470, ACM, Niigata, Japan, May 2012.

G. Sitarski, “On the small asteroid 1994 GV,” Acta Astronomica, vol. 50, no. 3, pp. 417–419, 2000.

J. E. Chambers, “A hybrid symplectic integrator that permits close encounters between massive bodies,” Monthly Notices of the Royal Astronomical Society, vol. 304, no. 4, pp. 793–799, 1999.

G. Sitarski, “Warsaw ephemeris of the solar system: DE405/WAW,” Acta Astronomica, vol. 52, no. 4, pp. 471–486, 2002.

A. Milani and Z. KnežEvić, “Asteroid proper elements and the dynamical structure of the asteroid main belt,” Icarus, vol. 107, no. 2, pp. 219–254, 1994.

A. M. Nobili, A. Milani, and M. Carpino, “Fundamental frequencies and small divisors in the orbits of the outer planets,” Astronomy and Astrophysics, vol. 210, no. 1-2, pp. 313–336, 1989.

Z. Knežević and A. Milani, “Synthetic proper elements for outer main belt asteroids,” Celestial Mechanics and Dynamical Astronomy, vol. 78, no. 1–4, pp. 17–46, 2000.

Z. Knežević and A. Milani, “Proper element catalogs and asteroid families,” Astronomy and Astrophysics, vol. 403, no. 3, pp. 1165–1173, 2003.