Mapping Long-Term Natural Orbits about Titania, a Satellite of Uranus

Close polar and circular orbits are of great interest for the exploration of natural satellites. There are still no studies in the literature investigating orbits around Titania, the largest satellite of Uranus. In this work, we present results of a set of numerical simulations carried out to obtain long-duration orbits for a probe around Titania. Through an expansion of the gravitational potential up to second order, the asymmetry of the gravitational field due to Titania's coefficient $C_{22}$, the zonal coefficient $J_2$, and the gravitational perturbation of Uranus is considered. The analysis of lifetime sensitivity due to possible errors in the values of $J_2$ and $C_{22}$ is investigated using multiple regression models. Simulations were performed for different eccentricity values, and lifetime maps were constructed. The results show that low-altitude and near-circular orbits have longer lifetimes due to the balance between the disturbance of Uranus and the gravitational coefficients of Titania. The results also show that non-zero values of the longitude of periapsis ($\omega$) and longitude of the ascending node ($\Omega$) are essential to increase the lifetime up to eight times compared to cases where $\omega= \Omega=0^\circ$. We also show that an orbit with eccentricity $10^{-3}$ is the most affected by errors in the values of $J_2$ and $C_{22}$.


Introduction
Uranus and its satellites are one of the least explored systems in our solar system. The only mission sent to visit Uranus took place in 1986 with Voyager 2. A mission to explore the Uranus system would be extremely important to understanding the formation and evolution of the solar system.
The main information about the Uranus system comes from observations made from Earth and data sent by space telescopes, but some important questions still remain unanswered. The distances of the ice giant planets are still an enigma, as some models suggest they should be closer to the Sun. Their high obliquity and low luminosity are intriguing. In addition to all these issues, information about Uranus' satellites would also be of great importance, as they have irregular surfaces and there is evidence of oceanic formation [1,2]. According to [3], NASA is working in partnership with ESA to plan possible missions to the ice giants during the period 2024-2037. These possible missions can help to understand the origin and evolution of this system. transformation between osculating elements and averaged elements in order to find initial osculating elements for frozen orbits. This transformation was done only for circular orbits. The results showed that it was necessary to transform the average elements into osculating elements in order to add the C 22 coefficient, as this term caused large variations in high-altitude orbits and can totally change the behaviour of the orbits.
Orbits that keep their eccentricity, inclination, and argument of periapsis constant over a period resemble frozen orbits. In [12] three families of periodic lunar orbits were found. The results showed that these orbits are very sensitive to the initial conditions and that these conditions have small differences with respect to the frozen orbits. Families of frozen orbits around the Moon were also found in the work of [13]. The authors eliminated the short period terms and considered the oblateness of the Moon to order 4 (J 4 ). For eccentric frozen orbits, all families had a critical inclination around 64 • . In the case of very low quasi-circular lunar orbits, several families with different critical inclinations were obtained, although they quickly disappeared due to small changes in the semi-major axis and eccentricity.
Despite causing large variations in some orbital elements, such as inclination, eccentricity, and argument of periapsis, the term C 22 has an important role in the dynamics of orbits around natural satellites. The combination of this term with the zonal term J 2 is able to attenuate the effects caused by the perturbation due to the third body acting as a "protection mechanism" to prolong/increase the lifetime of the desired orbit [14,15]. In [16] the combination of this coefficient with the J 2 term was also investigated in orbits around the Moon. The results showed that the values of J 2 and C 22 strongly affected the near-critical inclination and the libration amplitude of the argument of periapsis. This work also showed that these effects were greater when added to the contribution of the rotation rate.
In several works found in the literature, the perturbation due to the third body was considered when it comes to orbits around natural satellites. The effects of this perturbation were mainly observed in the eccentricity of the orbit around the observed satellite. These high eccentricity values achieved by the probe's orbit are related to the effect of the Kozai-Lidov mechanism. This mechanism was studied for the first time in [17,18]. This phenomenon causes the argument of periapsis to oscillate at a constant value. Whereas the argument of periapsis oscillates around a specific value, the eccentricity and inclination suffer periodic oscillations (e.g., see [19,20] and the references therein).
The structure of this work is organised as follows: in Section 2 the mathematical model adopted to analyze the system is described, as well as the methodology. In Section 3 we discuss the results obtained. In Section 4 we show the isolated effect of each perturbation through the analysis of the orbital elements of each specific case. In this same section, we also present the mechanism for obtaining long-time orbits through the evolution of the orbital elements. In Section 5 we present an estimate of the probe's lifetime when affected by possible errors in the values of Titania's gravitational coefficients. Finally, in Section 6, we present our conclusions.

Mathematical Model
The system addressed in this work consists of a central body (Titania), a space probe, and a perturber (Uranus). We assume that the perturber body is in a Keplerian orbit around Titania with a radius of 25.362 × 10 3 km and mass 8.68 × 10 25 kg, (https://ssd.jpl.nasa.gov/), accessed on April 05, 2021. The orbital parameters are given in Table 1. To make an accurate investigation of the problem, we consider the main gravity coefficients (J 2 and C 22 ) of Titania, regarding the oblateness and ellipticity of the body. These terms are the most important in the gravitational field, after the Keplerian term. Previous work [11,16,21,22] shows that these coefficients are the most relevant to the pertubative effects caused on a probe by a satellite and, in most cases, the most significant to model the irregular shape of a body. It is important to note that these two coefficients are the only coefficients available in the literature. As our focus is to analyze orbits very close to the surface of this satellite, these will be the main terms to be considered in this work. Titania has a mass of 35.27 × 10 20 kg and a radius of 788.9 km, (https://ssd.jpl.nasa.gov/), accessed on April 05, 2021. The other parameters of Titania are shown in Table 2. The equations of motion are described according to [22]: where P T and P U are the expansions of the gravitational potential up to second order for Titania and up to fourth order for Uranus, respectively; P Tx , P Ty , and P Tz are components of the P T vector and P Ux , P Uy , and P Uz the components of the P U vector; m is the mass of the space probe, M T the mass of Titania, and M U the mass of Uranus. The terms r and r U are the radius vector of the probe and Uranus, respectively. The investigation of orbits around natural satellites with a non-uniform mass distribution has been the subject of several studies in recent years. Studies [8,9,14,[22][23][24] emphasise the importance of considering the oblateness and ellipticity of the body to be orbited by the probe. However, it is understood that analyzing the influence of higher oblateness coefficients of the disturbing body can provide a much more accurate study of the orbits.
Uranus is very far from Earth, approximately 3 billion km, and the time to send signals from Earth is on the order of 2.8 h; therefore it takes approximately 5.6 h to send a command and receive a signal back. Because of this time lag, orbital maneuvering is a logistical problem. An important reason to search for long-term natural orbits in which to place a space probe is because they do not require frequent orbital maneuvers, thereby saving fuel and simplifying the orbit control. For the sake of comparison, we did a single test considering the highest order gravity coefficients J 2 and J 4 of the disturbing body. This implementation was made through an adaptation of the Mercury package [25], which was used in our numerical simulations.  Our main goal is to study circular orbits with high inclinations. Circular orbits are better for observations as eccentric orbits have a lower perigee and are more susceptible to collisions. Highly inclined orbits take advantage of the central body's rotation to observe the entire natural satellite. However, we understand that analyzing other values of eccentricities can contribute to a more complete approach to the problem, giving more choices for the mission designers. Therefore, we considered the following eccentricity values for the probe's orbit: e = 0, e = 10 −1 , e = 10 −2 , e = 10 −3 , and e = 10 −4 . In order to analyze those orbits very close to the surface of Titania, we adopted a semi-major axis interval (a) ranging from 810 km to 2500 km for each value of the eccentricity. The inclination is in the range between 75 • and 90 • .
We considered two cases regarding the argument of periapsis and the longitude of the ascending node (ω and Ω, respectively). First, we consider ω = Ω = 0 • and, after that, we choose specific values in regions with longer lifetimes, to explore best results.
All cases were simulated for 1000 days, equivalent to approximately 10 4 orbital periods for the closest orbit to Titania's surface and 10 3 revolutions for the farthest orbit. For our purposes, this time is adequate, sufficient to analyze the orbital evolution of a space probe and collect and send data for analysis. The lifetime of the orbit is defined by the time that the space probe remains in orbit around Titania without colliding with the central body. Therefore, the maps of a i versus I i show the lifetime of the probes numerically simulated up to 1000 days. Our main goal is looking for longer lived orbits around Titania under the effects of its gravitational coefficients and Uranus.
In order to analyze the influence of the gravity coefficients of Titania and the third body, we constructed life maps considering two cases. Point mass models are considered without taking into account the gravity coefficients of both bodies,with gravity coefficients for both bodies, Titania and Uranus, and including J 2 and C 22 of Titania. "Difference maps" (third body-J 2 , C 22 from Titania) are also presented. The case in which we analyzed the J 2 and J 4 coefficients showed insignificant changes in lifetimes compared to the case where these coefficients are not included. Therefore, we will adopt Uranus as a point of mass. Our purpose with these maps is to find regions where the orbits are affected by the gravity coefficients of Titania and the third body. These effects are in opposite directions, and the orbits can last longer without collision or escape. The time evolution of the orbital elements a, e, I, ω, and Ω of the orbits with a long lifetime are shown Section 4.

Analyzing the Results
The results show diagrams a × I considering: (i) only the contribution of the disturbing body (Uranus), modelled as a point of mass; (ii) the gravity coefficients of Titania (J 2 and C 22 ), and (iii) the more complete model, where both effects are considered. In the figures, the color bar represents the lifetimes of the orbits. We constructed "difference maps" by subtracting the lifetime of the orbits obtained when Uranus is considered. The lifetime of the orbits when considering both effects, the gravitational effects of Uranus, and the gravity coefficients of Titania, are presented in Section 3.
It should be noted that, in these figures showing the differences in the lifetimes, a positive sign indicates those values where only the presence of Uranus prolonged the probe's life, whereas negative values indicate that the gravity coefficients of Titania helps to naturally control the third-body perturbations Uranus, thus increasing the lifetime of the orbits. We show lifetime "difference maps" in Figure 1c. The maximum value of the semi-major axis is a = 1200 km; after that, there are no significant changes in the orbits' lifetime because the space probe is too far from Titania. Its irregular shape is unable to balance the third-body perturbations from Uranus. Figure 1c shows that the largest differences in lifetime appear for values of a = 810-1000 km and I = 75-83°(blue region), where the equilibrium between the gravity coefficients of Titania and Uranus occurs. We found increases in the lifetime on the order of 60 days, which is not negligible. In addition, those orbits are closer to Titania, which helps in the observation of this body. It is noted that the most negative numbers occur for lower values of inclinations and semi-major axes. In these regions, the third body perturbations are smaller, and the gravity coefficients dominate the dynamics and reduce the duration of the orbits.
For e = 10 −4 , the results present similar behaviour compared to e = 0. The simulations were performed using the same initial conditions used in Figure 1, except for the eccentricity, which now equals e = 10 −4 .
We found an equal maximum lifetime when the third body and the J 2 and C 22 coefficients of Titania are considered. However, when considering only Uranus as a perturber, the region with the longest orbit life ranges in the interval a = 810-900 km and I = 75-80°. In this region, the longest lifetime is approximately 450 days. When including the gravity coefficients of Titania, this region increases to a = 810-1000 km and I = 75-84°. Another fact is that the number of orbits with a lifetime between 300-350 days decreases when the gravity coefficients of Titania are considered.
The results showed that the orbits in the range a = 810-1050 km and I = 75-84°have increased by 60 and 20 days. Orbits with a = 810-900 km and I = 85-90°have their lifetime extended by approximately 40-80 days, thanks to the effect caused by the equilibrium of the perturbations.
For e = 10 −3 , shown in Figure 2, our results, given in the diagram a × I when only the gravitational effect of Uranus is considered, are similar to those observed for e = 0 and e = 10 −4 . As already noted, the terms J 2 and C 22 of Titania are more important to increasing the lifetime of the probe, confirming that the effects caused by the gravity coefficients of Titania compensate for the effects caused by Uranus' gravitational attraction, as expected.  The lifetime differences maps presented in Figure 2c have two nearby regions where the changes in the lifetimes are expressive (regions in purple, red, and orange squares). In the region with purple squares, the lifetime for those orbits only considering the third body is up to 250 days. At the red dots, the probe's lifetime can reach up to 200 days when J 2 , C 22 coefficients of Titania are considered. Figure 2b, including the gravity coefficients of Titania, presents a more extensive region with lifetimes of 350-500 days. This region is located very close to the surface of Titania (a = 810-1000 km), where the terms J 2 and C 22 cancel the effects of the third body. Figure 2a shows the numerical simulations where the gravitational effects of Uranus are considered. Orbits with a lifetime close to these values appear only on a small island, for all inclinations, and a semi-major axis closer to 1400 km.
In the regions where the probe remains around Titania for longer times, there is a balance between the perturbation caused by the third body and the perturbation due to the gravitational coefficients of Titania, causing the probe's lifetime to be extended. This balance is described in previous work [9,14], where it is shown that, when combined with J 2 and C 22 of the central body, this resultant perturbation acts as a "protection mechanism" responsible for softening the effects caused by the third body perturbation on the variation of the eccentricity.
The results presented in Figure 3 are obtained for e = 10 −2 . Note that there are no large differences between the results obtained when considering only the third body perturbation (a) and those when the J 2 and C 22 from Titania are added in the system (b). There is only a decrease in the region with lifetimes between 160 and 180 days, covering all the inclinations, when J 2 and C 22 of Titania are included. This region is located in the range a = 900-1050 km and I = 75-80°.
For e = 10 −1 , the orbital duration maps show that the results are very close when all the perturbations are considered, and the behaviour is similar to the one obtained for e = 10 −2 . The difference is that for e = 10 −1 , the probe's maximum lifetime is much shorter, 60 days, but the regions of greater and lesser duration are very similar. Previous studies [7,28] found similar results, where more eccentric orbits have shorter lifetimes, as expected, as these orbits have smaller periapsis and, therefore, collide with the central body in less time. In addition, they have larger apoapsis, which increases the effect of the third body.
Previous work [8,29,30] shows the relevance of the contribution of the argument of periapsis (ω) and longitude of the ascending node (Ω) in the duration of highly inclined orbits. In [8], the initial conditions capable of increasing the lifetime of the probe are ω = 180 • and Ω = 90 • . In [29], the best values found for these angles are ω = 140.5-148°, ω = 321-327°, and Ω = 165 • . Given the above, we analyzed the regions with longer lifetimes, shown in Figures 1-3, in a range of values of the periapsis argument and the node longitude. Within these regions, we chose the best values of a and I and built maps as a function of ω and Ω with the four values of eccentricities adopted earlier. For this analysis, we consider the perturbation due to the third body and also the gravity coefficients of Titania.
The first case to explore is for e = 10 −4 . The results provided by the simulations showed that the longest orbits occurred for a = 900-1058 km and I = 70-80°. The result of this analysis can be seen in Figure 4.
We fixed a semi-major axis value, 900 km, and varied the inclination to 70 • (Figure 4a) and 80 • (Figure 4b). We then fixed the inclination at 80 • and vary the semi-major axis to 1018 km, (Figure 4c) and 1058 km (Figure 4d). The angles ω and Ω range from 0 • to 360 • .    Figure 4a, where the considered inclination is 70 • , the longer duration orbits have a lifetime of up to 900 days. As we increase the inclination by 10 • , the lifetime increases by 100 days (Figure 4b). This is due to the action of the gravitational coefficients J 2 and C 22 , which are responsible for attenuating the effects caused by the third body for more inclined orbits.
The regions with longer duration orbits present in Figure 4c Figure 4a,b present orbits with longer lifetimes compared to the point chosen in Figure 4d. However, our main objective is to show that the adoption of these angles can increase the lifetime of the orbit even for regions with shorter orbits. In order to analyze the importance of these angles in the lifetimes, we made a new simulation with values of Ω = 155 • , Ω = 55 • for e = 10 −4 .
The insertion of the values ω and Ω in the new simulation for e = 10 −4 increase the probe's life by approximately 78%, a very large value. We extended the numerical integration for the two cases analyzed: third body and J 2 and C 22 from Titania. However, the new results show distinct regions with more extended life for each perturbation considered. In the case of the gravitational effects of Uranus, a small island of orbits with a lifetime of 600-800 days appears with a semi-major axis ranging from a = 810 km to 900 km and inclinations from 75 • to 78 • . If we include the gravity coefficients of Titania, the maximum lifetime is approximately 600 days. It is located in the region a = 810-1000 km and can be seen for all values of I, increasing from I = 78°to 90°. Orbits with a lifetime of 800 days are found for a value of a close to 900 km and I = 78-90°.
The differences indicate that after the inclusion of the angles ω and Ω, the gravity coefficients of Titania are more relevant for inclined orbits, I = 80-90°.
For a = 875-950 km, the orbits reach lifetimes between 200 and 400 days. This lifetime is also experienced for orbits with a semi-major axis ranging from 950 km to 1050 km and I = 75-80°, shown in the smaller island. From 810 km to 1100 km, the gravity coefficients of Titania can still extend the duration of the orbits up to 100 days for inclinations in the range of 80-90°. Before the inclusion of ω = 55 • and Ω = 155 • , the orbits that were most affected by the gravity coefficients of Titania had a maximum lifetime up to 60 days. Considering these specific values of those angles, this time has been extended to 400 days, a very large difference.
In the case of Uranus, the differences point out that, when ω and Ω angles are not zero, there is a deviation in the location of the orbits where their contribution was most significant.
When ω = Ω = 0 • , the most inclined orbits were the most affected by the third body. After that there are orbits with inclinations between 75 and 76°with lifetimes of 100-300 days.
For the case e = 10 −3 , the best initial conditions presented in Figure 2 point to a semi-major axis in the range of 834-1265 km and an inclination of 70-80°. In Figure 5 we present a diagram ω × Ω for some specific cases. In Figure 5a, a = 1265 km and I = 70 • , whereas in Figure 5b, the initial conditions are a = 1265 km and I = 80 • . Figure 5a shows large islands with orbits that reach 450 days. These islands are located at ω = 100-200°, ω = 300-360°, Ω = 50 • − 200 • , and Ω = 250-360°. In the case of Figure 5b, where the inclination is 80 • , the "protection mechanism" caused by the combination of the terms J 2 and C 22 increases the orbital duration to 550 days. However, these islands are much smaller and can be seen at ω = 125-180°, ω = 300 • − 360 • and Ω = 50 • and Ω = 250 • .  Figure 5c,d have the same inclination 80 • , but different semi-major axes. In Figure 5c, with a = 834 km, the maximum lifetime is almost 3 years. However, the regions are smaller compared to those seen in Figure 5a,b. They are found with ω around 150 • and Ω close to 325 • . More distant orbits, with a semi-major axis equal to 1000 km, (Figure 5d), it is noted that the lifetime decreases, and the values of ω and Ω for orbits with the longest duration are 100 • and 275 • , respectively.
An initial condition (a point) with values of ω = 165 • and Ω = 205 • is highlighted (black circle) in Figure 5b. The values of these angles were added in the new simulations for the initial conditions of Figure 2, in order to investigate how they would affect the lifetime of the probe. The new investigation was carried out by considering only the perturbation of the third body and with the complete system, third body plus the J 2 , C 22 coefficients of Titania. Shown in Figure 6, it is observed that the use of the non-zero values for these angles in e = 10 −3 increases by 100% of the lifetime of the probe in both cases. When both effects are included, the gravitational attraction of Uranus and the gravity coefficients of Titania, the use of non-zero values for ω and Ω also changed those regions where the lifetimes are longer. Orbits with a lifetime up to 600 days are found for inclinations ranging from 84 • to 90 • and a = 810 km to a = 1000 km (third-body effects). When considering the gravity coefficients of Titania, orbits with this lifetime appear for an interval of semi-major axes equal to a = 810-1050 km and I in the range 75-90°(light green squares). The location of the orbits with a lifetime of approximately 400 days is practically the same for both cases, a = 810-1000 km and I = 75-90°(green squares). The blue dot region has lifetimes in the range 200-300 days for all inclinations and with a in the range 1000-1600 km. Figure 6c shows that Uranus affects the lifetime of the more inclined orbits. This effect appears in a small orange island that goes from I = 84-90°for a = 810-1300 km. In this region, the duration of the orbits is approximately 400 days, considering the lifetime of the probe affected by the gravity coefficients of Titania and Uranus. The effect of J 2 and C 22 of Titania and the third-body acts on orbits with inclinations of 75-90°and semi-major axes of 810-1000 km. The lifetimes of these orbits range from 100 to 600 days compared to the case where only the effect of Uranus is considered; there are no significant changes for a larger than 1400 km.
For e = 10 −2 , the longest-lived orbits are located in the range 900-1160 km and I = 70-80°. Therefore, the initial conditions used to build the ω × Ω graphs are within this range of values. In Figure 7 we present four scenarios for analysis of these angles. In Figure 7a, the semi-major axis considered is 1000 km and the inclination is 70 • , whereas in Figure 7b, a also assumed the value 1000 km and I = 70 • . There is a similarity between the regions of orbits with longer orbital duration shown in Figure 7a,b. In both cases, the islands with the longest lifetime are located for values of Ω close to 150 • , whereas for ω the best values are around 140 • and 325 • . Despite the high value of the eccentricity, with the combination of Titania's gravity terms with the third body perturbation, long-lived orbits can be found. In Figure 7a the maximum lifetime is 550 days, whereas in Figure 7b, where the inclination is greater, the duration of the orbits is up to 650 days.
The best orbits for I = 80 • (Figure 7c,d), are for Ω = 150 • and for ω close to 140 • and 325 • . Due to the proximity to the surface of Titania, the lifetime for the cases analyzed in Figure 7d is longer, 700 days. New numerical simulations for e = 10 −2 were made according to the initial condition highlighted in Figure 7c (black circle), where the orbital elements are equal to a = 1162 km, I = 80 • , ω = 145 • , and Ω = 95 • . These values for the argument of periapsis and the longitude of the ascending node were used in the new simulations with the same values of a, e, and I used in Figure 3. The results of the inclusion of these new values of ω and Ω are shown in Figure 8.
When considering ω and Ω non-zero, the increase of the lifetime of the probe is 2.5 times greater in comparison with ω = Ω = 0 • . The lifetimes increased from 180 to 450 days.
In Figure 8a, the orbits closer to the surface of Titania have lifetimes of 250-300 days with a semi-major axis in the range 810-1000 km. It occurs for all inclinations when we considered the perturbation of Uranus. When combining the effects of Uranus and the gravity coefficients of Titania, Figure 8b, orbits with this lifetime occupy a slightly larger region, a = 810-1350 km, also for all analyzed values of I. The islands with lifetimes up to 450 days are seen in both cases: when only the perturbation of Uranus is considered, and when the gravity coefficients of Titania are also taken into account. However, these orbits change locations according to each perturbation investigated. In Figure 8a, these orbits can be seen from I = 75-90°and a = 1300-1500 km, whereas in Figure 8b, this lifetime is observed for orbits in the range a = 1100-1400 km and I =75-76°(orange and red dots).
The orbits with a = 900-1450 km and I =75-82°have longer lifetimes when subjected to the effect of the gravity coefficients of Titania. The difference in the lifetimes is up to 50 days, in comparison with the situation where only the influence of the third body is considered. These orbits can be seen in Figure 8c (green dots).
Numerical integration for e = 10 −1 showed that the best values for a and I would be 900 km and I = 80 • , respectively. We use these values and plot the lifetimes as a function of ω × Ω varying these angles over the range 0-360°. The results of this analysis are not shown here in this work, as they are similar to the case e = 10 −2 . However, it is important to note that the adoption of specific values for ω and Ω also increased the orbit lifetime.
The results shown in Figure 6b point out that the angles ω and Ω contribute to the increase in the lifetimes for orbits with high eccentricities, as in the case e = 10 −1 . With the values found for these angles, the duration of the orbits around Titania increased by two times when considering only the third body perturbation and more than eight times when the gravity coefficients of Titania were also included. In the analysis involving the contribution from Uranus gravitational attraction, the regions of greater and smaller duration remained similar, increasing only the maximum lifetime.
However, exploring the scenario where J 2 and C 22 of Titania were added to the model, the orbits with zero and non-zero values of ω and Ω had very different maximum lifetimes as well as larger regions with longer lifetimes. Without the inclusion of specific values for these angles, the orbits with the longest lifetime were located at I = 75-90°and a = 900-1200 km. In this interval, the maximum times reached by these orbits were 60 days. With the adoption of better values for ω and Ω in the new simulations, these orbits started to have lifetimes from 200 to 500 days, for a = 900-1000 km and I = 75-90°. The analysis carried out in this section shows the importance of a different approach, first including variations in ω and Ω, and then mapping a × I for the best values of ω and Ω. In all cases investigated, the best values for ω can be found in two regions, in the ranges 100-150°and 27-360°. In the case of the ascending node longitude, the best values are close to 50-100°and 250-360°. The orbits that had their lifetimes extended were those orbits closer to Titania's surface and with the greatest inclination. The reason for this is associated with the fact that lower and more inclined orbits are subject to the action of the zonal oblateness coefficient J 2 , which, added to the term C 22 and the perturbation caused by the third body, promotes an equilibrium capable of increasing the duration of the orbit.

Long-Duration Orbits
Previous studies of systems similar to the one presented here often use the double averaging method to find frozen orbits. According to the work of [9], in order to find a frozen orbit, the relations de dt = 0, di dt = 0, dω dt = 0 must be satisfied. For such relationships to be met, it is common in the literature to define ω = π/2 or ω = 3π/2, because, according to the Kozai-Lidov mechanism, under the presence of the third body, the argument of periapsis librate around these values as the eccentricity and inclination oscillate [17,18]. However, when these values are fixed, di dt is non-zero due to the term C 22 , which causes the inclination to increase significantly and hence the eccentricity. In other work, where the double mean model is used, it is possible to eliminate the short period terms that contain the C 22 coefficient and then obtain the necessary conditions to find the frozen orbits.
In this work, we are not interested in finding frozen orbits, whose noted derivatives are equal to zero. Despite showing similarity to a frozen orbit, orbits that last longer are our main interest, independently of the value of these derivatives. This is because our goal is to analyze more realistic models and investigate the isolated effect of the term C 22 on the desired orbits.
In this section, we present a study on the isolated effect of each perturbation on the evolution of some orbital elements of a probe for some specific cases. The perturbations considered here are: the gravitational effect of Uranus, the gravity coefficient J 2 from Titania, the ellipticity of Titania (C 22 ), Uranus + J 2 from Titania, Uranus + C 22 of Titania, and Uranus + J 2 + C 22 of Titania. This analysis is presented in Figures 9 and 10, where the temporal evolution of the orbital elements e, I, and ω are shown. In all cases analyzed, the numerical integration stops when the probe collides with the surface of Titania. Figure 9a shows the temporal evolution of the eccentricity under the effect of these perturbations. It can be observed that the effect of the J 2 and C 22 gravity coefficients of Titania are small compared to the effect caused by the third body. This is seen even when the other perturbations are added. When J 2 and C 22 are included in the system, the lifetime increases, as already noted; however, the eccentricity continues to increase due to the presence of the third body.
In the circular restricted three body problem, the Kozai-Lidov effect is common in some cases. This phenomenon is directly related to the presence of the third body, which, as shown in Figure 9, causes the eccentricity of the probe to reach high values in a short period of time. The Kozai-Lidov mechanism is still capable of making the argument of periapsis librate around a constant value. These values are usually 90 • or 270 • . However, in cases where the eccentricity grows rapidly, as is the case of the present work, a region of libration is not found for the diagram e × ω (e.g., see [9]). When J 2 and C 22 of Titania are included, there is a balance between the effects caused by the gravity coefficients of Titania and those caused by the Kozai-Lidov effect. This balance causes a slow increase in the eccentricity and, consequently, the lifetime of the orbit increases.
Regarding the evolution of the inclination, shown in Figure 9b, the third body and the J 2 coefficient do not cause large variations in the probe inclination. The variation caused by J 2 is small and does not exceed 1 • , although the lifetime is 1000 days. When Uranus is considered, the variation is small and the lifetime is just 300 days. When only Uranus and J 2 are taken into account, it is noted that the presence of the third body reduces lifetime compared to the case when only the J 2 coefficient was considered.
The effect responsible for causing large variations in the inclination is the tesseral coefficient C 22 . When we analyze the dynamics only assuming C 22 , the inclination increases from 80 • to 100 • degrees and the lifetime is 1000 days. In all cases in which C 22 is considered, the variation of I is large. This can be seen in the cases Uranus +C 22 and Uranus +J 2 + C 22 , in which the variation of I can reach 30 • . This effect caused by C 22 is expected, as shown in previous work [9,15,16,21].
In the analysis of the argument of periapsis, presented in Figure 10, it is possible to note that, as well as in the evolution of the eccentricity shown in Figure 9a, the third body is responsible for causing a greater perturbation in this element. The behaviour of ω is similar in all cases where Uranus is added. The argument undergoes a great variation in the first days, with oscillations ranging from 0 • to 360 • . After a certain time, the argument of periapsis begins to librate around ∼50°for a short period of time. An expected behaviour for ω would be its libration around π/2 or 3π/2, as described by 18]. However, as the short period terms that contain the C 22 coefficient were not eliminated, the inclination grows very quickly and, as the eccentricity depends on the slope, it also reaches high values in a few days.
When only the effects of J 2 of Titania are assumed, the argument of periapsis reaches 360 • and then decays to approximately 240 • until a collision occurs in 1000 days. When only the perturbation due to C 22 is considered, the lifetime of the argument of periapsis is also 1000 days; however, the argument of periapsis increases to approximately 360 • before the collision.
The relevance of the ω and Ω angles have already been evident in the lifetime maps shown in Figures 1-8. In Figure 11, an analysis of the influence of these angles on the orbital elements of two orbits is presented. The orbital elements analyzed in Figure 11 have different values of ω. The first adopted value of ω and the other values of the orbital elements were taken from the point highlighted in Figure 7 (black circle).
The second orbit has the same initial conditions as the highlighted initial condition, except ω. The value of ω for the highlighted point is 145 • . The use of this value is able to significantly increase the lifetimes of the orbits. For the sake of comparison, a value of ω just below the circled point was chosen, ω = 100 • .
In all cases shown in Figure 11, the value of ω = 145 • causes an increase in the lifetimes of the probe. The variation of the semi-major axis, shown in Figure 11a, is very small in both cases. For ω = 100 • and ω = 145 • , the variations are no more than 1 km. For the eccentricity, shown in Figure 11b, the value of ω = 145 • attenuates its rapid growth. It remains almost constant for 300 days and only then begins to grow, reaching a maximum of 0.3 in 430 days.
Using ω = 100 • , the eccentricity increases much faster, and in just 135 days, the collision with the surface of Titania occurs. The evolution of the inclination, shown in Figure 11c, presents a variation of 30 • when ω = 145 • , but has a long lifetime, 430 days. Using ω = 100 • , the variation is only 5 • , but the collision happens in a few days.
In the case of the evolution of the argument of periapsis, Figure 11d, the behaviour is similar to the one shown in Figure 11 in the cases where the third body was considered. For both ω = 100 • and ω = 145 • , the argument of periapsis oscillates around a constant value for some time. When we consider ω = 100 • , the argument decays and librates around 40 • for approximately 80 days. For ω = 145 • , the argument of periapsis circulates for about 200 days and then proceeds to librate around 238 • for almost 230 days. This behaviour is expected due to the aforementioned Kozai-Lidov mechanism.
In Figure 11e, the evolution of the ascending node is presented for ω = 100 • and ω = 145 • . For ω = 100 • , the ascending node decays almost 25 • in 140 days, whereas for ω = 145 • this decay is almost 45 • during the first 250 days of integration. After that, the node starts to describe a periodic behaviour until the collision occurs.
In order to gather all the best orbits investigated, a summary of the best lifetimes for each initial condition and for each case analyzed is presented in Table 3. This table presents the best initial conditions for the five values of the eccentricity studied in this work. In addition to the results for ω = Ω = 0 • , the results considering ω and Ω different from 0 • are also presented. The table shows the best initial conditions obtained for long-duration orbits, with these conditions given by: a = 826 km, e = 10 −3 , I = 86.2 • , ω = 165 • , and Ω = 205 • considering only the third body perturbation. For the system, third-body +J 2 , C 22 of Titania, we have a = 874 km, e = 10 −3 , I = 80 • , ω = 165 • , and Ω = 205 • . In both cases, the orbits have a lifetime of 1000 days. In Figure 12 we present different projections of one of the orbits with the best initial conditions found in this work. The trajectory was numerically simulated for 991 days, considering the complete system, the gravitational effects of the third body (Uranus), and Titania's gravi-tational coefficients (J 2 and C 22 ). The initial conditions are: a = 898 km, e = 10 −3 , I = 81.4 • , ω = 140 • , and Ω = 225 • .

Analyzing the Gravitational Coefficients of Titania
In order to analyze how possible errors in the values of Titania's gravitational coefficients can affect the probe's lifetime, our next results present a detailed analysis of possible error values in both coefficients. We also present multiple regression models in order to find a relationship as a function of the coefficients J 2 and C 22 capable of predicting error in the probe's lifetime. Figure 13 shows the response surface found for the more complete regression model that matches the equation.
We investigate a range from −10% to 10% for J 2 and C 22 . To find the noted function, we consider only the case where the lifetime variation is significant compared to the cases obtained from [27]. The initial conditions for the regression were taken from the percentages analyzed for e = 10 −3 , as it shows the largest increasing in the lifetime. We present four multiple regression models for different values of the coefficient of determination R 2 . The models and their fit coefficients are presented in Equations (8)- (11).
where Y is the lifetime dependent variable, and J 2 and C 22 are independent variables. It is important to highlight that in the values of J 2 and C 22 were added percentages ranging from −10% to +10% with respect to the nominal value. Equation (8) presents a relationship for the lifetime as a function the percentages of J 2 and C 22 . For this model (Model 1), the coefficient of determination R 2 has a value of 0.84. Although the value of R 2 is not as close as to 1 and the model does not present an interaction term between J 2 and C 22 , the function provides a good approximation of the lifetime as a function of the percentage values of Titania's gravitational coefficients. For Equations (9)- (11), the interactions between the percentages of the coefficients are analyzed. With a coefficient R 2 = 0.92, Equation (11) (Model 4) presents a better precision for estimating the lifetime as a function of J 2 and C 22 . The graphical representation of this model is shown in Figure 13.

Final Comments
In this work we study the dynamics of orbits around the natural satellite Titania, the largest moon in the Uranus system. The gravitational attraction due to Uranus is considered, as well as the gravity coefficients J 2 and C 22 of Titania. Titania is at the center of the system, and Uranus is in an elliptical orbit and acts as the third body to perturb the motion of the space probe around Titania.
Through a set of numerical simulatons of the equations of motion, we build lifetime maps for different values of eccentricities of the orbit of the space probe. The construction of these maps allowed us to find values of a, e, I, ω, and Ω capable of increasing the lifetime of the probe. We also constructed lifetime "difference maps" in order to analyze the regions where the effect of the third body and the gravity coefficients of Titania were in equilibrium, allowing longer life for the space probe.
We also present an analysis of the importance of the angles ω and Ω in the lifetime of the probe. In addition, we studied the response of each effect on the orbital element.
Our results showed that, for the altitudes adopted in this work, the oblateness coefficients J 2 and J 4 of Uranus do not affect the orbits of the space probe. We also showed that the orbits with a longer lifetime had a semi-major axis close to 900 km, whereas the best eccentricity is between 10 −3 and 10 −4 when considering ω = Ω = 0 • . The results also point out that the inclusion of specific values of ω and Ω are extremely important in increasing the lifetime of the probe. In all cases, the increase in the lifetime was more than 50% after adopting specific values for these two angles. The best values for ω are around 145-165°and for Ω around 55-205°. In some cases, including these values resulted in the probe lifetime increasing up to eight times compared to cases with ω = Ω = 0 • .
Our main objective is to find long-duration orbits around Titania. After analyzing the results including the values of ω and Ω, we found ideal initial conditions for an orbit lasting up to 1000 days. This orbit has initial orbital elements given by: a = 874 km, e = 10 −3 , I = 80 • , ω = 165 • , and Ω = 205 • .
Analyzing the effect of each perturbation, we show that the third body is responsible for causing a greater perturbation in the eccentricity and argument of periapsis, whereas the ellipticity coefficient C 22 of Titania causes great oscillations in the inclination. The zonal coefficient J 2 , when acting alone, has little influence on the orbital elements; however, when added to the other perturbations, it helps to promote a balance between third body and gravity coefficients, and then we have the so-called "protection mechanism" prolonging the lifetimes of the orbits.
The gravitational coefficients of Titania have not yet been used to investigate orbits close to its surface. Thus, we present an original study on the sensitivity of the lifetime as a function of some possible errors attributed to the values of J 2 and C 22 . Our results show that the hypothetical errors in the values of Titania's gravitational coefficients can greatly increase or decrease the lifetime of the orbit, especially when the eccentricity is equal to 10 −3 . To try to predict the lifetime as a function of these errors, we present multiple regression models for the lifetime as a function of the coefficients J 2 and C 22 . We found a function with a value of R 2 equals to 0.92, capable of providing a good approximation of the lifetime for values of J 2 and C 22 between −10% and +10% of the values given in [27].
The results presented in this paper can help plan missions to Titania, which can bring important scientific data, in particular considering that the Uranus system has been studied mostly from Earth.

Conflicts of Interest:
The authors declare that they have no conflicts of interest with research institutions, professionals, researchers, and/or financial supports.

Ethical Approval:
The submitted work is original and has not been published elsewhere in any form or language. The work presents the results of a single study. The results are presented clearly, honestly, and without fabrication, falsification, or inappropriate data manipulation. No data, text, or theories by others are presented as if they were the authors' own.