AbstractAccurate trajectories for the orbits of the planets in our solar system depends on obtaining an accurate solution for the elliptical Kepler equation. This equation is solved in this article using the homotopy perturbation method. Several properties of the periodicity of the obtained approximate solutions are introduced through some lemmas. Numerically, our calculations demonstrated the applicability of the obtained approximate solutions for all the planets in the solar system and also in the whole domain of eccentricity and mean anomaly. In the whole domain of the mean anomaly, 0≤M≤2π, and by using the different approximate solutions, the residuals were less than 4×10−17 for e∈[0, 0.06], 4×10−9 for e∈[0.06, 0.25], 3×10−8 for e∈[0.25, 0.40], 3×10−7 for e∈[0.40, 0.50], and 10−6 for e∈[0.50, 1.0]. Also, the approximate solutions were compared with the Bessel–Fourier series solution in the literature. In addition, the approximate homotopy solutions for the eccentric anomaly are used to show the convergence and periodicity of the approximate radial distances of Mercury and Pluto for three and five periods, respectively, as confirmation for some given lemmas. It has also been shown that the present analysis can be successfully applied to the orbit of Halley’s comet with a significant eccentricity.