Combining the 11 resonant returns of Fig. 1 and the 14 nonresonant ones of Fig. 2, our theory predicts 25 close approach solutions; for each of these we could perform a detailed close approach analysis to determine the minimum distance possible. However, this is not necessary because the minimum distance is essentially the local MOID near the relevant node; this is also not sufficient to identify all the possible returns, because a secondary return from a previous one is possible, and so on. For this reason we have devised a global method to find returns.
We started from the same catalog of alternate orbital solutions used for the figures. Each solution was propagated forward from the 2027 encounter, recording the position and the nodal distances every time the Earth is passing at the nodes. We determine if there was a crossing of the relevant node (near that time) by the changes of sign of the z coordinate in an ecliptic reference frame. We interpolate between these adjacent solutions to find the value corresponding to the node crossing at the time when the Earth is there. We similarly obtain the minimum distance between the orbit of the Earth and of the asteroid around the relevant node. By a continuity argument, if z changes sign between two solutions at and , there is an intermediate value of for which z is zero, that is at least one solution along the LOV always exists that passes at a distance from the Earth as low as the local MOID (even slightly less, due to gravitational focusing).
This argument cannot be applied for values of z too large, otherwise the two consecutive solutions could be out of phase by more than one period. Thus the limit of the method is the stretching , which is the ratio between the distance in physical space of two orbit solutions at some time and the distance of the corresponding values of , which parametrises the LOV. For , as in our solutions catalog, AU would result in two consecutive orbits being out of phase by 1 revolution; we can reliably detect a close approach only up to . (P. Chodas, private communication, has found another return in August 2039 which has escaped our search because it has .) After a very close approach such values of do occur, and even more after a sequence of close approaches. For this reason we have densified our sampling of the LOV in the region of high stretching around the solution with the closest approach in 2027, namely for , by computing another alternate orbits. With , even returns with AU can be detected.

Table I presents all the returns up to August 2040 that we have found with this method, using both the catalog and the denser catalog. The stretching in the Table is not , computed with distances in the 3D space, but its projection upon the MTP, which is in a fixed ratio to . That is, we use the product of the time difference in the node crossing and the relative encounter velocity divided by . allows one to compute the size of the interval along the LOV, in units, where approaches within a given distance occur. Given a probability density function on the LOV, the probability of such an event can be determined. But, there is no such thing as a unique probability of an event involving an orbit obtained by a least squares fit: it depends upon assumptions on the statistical distribution of the observational errors. In the Table we have used a uniform probability density along the LOV for , to estimate the probability of an encounter within the mean distance of the Moon. Note that the lower the stretching, the higher the probability of an encounter within a given distance; thus shallow encounters can be more effective in generating likely returns than the deep ones.
Each of the 25 returns predicted by our theory appear in the Table, with . 6 solutions not predicted by the figures appear; they can all be interpreted as secondary returns. Both the 5/3 and 2/1 returns become possible after the 2034 encounter. Among these secondary returns there is one in August 2039 for which the interpolated MOID is less than the radius of the Earth. Since the stretching is extreme, we have checked by performing close approach analysis: a collision solution does exist. But appears as divisor in the formula for the probability, so the probability for this impact is of the order of 10^{9}. If the probability of an impact by an undiscovered 1 km asteroid is of the order of 10^{5} per year [Chapman & Morrison 1994], the probability of impact in 2039 is less than the probability of being hit by an unknown asteroid of this size within the next few hours. In any case the asteroid orbit will soon be refined by further observations and this possible solution may be ruled out.
The stretching coefficients used here are related to the dimensionless stretching used in the computations of the Lyapounov characteristic exponents: they differ only by a constant factor. Thus the data in the Table indicate the level of chaos of each return orbit. The cascade of successive returns could be described by a symbolic dynamics, as in other chaotic celestial mechanics problems [Zare & Chesley 1998].