## 1 Introduction  
One of the main problems in celestial mechanics is the two-body problem, especially in its elliptical case. This problem can be formulated as follows: given two masses \(m_{1}\), \(m_{2}\) called primary and secondary, respectively, establish the position of the secondary relative to the primary, denoted by \(\vec{r}\), as a function of time. The second-order differential equations define the relative motion of the secondary with respect to the primary:  
\[\frac{d^{2}\vec{r}}{dt^{2}}=-\mu\frac{\vec{r}}{r^{3}} \tag{1}\]  
where \(\vec{r}\) is the vector radius of the secondary and \(\mu=G(m_{1}+m_{2})\), where \(G\) is the gravitational constant and \(m_{1}\), \(m_{2}\) the masses of the primary and the secondary, respectively. The two-body problem in celestial mechanics is a well-known classical integrable problem whose analytical solution is well-known. It has been studied by numerous authors, e.g., [1; 2; 3; 4; 5]. This problem is also appropriate for testing numerical methods because we can compare the numerical and the analytical solutions. The analytical solution is described by the set III of elements of Brouwer and Clemence \((a,e,i,\Omega,\omega,M)\)[1], where \(a\) is the semimajor axis of the ellipse, \(e\) the eccentricity, \(i\) the inclination, \(\Omega\) the argument of ascending node, \(\omega\) the argument of the apoapsis, and \(M\) the mean anomaly for the epoch. The performance of the numerical methods is good. However, the distribution of points provided by the natural time is contrary to the dynamics by accumulating more points in the apoapsis region than in the periapsis region if constant time intervals are taken as base points of integration when the opposite would be desirable. One way to achieve this goal is to find a temporal variable for which the distribution of points at constant intervals in the new variables is more significant in the areas of higher speed of the secondary and those with greater curvature. This is particularly interesting in the case of very high eccentricities, like most periodic comets, in whose study numerical integrators can and do present problems. There are several ways to solve this problem. In this paper, we use the analytical regularization of step size to reparametrize the orbit to obtain a more appropriate point distribution in the orbit. This method can be combined with variable step-size integrators, symplectic integrators, and other techniques. The reparametrization technique is a very interesting method to solve several problems. A review of these methods can be seen in [6]. In 1912, [PERSON] [7] introduced a change in the temporal variable \(\tau\) using the transformation \(r\,d\tau=dt\), a temporal reparametrization of the motion known as the Sundman transformation. In this transformation, we can change the time by the mean anomaly; multiplying by the mean motion \(n\) and including in \(\tau\) a normalization factor \(\frac{1}{a}\), the Sundman transformation can be rewritten as \(dg=\frac{\tau}{a}dM\), where \(g\) is the eccentric anomaly. To solve this problem, several authors have tried to use a time transformation to obtain a point distribution more according to the dynamics. * [PERSON] [8] pioneered in this field, introducing the relative anomalies. * [PERSON], in 1912 [7], introduced a new temporal variable \(s\) related to \(t\) by means of \(r\,ds=dt\). * [PERSON] [9] introduced a new anomaly \(\tau\) related by \(r^{\frac{3}{2}}d\tau=dt\). * [PERSON] and [PERSON] [10; 11] extended this transformation to \(\Psi_{\alpha}\), defined by \(r^{\alpha}d\Psi_{\alpha}=dt\). * [PERSON] [12] introduced the regularized length of arc \(s^{*}\) by \(vds^{*}=dt\), where \(v\) is the velocity of the secondary. * [PERSON] and [PERSON] [13] introduced the elliptic anomaly \(w\) as \(w=\frac{\pi\varepsilon}{2K(\varepsilon)}-\frac{\pi}{2}\), where \(\,\mathrm{am}\,z=g+\frac{\pi}{2}\).