## 3 Numerical Experiments with Real Debris Orbital Data  
Now we show the results obtained from the numerical experiments. The data for the space debris used in the simulations is taken from the \"stuffin.space\" website, a platform that receives daily updates with orbital data from \"Space-Track.org\". The orbital data of the debris considered in the dynamics for each result presented here is presented in the legend of each figure. The disturbing potential given by Equation (5) is replaced into the system of the Lagrange planetary equations. This system of equations is numerically integrated using Maple software. For the numerical integration of the dynamical system, we use Maple's dsolve routine with the options \"numeric\" and \"method = rkf45\", which represents an integration using the Runge-Kutta-Fehlberg method. In the mentioned website, [[http://stuffin.space/](http://stuffin.space/)]([http://stuffin.space/](http://stuffin.space/)) (accessed on 7 September 2022), we obtain the initial conditions of real debris, called BLOCK DM-SL-DEB (2019-095E), which we will call Case 1. Here, we investigate the mitigation process by removing the debris into a graveyard orbit. Thus, we consider raising the debris to a new orbit 1000 km above the geostationary ring. To perform this maneuver, we assume that we have a solar sail coupled to the satellite with different values of \(A/m\), namely \(A/m=[1,2,3,4,5,6]\) m\({}^{2}\)/kg and natural disturbances such as SRP, the flattening of Earth (\(I_{2}\)), and the gravitational perturbation of the Moon and the Sun. Figure 1a shows that the eccentricity increases proportionally to the increase in the \(A/m\) ratio in a short period of time. That is, the greater the value of A/m, the greater the variation in the eccentricity of the debris in a short period of time. We see that, in less than one year, the eccentricity reaches a considerable variation for all cases analyzed for the values of A/m. Figure 1b shows the increase in inclination also for different values of \(A/m\). Note that the time to increase the variation of the inclination is much longer than the eccentricity case, but there is a large inclination variation anyway. This can be used to change the orbital plane of the debris. Note that we use zero for the argument of perigee and longitude of the ascending node, because the site does not provide the values of these orbital elements. Figure 1c shows the apogee position. The horizontal lines in yellow and red represent the protected GEO region. The black horizontal line represents the altitude of the chosen graveyard orbit, which in this case is 1000 km above \(r_{GEO}\) (35,786 km). Note that, with just a solar sail with \(A/m=2\) m\({}^{2}\)/kg, the debris is directed to the desired orbit. However, note also that, if no other maneuvers are performed, the debris returns to its initial orbit. Then, when it reaches the desired altitude, a new maneuver must be performed so that the debris remains in the orbit 1000 km above the protected region. This maneuver, which must be performed when the debris reaches the desired altitude, is not the focus of this work. Here, our main objective is to choose the \(A/m\) parameter for the debris to reach the desired altitude. Furthermore, another objective is to choose convenient values for the argument of perigee (\(g\)) and the longitude of the ascending node (\(h\)) that amplify the growth of the eccentricity. To verify the contributions of these orbital elements, we plotted dynamical (color) maps (see Figure 2a) varying the initial conditions for the argument of perigee and longitude of the ascending node from 0 to 360 degrees with steps of 5 degrees. In the color scale we show the amplitude of the eccentricity (\(Amp(e)=e_{max}-e_{min}\)). As shown in Figure 2a, the amplitude of the eccentricity depends on the initial values of \((g,h)\). In this way, we can direct the debris to one of these values of \((g,h)\) that will help to amplify the orbital eccentricity of the debris.