# Mapping Orbits regarding Perturbations due to the Gravitational Field of a Cube

Authors: Flaviane C. F. Venditti; Antonio F. B. A. Prado
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

The orbital dynamics around irregular shaped bodies is an actual topic in astrodynamics, because celestial bodies are not perfect spheres.  When it comes to small celestial bodies, like asteroids and comets, it is even more import to consider the nonspherical shape. The gravitational field around them may generate trajectories that are different from Keplerian orbits. Modeling an irregular body can be a hard task, especially because it is difficult to know the exact shape when observing it from the Earth, due to their small sizes and long distances. Some asteroids have been observed, but it is still a small amount compared to all existing asteroids in the Solar System. An approximation of their shape can be made as a sum of several known geometric shapes. Some three-dimensional figures have closed equations for the potential and, in this work, the formulation of a cube is considered. The results give the mappings showing the orbits that are less perturbed and then have a good potential to be used by spacecrafts that need to minimize station-keeping maneuvers. Points in the orbit that minimizes the perturbations are found and they can be used for constellations of nanosatellites.


---

## Body

## 1. Introduction

The main reason to study the properties of irregular bodies in space science is because there are many irregular small bodies in the Solar System. Most of them are located between the orbits of Mars and Jupiter, in the main asteroid belt. According to their orbits, comets and asteroids are also called NEOs, abbreviation for near Earth objects. For asteroids, they can also be called NEAs, meaning near Earth asteroids. These are objects with orbits passing near the orbit of the Earth, as the name already says. The NEAs are divided in four groups, depending on the perihelion distance, aphelion distance and the semimajor axis of the orbit: Atiras, Athens, Apollos, and Amors [1].Geometric figures with nonspherical shapes can be used as a model for the irregular bodies, and this technique can be very useful in the preliminary study of trajectories around asteroids. Initial works about the formulation for the potential of rectangular objects were made by MacMillan [2], and this type of research continued being developed on the next decades [3–7]. Studies of the property of orbits around elongated rectangular objects were also developed by Broucke and Prado [8].There are some works related to the potential of a cube. Michalodimitrakis and Bozis [9] studied the two-body problem with a massive point particle and a homogeneous cube and proved the existence of a ring-type bounded motion. The gravitational field around a cube was also analyzed by Chappell et al. [10]. The search for periodic orbits around a fixed homogeneous cube and the potential based on the polyhedron method for a cube given by a volume integral were studied by Liu et al. [11, 12].In this work, the gravitational perturbation due to the shape of a cube is analyzed. The formulation of the potential for the cube is given by closed equations [2, 13]. The perturbation is studied with a method that measures the amount of the change made in the velocity over the time, which can be seen in more details in previous works by Prado [14, 15], Sanchez et al. [16], and Oliveira and Prado [17]. This methodology consists in calculating the integral of the acceleration over the time, which will be referred here as the perturbation integral and will be explained in the next sessions.


## 2. Problem Formulation

The potential of a cube can be found in the literature in terms of spherical harmonics, and the expansion limited to order four is shown here in (1). An expansion with higher degrees can be found in Broucke and Prado [8]
(1)V′=Gmr-edr4C40P4sin⁡ϕ+C44P44cos⁡4λ,
where ϕ is the latitude; λ is the longitude; P4sin⁡ϕ represents a Legendre polynomial in sin(ϕ)=z/r; P44 is the associate Legendre polynomial, where P44=105cos⁡(ϕ); and the two stokes coefficients are C40=-7/30 and C44=-1/720.Approximations of the potential for a cube can be found in the literature as a closed equation [2] and is given by
(2)V′=Gmr-7ed4Gmx4+y4+z4-3x2y2+x2z2+y2z230r9,
where G is the gravitational constant 6.67384 × 10−11 (m3/kg·s2); m is the mass of the cube; ed is the length of the edges of the cube; and r is the distance between the spacecraft and the center of the body.Another closed expression for the potential of a cube can be seen in Kellog [13], in the form of the equation shown below:
(3)V′=Gmr+7ed4Gm3r4-5x4+y4+z460r9.The termr is replaced by
(4)r=x2+y2+z2.Having the equation for the potential, the first step is to obtain the acceleration by calculating the partial derivatives of the potential. Neglecting the Keplerian term, to consider only the perturbative potential, the derivatives with respect tox, y, and z of the disturbing potential are shown in (5). They are obtained using the potential shown in (2). The potential without the Keplerian term (perturbing term) is indicated by V instead of V′(5)∂V∂x=7ed4Gmxx4-5x2y2+z2+3y4-y2z2+z46x2+y2+z211/2,∂V∂y=7ed4Gmy3x4+y4-5y2z2+3z4-x25y2+3z26x2+y2+z211/2,∂V∂z=7ed4Gmz3x4+3y4-5y2z2+z4-x23y2+5z26x2+y2+z211/2.The acceleration for the cube is obtained in order to study the effects of the shape of the body. By calculating the integral of the acceleration for one orbital period, the perturbation can be analyzed in terms of the variation of the velocity caused in the spacecraft. This integral, that is given by (6), will be called the perturbation integral (PI), which can be seen in previous works by Prado [14, 15], and Oliveira and Prado [17]
(6)PI=∫0TGradVdt,
where GradV is the gradient of the potential.The expression for the mean anomaly is shown by (7). Then, after some arrangements, (8) is replaced in (6), so that the integral PI can be evaluated in terms of the eccentric anomaly [14, 15, 18]
(7)M=M0+nt-t0,M=E-esin⁡E,(8)dM=ndt=1-ecos⁡EdE,(9)PI=1n∫02πGradV1-ecos⁡EdE,
where M is the mean anomaly, n is the mean motion, e is the eccentricity, and E is the eccentric anomaly of the spacecraft orbiting the cube.After the substitutions are made, the expression for the perturbation integral of a cube may now be obtained by calculating the integral given by(10)PI=1n·∫02π∂V∂x2+∂V∂y2+∂V∂z21-ecos⁡EdE.The perturbation integral depends on the initial values of the orbital elements. The results will be shown for several different conditions, and the search for the least perturbed orbits can be obtained from those mappings.The next step is to consider that small bodies, in the same way as the planets, are spinning [19]. Most of them are rotating around one principal axis, but some asteroids were found to be spinning in more than one axis like, for example, the asteroid Toutatis, said to be in tumbling rotation [20–22]. To take into account that rotation is an important factor, even more in the case of gravitational perturbation of non-spherical bodies. Therefore, the results will also consider the cases of a rotating cube.


## 3. Results

The next sessions will show some results for the PI under different conditions, in order to map the orbits to find the ones that are less perturbed.


### 3.1. Varying Initial Configurations

The results presented next intend to show the differences of having a body that cannot be treated as a point mass. An orbit around an irregular body is different from a Keplerian orbit, and the perturbation generated due to the asymmetric shape of the orbited body must be considered.Figure1 shows the acceleration using the potential for one revolution around the cube, considering inclinations of 0°, 45°, and 90° for the orbit of the spacecraft. The green line shows a Keplerian orbit around a sphere that has the same mass of the cube, as a comparison. The following results were obtained for a circular orbit using the value of 4 × 1014 kg for the mass of the cube, length of the edge of the cube of 5 km, and orbital radius of 10 km.Figure 1
Acceleration for one revolution around a cube, for equatorial, polar orbit,i=45°, and a Keplerian orbit.It can be noticed from Figures1 and 2 that the acceleration for equatorial and polar orbits is the same, as expected, due to the symmetry of the cube. Equatorial and polar orbits make the spacecraft pass exactly by the same mass distribution. The spacecraft starts its motion in front of the middle point of the face, which is the position with the smallest potential, due to the largest distance between the subspacecraft point and the spacecraft. The potential increase with the spacecraft going to the position near the edge of the cube reaches its maximum at this point and then decreases again until the middle of the next face is reached. This geometry is repeated until a complete revolution around the cube is completed. For an orbit with inclination of 45°, the spacecraft also starts its motion near the face, in the position of the minimum potential. Then it moves in the direction of the vertex, where the potential reaches its maximum value, which is higher than its equivalent for equatorial and polar orbits. After that, the spacecraft passes near the middle point of the edge, where it reaches an intermediate minimum (between both maximum). The acceleration at this point is smaller than the one achieved by the spacecraft at the vertex but higher than the values when it is crossing the equatorial plane. The motion is repeated again for a full revolution. It is also visible from Figure 1(a) that the potential affecting the spacecraft is smaller than the Keplerian one when the spacecraft passes near the face of the cube, which is explained by the largest distance of the face to the spacecraft when compared to the distance from a sphere of the same mass. The opposite occurs when the spacecraft passes near the edges, and near this position the potential is larger than the equivalent one from the sphere. It is important to note that the difference in the potential of the cube when compared to the sphere is larger when passing near the face, which means that the potential is smaller near the faces but the perturbation is larger. The points where the green line (Keplerian orbit) crosses with the potential of the cube represent positions where the magnitude of the potential of the sphere is the same as the potential from the cube, but they are not equilibrium points, because the three components of the acceleration are not the same. But they are points where the difference of the potentials is interesting and deserves to be investigated. Figure 2 shows the acceleration due to the perturbation term only, after subtracting the Keplerian term, which also represents the difference between the potential of a cube and a sphere. There are the points where this difference is at minimum. Those points are the best ones to locate a satellite that is moving in a synchronous orbit with the cube (same period of rotation), because the efforts to keep this geometry are at minimum. Those points are interesting to place spacecrafts that have the goal of staying in near stable orbits. Since most of the efforts in interplanetary research will probably be done by constellations of nanosatellites instead of single larger spacecrafts, one or more of the satellites could stay in one of those points to be used as a communication satellite to transfer information from the others, in different orbits to explore more the body, to the Earth. They would be similar to geostationary satellites in orbit around the Earth. For equatorial and polar orbits there are periodic oscillations, for the physical reasons already explained, based in the regions of the cube (face or edge) that the spacecraft is passing by. The inclined orbit shows more oscillations, with six maxima and minima for each period, due to the passages by the faces, edges, and vertexes.Figure 2
Acceleration considering the disturbing term for one revolution around a cube, for equatorial and polar orbit,i=45°.The perturbation integral, or PI, gives the variation of velocity over the time after one revolution of the spacecraft. Therefore, by analyzing the variation of the values of PI, it is possible to obtain the most or the least perturbed orbits, depending on the initial orbital elements chosen. This integral is obtained by the integration of the equation for the acceleration of the cube (showed by (5)) for one orbital period of the spacecraft. For a circular orbit with argument of periapsis and right ascension of ascending node equal zero, the PI is given by
(11)PI=49ed8G2M2π26865+372cos⁡4i+315cos⁡8i147456na12.This analytical equation shows the dependence of the PI with the semimajor axis and inclination of the orbit, and it allows us to find the orbits that are less perturbed. The dependence with the semimajor axis quantifies the expected fast decrease with increase of this variable. A simple analysis shows that when the spacecraft is far from the cube, it feels that the gravitational potential is like the one generated by a sphere, so the perturbations go to zero. Although expected, this integral allows a better comparison of values by presenting an objective form to quantify the different effects of the geometric shape. A more interesting result is the study with respect to the inclination. The values for the perturbation integral that are obtained by (11) shown in a graph as a function of the inclination. Figure 3 shows the results. The existence of four points of global minimum, two global maxima, and two more local maxima is noticed. Those extreme points represent orbits that would require maximum and minimum fuel consumptions for a spacecraft that needs to make station keeping maneuvers to stay all the time in a Keplerian orbit. This hypothesis is usually too strong and not common in real missions, but it is used here to map trajectories, expecting that the more perturbed orbits present faster and larger variations of the Keplerian elements, with a good potential to require more fuel consumption for the station keeping maneuvers. This is the reason for justifying this type of mappings for the orbits. Figure 3 shows that the minima are located at the inclinations of 0.45, 1.10, 2.05, and 2.65 radians and the maxima are located at the inclinations of zero and π radians. The difference in the values of the PI is near 0.000006 m/s, which represents about 14%, not a negligible amount. For orbits with smaller semimajor axis, this difference increases in absolute values. It means that a graph similar to Figure 3 is important for mission designers planning orbits around a body that has a geometrical form similar to a cube. A deep physical analysis of the perturbations can explain that the maximum PI for the polar and equatorial orbits is due to the larger distances from the spacecraft to the cube, when compared to a sphere with the same mass, because it passes near the faces of the cube. So the perturbation results from a gravity field that is weaker than the one generated by a sphere or point of mass. The spacecraft tends to move away from the body. In the opposite side, the 45-degree inclined orbit has its maximum value for the PI due to the smaller distances from the spacecraft to the cube, when compared to a sphere with the same mass, because it is passing near the edges and the vertex of the cube. So the perturbation results from a gravity field that is stronger than the one generated by a sphere or point of mass and the spacecraft is more attracted by the body, tending to reduce the distance spacecraft cube. Figure 3 shows that the value of the PI is larger for the polar and equatorial orbits, when compared to the 45-degree inclined orbit, which means that the excess of distances experienced by the polar and equatorial orbits is larger, in magnitude, than the smaller distances experienced by the 45-degree inclined orbit.Figure 3
Perturbation integral versus inclination for circular orbit around a cube.The previous results were obtained considering circular orbits. It is now interesting to study eccentric orbits. A graph of the variation of the PI with respect to the inclination is shown in Figure4, now considering different values for the eccentricity. The initial values used for the eccentricity were zero (circular orbit), 0.1, 0.2, and 0.3. The effect of the eccentricity for equatorial and polar orbits are the same, which can be seen better in Figure 4(b) that plotted the variation of PI with the eccentricity for orbits with inclination of 0°, 90°, and 45°. Figure 4(a) shows the same values, but the eccentricity goes up to 0.7, which generate results that are too similar to allow a good comparison among the orbits. It is demonstrated that, as the eccentricity increases, the higher are the values for the PI and consequently eccentric orbits are more perturbed than circular orbits. This happens slightly faster for orbits with inclination of 45°. The physical reason for this effect is that eccentric orbits cause closer approximations of the spacecraft to the cube, during the periapsis passages, so it feels stronger perturbations during these passages and these passages contribute to the accumulated effects of the perturbation during one orbital period. Figure 4(b) shows the same information shown in Figure 4(a), but reducing the interval of eccentricities to a maximum of 0.3. It is possible to see an inversion in the values of the PI at e=0.15. The polar and the equatorial orbits are less perturbed than the inclined 45° after this limit and less perturbed before that.(a) Perturbation integral versus eccentricity. Solid line is for equatorial and polar orbits, and dashed line is fori=45°, for eccentricities up to 0.7. (b) Perturbation integral versus eccentricity. Solid line is for equatorial and polar orbits, and dashed line is for i=45°, for eccentricities up to 0.3.
(a)
(b)Figure5(a) shows the value of the PI (m/s) as a function of the inclination for different values of the eccentricity. The same oscillations are noticed for each value of the eccentricity. The results also emphasize the increase of the PI with the eccentricity. It is visible that the eccentricity has stronger effects when compared to the inclination. Figures 5(b), 5(c), and 5(d) show the same information of Figure 5(a), but for each eccentricity in order to give more details about each curve. From Figure 5(a) it is visible that the amplitude of oscillation of the curves increases with the eccentricity of the orbit. The combination of Figures 5(a), 5(b), 5(c), and 5(d) show that the increase of the eccentricity makes the values of the local maximum (inclination near 45 degrees) increase until a point (near e=0.2) where they become the global maximum. The global maxima of the circular orbits (equatorial and polar orbits) decrease in magnitude. When the eccentricity reaches the value 0.3, the figure shows only two maxima. The explanation for this phenomenon comes from the analysis already made regarding the physical reasons to have maximum at the polar, equatorial, and 45-degree inclined orbits. The polar and the equatorial orbits have their maximum based on the larger distances spacecraft cube, when compared to a sphere of the same mass, so the increase of the eccentricity that causes an approximation spacecraft cube during the passages of the spacecraft by the periapsis of the orbit works against the physical reason of the perturbation. In the opposite side, the 45-degree inclined orbits have their maximum based on the smaller distances spacecraft cube, when compared to a sphere of the same mass, so the increase of the eccentricity works in favor of the physical reason of the perturbation. So the net result is an increase of the PI for all orbits, because the spacecraft gets closer to the cube and it increases the perturbations, but this increase is different for each inclination. It is stronger for the 45-degree inclined orbit, causing an inversion of the large maximum. The minima are located in about the same positions, with a small shift to the left of the scale, which means that the increase of the eccentricity reduces the value of the inclination that generates the minimum perturbation. The physical reason is the same as explained above. The increase of the eccentricity makes the spacecraft pass closer to the cube, so the minimum perturbed orbits have smaller values for the inclination, to be far away from the edges and vertexes, that have their perturbations increased by the closer distances generated by larger eccentricities.(a) Perturbation integral versus inclination. Blue,e=0; red, e=0.1; green, e=0.2; orange, e=0.3. (b) Perturbation integral versus inclination for e=0.1. (c) Perturbation integral versus inclination for e=0.2. (d) Perturbation integral versus inclination for e=0.3.
(a)
(b)
(c)
(d)Another relevant analysis to be made is the influence of the semimajor axis of the orbit. Figure6 shows, for a circular equatorial orbit, the PI (m/s) as a function of the semimajor axis from 7.5 km to 25 km, for inclinations of 0, 45°, and 90°. The dependence of the perturbation relative to the semimajor axis is very clear, and when moving away from the body the value for PI decreases fast. The perturbation is larger for smaller distances, which is an expected result, because the perturbation in the potential equation for a cube changes at a rate of 1/r9. The interesting aspect of the results shown here is that they can quantify the values involved using a new approach.Figure 6
Perturbation integral versus semimajor axis fori=0=90° (blue line), and i=45° (red line), for a circular orbit.From the results presented here, there are two main things that need to be emphasized: first, the smaller the semimajor axis is, the larger the values of the PI are; second, as shown in Figure4, eccentric orbits have a higher value for the PI. Therefore, close eccentric orbits are the situations where the effects are larger, which means that the orbits are more perturbed.Besides the semimajor axis and the eccentricity, from the results obtained studying the inclination, it was proved that orbits passing near the vertex and edges of the cube are more affected and that there are maximum and minimum values for the perturbation as a function of the inclination. Another important consideration to take into account when studying nonspherical objects is that these objects are rotating. Therefore, the results for a rotating cube will be analyzed and showed next.


### 3.2. Adding Rotation

The results presented up to now considered that the cube is standing still. The next step is to consider that the cube is rotating around one of the principal axes (x, y, and z). Since the small celestial bodies are spinning, adding this fact will enrich the analysis of the gravitational perturbation due to the shape of the object.The rotation of an object may be obtained with a rotation matrix considering the direction of the rotation, which is usually represented by three matrices with respect to the axes of symmetry. If the object is rotating only in one principal axis, one of the rotation matrices is applied; therefore, knowing the direction in which the body is spinning is important. There can also be rotation in more than one axis, as already found in asteroids, and thus the sequence of the rotations must be chosen properly. The coordinates for rotation inx, y, and z axes, obtained by multiplying each rotation matrix by the position coordinates with respect to the orbital elements, are represented by the following equations:
(12)x′y′z′x=1-e2acos⁡ΩCos⁡v+o-cos⁡isinΩsinf+o1+ecos⁡f1-e2acos⁡ωtcos⁡f+osinΩ+cos⁡icos⁡Ωsinf+o1+ecos⁡f-1-e2asinisinf+osinωt1+ecos⁡f1-e2acos⁡ωtsinisinf+o1+ecos⁡f+1-e2acos⁡f+osinΩ+cos⁡icos⁡Ωsinf+osinωt1+ecos⁡f,x′y′z′y=1-e2acos⁡ωtcos⁡Ωcos⁡f+o-cos⁡isinΩsinf+o1+ecos⁡f+1-e2asinisinf+osinωt1+ecos⁡f1-e2acos⁡f+osinΩ+cos⁡icos⁡Ωsinf+o1+ecos⁡f1-e2acos⁡ωtsinisinf+o1+ecos⁡f-1-e2acos⁡Ωcos⁡f+o-cos⁡isinΩsinf+osinωt1+ecos⁡f,x′y′z′=-1-e2asinωtcos⁡o+fsinΩ+cos⁡icos⁡Ωsino+f1+ecos⁡f+1-e2acos⁡⁡ωtcos⁡Ωcos⁡o+f-cos⁡isinΩsino+f1+ecos⁡f1-e2acos⁡ωtcos⁡o+fsinΩ+cos⁡icos⁡Ωsino+f1+ecos⁡f+1-e2asinωtcos⁡Ωcos⁡o+f-cos⁡isinΩsino+f1+ecos⁡f1-e2asinisino+f1+ecos⁡f
where a is the semimajor axis, e is the eccentricity, i is the inclination, ω is the angular velocity, f is the true anomaly, o is the argument of periapsis, and Ω is the right ascension of ascending node.Several rotation periods were tested in order to compare the differences from a nonrotating situation, and among the rotations with different velocities. Reminding that the mass used was 4 × 1014 kg and the edge length of the cube is 5 km. The cube is assumed to have homogeneous density of 3.2 g/cm3, which is near the values found for the density of many asteroids [23, 24]. For all the examples that will be presented in the case of a rotating cube, the orbit is circular. Figure 7 gives the PI (m/s) as a function of the inclination for different angular velocities about the z axis.Figure 7
Perturbation integral versus inclination for different rotation periods around thez axis for circular orbits with semimajor axis of 10 km.The letterR in the figures stands for the rotation period, where R=0 means that the cube is not spinning. The physical reason to have different values when the cube is rotating is that the spacecraft passes by different parts of the cube. If it spends more time near the faces, the total perturbation increases. If it spends more time near the equilibrium position, which is a position where the potentials of a sphere and a cube are the same, the total value of the perturbation decreases. From Figure 7, it is noticeable that the results are not symmetric around the polar orbit, since prograde and retrograde orbits have different rotations with respect to the cubic body and so the regions where the subsatellite point passes are different. The equatorial orbits are the ones where the rotation generates smaller effects. The polar orbits have stronger effects when compared to the equatorial orbits, showing a decrease in the total perturbation effects. The largest reduction of the perturbation occurs near the inclination of 1 rad. In general, the rotations decreased the value of the PI. Retrograde orbits are different from direct orbits. Since they are rotating in the opposite direction of the cube, the edges and vertices spend more times by the spacecraft and it increases the perturbation. It is important to note the large reduction of the PI, going from near 0.000042 m/s for an equatorial orbit to near 0.000026 m/s in the minimum situation (inclination near 1 rad), a reduction of about 38%. The period of the rotation of the cube also has some effects in the perturbation, with the largest one simulated (24 hours) being the one that generates the smallest perturbations. Anyway, those differences are small when compared to the reduction generated by the inclusion of the rotation.Figure8 shows the perturbation integral as a function of the inclination for a cube rotating around the y axis. Each line color represents a different angular velocity, and the blue line shows the case of a fixed cube. The same representation is used in all similar figures. The figure is now symmetric with respect to the polar orbit and the minimum is near 0.55 rad. The reduction of the PI is again very large when compared to the fixed case, reaching a value near 33% at the minimum perturbation point. For the polar orbits, the perturbation is almost the same with and without rotation of the cube. When the period of rotation is 10 hours, the rotating situation generates even more perturbation. All those effects are a result of the averaged distance between the spacecraft and the mass distribution of the cube.Figure 8
Perturbation integral versus inclination for different rotation periods iny axis for circular orbits with semimajor axis of 10 km.The last case of rotation is in thex axis, showed in Figure 9. It is noticed that, now, some periods of rotation decrease too much the PI, generating orbits with lower perturbation. The results are now more dependent of the period of the rotation. The results have symmetry around the polar orbit again. In general, polar and equatorial orbits are more perturbed in the nonrotating case. When the body is rotating, the inclination of the orbit plays an important role in the stability of the orbit. For example, there are studies that suggest that retrograde orbits are more stable against the effects of perturbation around nonuniform gravitational fields [25–27]. Also, depending on the rotation axis and the period of rotation, the orbit could be synchronous or resonant. It is clear that the minimum perturbation situation changes with the inclination of the orbit, so the decision of the best orbit to place a satellite needs to take into account the period of rotation of the body.Figure 9
Perturbation integral versus inclination for different rotation periods around thex axis.The acceleration that the spacecraft feels due to a rotating cube is very dependent on the inclination of the orbit. Figure10 shows the acceleration for an orbit with different values for the inclination, with a rotation period of 5 hours around the z axis.Figure 10
Acceleration for one revolution, for different inclinations, with rotation around thez axis and period of 5 hours.Still considering a rotation period of 5 hours, Figure11 shows the acceleration for different inclinations for a cube with rotation about the y axis. In this case, the results for retrograde orbits are the same as for the prograde orbits, because the rotation axis in this case influences the orbits on the equatorial plane in the same way.Figure 11
Acceleration for one revolution, for different rotation inclinations with rotation around they axis, with rotation period of 5 hours.The last graph showing the acceleration for a rotating cube is around thex axis and can be seen in Figure 12. For this direction of rotation, the inclination of 45° is affected differently compared to the rest of the initial inclinations used.Figure 12
Acceleration for one revolution, for different rotation inclinations with rotation around thex axis, with rotation period of 5 hours.


### 3.3. Lagrange’s Planetary Equations

Considering that there is a gravitational perturbation, the variation of each orbital element can be obtained using Lagrange’s planetary equations. The disturbing potential of the cube with respect to a sphere is used, and six differential equations are obtained to express the evolution of the Keplerian elements of the orbit of the spacecraft. These equations are given by(13)dadt=2na∂V∂M,dedt=1-e2na2e∂V∂M-1-e2na2e∂V∂o,didt=-1na21-e2sin⁡i∂V∂Ω+cot⁡ina21-e2sin⁡i∂V∂o,dωdt=1-e2na2e∂V∂e-cot⁡ina21-e2sin⁡i∂V∂i,dΩdt=1na21-e2sin⁡i∂V∂i,dMdt=n-2na∂V∂a-1-e2na2e∂V∂e.The variation rate of the semimajor axis, eccentricity, and inclination, as a function of the semimajor axis, is shown in Figures13, 14, and 15, respectively. The results are for one revolution considering an equatorial orbit, a polar orbit, and an orbit with inclination of 45°. For the integration of the Lagrange equations, the magnitude was used to avoid compensations that leads to near zero integrals, but that still disturbs the orbits. This view of the mapping of the orbits is different from the view that comes from the integral approach, and it complements the study. The integral approach shows the magnitude of the total perturbation in the orbit of the spacecraft, independent of which orbital element is affected. The approach using Lagrange’s planetary equations studies the perturbation of each orbital element individually, so the relation between the results has to be seen under this perspective. In general, the graphs showed that the evolution of the orbital elements as functions of the semimajor axis, using Lagrange’s planetary equations, have a good correlation with the results presented for the perturbation integral method. It shows that the further from the body the spacecraft is, the smaller the change in the orbital elements is. It is noticed that there are some cases where the rate of variation is zero, like the polar orbit for the variation of the eccentricity and the equatorial orbit for the variation of the inclination. This point is not against the integral approach, which shows that the orbit is perturbed, because the rates of variation of the other elements are not zero. This situation just shows that the existing perturbation is disturbing another orbital element.Figure 13
d
a
/
d
t as a function of the semimajor axis around a cube.Figure 14
d
e
/
d
t as a function of the semimajor axis around a cube.Figure 15
d
i
/
d
t as a function of the semimajor axis around a cube.A complete comparison between both approaches is not possible due to the singularities of Lagrange’s planetary equations in the eccentricity and inclination. A future work can extend this comparison in further details by using nonsingular variables. One more example is possible and is shown in Figure16 that illustrates another example of the connection of the PI with the Lagrange equations, showing the variation rate of the semimajor axis as a function of the eccentricity. It shows a near zero value at the beginning, with only a slight motion for all the curves, increasing drastically after eccentricity 0.5. Similarly to what happens for the PI, the fast increase of the time derivatives of the Keplerian elements starts near e=0.5.Figure 16
d
a
/
d
t as a function of the eccentricity around a cube.


## 3.1. Varying Initial Configurations

The results presented next intend to show the differences of having a body that cannot be treated as a point mass. An orbit around an irregular body is different from a Keplerian orbit, and the perturbation generated due to the asymmetric shape of the orbited body must be considered.Figure1 shows the acceleration using the potential for one revolution around the cube, considering inclinations of 0°, 45°, and 90° for the orbit of the spacecraft. The green line shows a Keplerian orbit around a sphere that has the same mass of the cube, as a comparison. The following results were obtained for a circular orbit using the value of 4 × 1014 kg for the mass of the cube, length of the edge of the cube of 5 km, and orbital radius of 10 km.Figure 1
Acceleration for one revolution around a cube, for equatorial, polar orbit,i=45°, and a Keplerian orbit.It can be noticed from Figures1 and 2 that the acceleration for equatorial and polar orbits is the same, as expected, due to the symmetry of the cube. Equatorial and polar orbits make the spacecraft pass exactly by the same mass distribution. The spacecraft starts its motion in front of the middle point of the face, which is the position with the smallest potential, due to the largest distance between the subspacecraft point and the spacecraft. The potential increase with the spacecraft going to the position near the edge of the cube reaches its maximum at this point and then decreases again until the middle of the next face is reached. This geometry is repeated until a complete revolution around the cube is completed. For an orbit with inclination of 45°, the spacecraft also starts its motion near the face, in the position of the minimum potential. Then it moves in the direction of the vertex, where the potential reaches its maximum value, which is higher than its equivalent for equatorial and polar orbits. After that, the spacecraft passes near the middle point of the edge, where it reaches an intermediate minimum (between both maximum). The acceleration at this point is smaller than the one achieved by the spacecraft at the vertex but higher than the values when it is crossing the equatorial plane. The motion is repeated again for a full revolution. It is also visible from Figure 1(a) that the potential affecting the spacecraft is smaller than the Keplerian one when the spacecraft passes near the face of the cube, which is explained by the largest distance of the face to the spacecraft when compared to the distance from a sphere of the same mass. The opposite occurs when the spacecraft passes near the edges, and near this position the potential is larger than the equivalent one from the sphere. It is important to note that the difference in the potential of the cube when compared to the sphere is larger when passing near the face, which means that the potential is smaller near the faces but the perturbation is larger. The points where the green line (Keplerian orbit) crosses with the potential of the cube represent positions where the magnitude of the potential of the sphere is the same as the potential from the cube, but they are not equilibrium points, because the three components of the acceleration are not the same. But they are points where the difference of the potentials is interesting and deserves to be investigated. Figure 2 shows the acceleration due to the perturbation term only, after subtracting the Keplerian term, which also represents the difference between the potential of a cube and a sphere. There are the points where this difference is at minimum. Those points are the best ones to locate a satellite that is moving in a synchronous orbit with the cube (same period of rotation), because the efforts to keep this geometry are at minimum. Those points are interesting to place spacecrafts that have the goal of staying in near stable orbits. Since most of the efforts in interplanetary research will probably be done by constellations of nanosatellites instead of single larger spacecrafts, one or more of the satellites could stay in one of those points to be used as a communication satellite to transfer information from the others, in different orbits to explore more the body, to the Earth. They would be similar to geostationary satellites in orbit around the Earth. For equatorial and polar orbits there are periodic oscillations, for the physical reasons already explained, based in the regions of the cube (face or edge) that the spacecraft is passing by. The inclined orbit shows more oscillations, with six maxima and minima for each period, due to the passages by the faces, edges, and vertexes.Figure 2
Acceleration considering the disturbing term for one revolution around a cube, for equatorial and polar orbit,i=45°.The perturbation integral, or PI, gives the variation of velocity over the time after one revolution of the spacecraft. Therefore, by analyzing the variation of the values of PI, it is possible to obtain the most or the least perturbed orbits, depending on the initial orbital elements chosen. This integral is obtained by the integration of the equation for the acceleration of the cube (showed by (5)) for one orbital period of the spacecraft. For a circular orbit with argument of periapsis and right ascension of ascending node equal zero, the PI is given by
(11)PI=49ed8G2M2π26865+372cos⁡4i+315cos⁡8i147456na12.This analytical equation shows the dependence of the PI with the semimajor axis and inclination of the orbit, and it allows us to find the orbits that are less perturbed. The dependence with the semimajor axis quantifies the expected fast decrease with increase of this variable. A simple analysis shows that when the spacecraft is far from the cube, it feels that the gravitational potential is like the one generated by a sphere, so the perturbations go to zero. Although expected, this integral allows a better comparison of values by presenting an objective form to quantify the different effects of the geometric shape. A more interesting result is the study with respect to the inclination. The values for the perturbation integral that are obtained by (11) shown in a graph as a function of the inclination. Figure 3 shows the results. The existence of four points of global minimum, two global maxima, and two more local maxima is noticed. Those extreme points represent orbits that would require maximum and minimum fuel consumptions for a spacecraft that needs to make station keeping maneuvers to stay all the time in a Keplerian orbit. This hypothesis is usually too strong and not common in real missions, but it is used here to map trajectories, expecting that the more perturbed orbits present faster and larger variations of the Keplerian elements, with a good potential to require more fuel consumption for the station keeping maneuvers. This is the reason for justifying this type of mappings for the orbits. Figure 3 shows that the minima are located at the inclinations of 0.45, 1.10, 2.05, and 2.65 radians and the maxima are located at the inclinations of zero and π radians. The difference in the values of the PI is near 0.000006 m/s, which represents about 14%, not a negligible amount. For orbits with smaller semimajor axis, this difference increases in absolute values. It means that a graph similar to Figure 3 is important for mission designers planning orbits around a body that has a geometrical form similar to a cube. A deep physical analysis of the perturbations can explain that the maximum PI for the polar and equatorial orbits is due to the larger distances from the spacecraft to the cube, when compared to a sphere with the same mass, because it passes near the faces of the cube. So the perturbation results from a gravity field that is weaker than the one generated by a sphere or point of mass. The spacecraft tends to move away from the body. In the opposite side, the 45-degree inclined orbit has its maximum value for the PI due to the smaller distances from the spacecraft to the cube, when compared to a sphere with the same mass, because it is passing near the edges and the vertex of the cube. So the perturbation results from a gravity field that is stronger than the one generated by a sphere or point of mass and the spacecraft is more attracted by the body, tending to reduce the distance spacecraft cube. Figure 3 shows that the value of the PI is larger for the polar and equatorial orbits, when compared to the 45-degree inclined orbit, which means that the excess of distances experienced by the polar and equatorial orbits is larger, in magnitude, than the smaller distances experienced by the 45-degree inclined orbit.Figure 3
Perturbation integral versus inclination for circular orbit around a cube.The previous results were obtained considering circular orbits. It is now interesting to study eccentric orbits. A graph of the variation of the PI with respect to the inclination is shown in Figure4, now considering different values for the eccentricity. The initial values used for the eccentricity were zero (circular orbit), 0.1, 0.2, and 0.3. The effect of the eccentricity for equatorial and polar orbits are the same, which can be seen better in Figure 4(b) that plotted the variation of PI with the eccentricity for orbits with inclination of 0°, 90°, and 45°. Figure 4(a) shows the same values, but the eccentricity goes up to 0.7, which generate results that are too similar to allow a good comparison among the orbits. It is demonstrated that, as the eccentricity increases, the higher are the values for the PI and consequently eccentric orbits are more perturbed than circular orbits. This happens slightly faster for orbits with inclination of 45°. The physical reason for this effect is that eccentric orbits cause closer approximations of the spacecraft to the cube, during the periapsis passages, so it feels stronger perturbations during these passages and these passages contribute to the accumulated effects of the perturbation during one orbital period. Figure 4(b) shows the same information shown in Figure 4(a), but reducing the interval of eccentricities to a maximum of 0.3. It is possible to see an inversion in the values of the PI at e=0.15. The polar and the equatorial orbits are less perturbed than the inclined 45° after this limit and less perturbed before that.(a) Perturbation integral versus eccentricity. Solid line is for equatorial and polar orbits, and dashed line is fori=45°, for eccentricities up to 0.7. (b) Perturbation integral versus eccentricity. Solid line is for equatorial and polar orbits, and dashed line is for i=45°, for eccentricities up to 0.3.
(a)
(b)Figure5(a) shows the value of the PI (m/s) as a function of the inclination for different values of the eccentricity. The same oscillations are noticed for each value of the eccentricity. The results also emphasize the increase of the PI with the eccentricity. It is visible that the eccentricity has stronger effects when compared to the inclination. Figures 5(b), 5(c), and 5(d) show the same information of Figure 5(a), but for each eccentricity in order to give more details about each curve. From Figure 5(a) it is visible that the amplitude of oscillation of the curves increases with the eccentricity of the orbit. The combination of Figures 5(a), 5(b), 5(c), and 5(d) show that the increase of the eccentricity makes the values of the local maximum (inclination near 45 degrees) increase until a point (near e=0.2) where they become the global maximum. The global maxima of the circular orbits (equatorial and polar orbits) decrease in magnitude. When the eccentricity reaches the value 0.3, the figure shows only two maxima. The explanation for this phenomenon comes from the analysis already made regarding the physical reasons to have maximum at the polar, equatorial, and 45-degree inclined orbits. The polar and the equatorial orbits have their maximum based on the larger distances spacecraft cube, when compared to a sphere of the same mass, so the increase of the eccentricity that causes an approximation spacecraft cube during the passages of the spacecraft by the periapsis of the orbit works against the physical reason of the perturbation. In the opposite side, the 45-degree inclined orbits have their maximum based on the smaller distances spacecraft cube, when compared to a sphere of the same mass, so the increase of the eccentricity works in favor of the physical reason of the perturbation. So the net result is an increase of the PI for all orbits, because the spacecraft gets closer to the cube and it increases the perturbations, but this increase is different for each inclination. It is stronger for the 45-degree inclined orbit, causing an inversion of the large maximum. The minima are located in about the same positions, with a small shift to the left of the scale, which means that the increase of the eccentricity reduces the value of the inclination that generates the minimum perturbation. The physical reason is the same as explained above. The increase of the eccentricity makes the spacecraft pass closer to the cube, so the minimum perturbed orbits have smaller values for the inclination, to be far away from the edges and vertexes, that have their perturbations increased by the closer distances generated by larger eccentricities.(a) Perturbation integral versus inclination. Blue,e=0; red, e=0.1; green, e=0.2; orange, e=0.3. (b) Perturbation integral versus inclination for e=0.1. (c) Perturbation integral versus inclination for e=0.2. (d) Perturbation integral versus inclination for e=0.3.
(a)
(b)
(c)
(d)Another relevant analysis to be made is the influence of the semimajor axis of the orbit. Figure6 shows, for a circular equatorial orbit, the PI (m/s) as a function of the semimajor axis from 7.5 km to 25 km, for inclinations of 0, 45°, and 90°. The dependence of the perturbation relative to the semimajor axis is very clear, and when moving away from the body the value for PI decreases fast. The perturbation is larger for smaller distances, which is an expected result, because the perturbation in the potential equation for a cube changes at a rate of 1/r9. The interesting aspect of the results shown here is that they can quantify the values involved using a new approach.Figure 6
Perturbation integral versus semimajor axis fori=0=90° (blue line), and i=45° (red line), for a circular orbit.From the results presented here, there are two main things that need to be emphasized: first, the smaller the semimajor axis is, the larger the values of the PI are; second, as shown in Figure4, eccentric orbits have a higher value for the PI. Therefore, close eccentric orbits are the situations where the effects are larger, which means that the orbits are more perturbed.Besides the semimajor axis and the eccentricity, from the results obtained studying the inclination, it was proved that orbits passing near the vertex and edges of the cube are more affected and that there are maximum and minimum values for the perturbation as a function of the inclination. Another important consideration to take into account when studying nonspherical objects is that these objects are rotating. Therefore, the results for a rotating cube will be analyzed and showed next.


## 3.2. Adding Rotation

The results presented up to now considered that the cube is standing still. The next step is to consider that the cube is rotating around one of the principal axes (x, y, and z). Since the small celestial bodies are spinning, adding this fact will enrich the analysis of the gravitational perturbation due to the shape of the object.The rotation of an object may be obtained with a rotation matrix considering the direction of the rotation, which is usually represented by three matrices with respect to the axes of symmetry. If the object is rotating only in one principal axis, one of the rotation matrices is applied; therefore, knowing the direction in which the body is spinning is important. There can also be rotation in more than one axis, as already found in asteroids, and thus the sequence of the rotations must be chosen properly. The coordinates for rotation inx, y, and z axes, obtained by multiplying each rotation matrix by the position coordinates with respect to the orbital elements, are represented by the following equations:
(12)x′y′z′x=1-e2acos⁡ΩCos⁡v+o-cos⁡isinΩsinf+o1+ecos⁡f1-e2acos⁡ωtcos⁡f+osinΩ+cos⁡icos⁡Ωsinf+o1+ecos⁡f-1-e2asinisinf+osinωt1+ecos⁡f1-e2acos⁡ωtsinisinf+o1+ecos⁡f+1-e2acos⁡f+osinΩ+cos⁡icos⁡Ωsinf+osinωt1+ecos⁡f,x′y′z′y=1-e2acos⁡ωtcos⁡Ωcos⁡f+o-cos⁡isinΩsinf+o1+ecos⁡f+1-e2asinisinf+osinωt1+ecos⁡f1-e2acos⁡f+osinΩ+cos⁡icos⁡Ωsinf+o1+ecos⁡f1-e2acos⁡ωtsinisinf+o1+ecos⁡f-1-e2acos⁡Ωcos⁡f+o-cos⁡isinΩsinf+osinωt1+ecos⁡f,x′y′z′=-1-e2asinωtcos⁡o+fsinΩ+cos⁡icos⁡Ωsino+f1+ecos⁡f+1-e2acos⁡⁡ωtcos⁡Ωcos⁡o+f-cos⁡isinΩsino+f1+ecos⁡f1-e2acos⁡ωtcos⁡o+fsinΩ+cos⁡icos⁡Ωsino+f1+ecos⁡f+1-e2asinωtcos⁡Ωcos⁡o+f-cos⁡isinΩsino+f1+ecos⁡f1-e2asinisino+f1+ecos⁡f
where a is the semimajor axis, e is the eccentricity, i is the inclination, ω is the angular velocity, f is the true anomaly, o is the argument of periapsis, and Ω is the right ascension of ascending node.Several rotation periods were tested in order to compare the differences from a nonrotating situation, and among the rotations with different velocities. Reminding that the mass used was 4 × 1014 kg and the edge length of the cube is 5 km. The cube is assumed to have homogeneous density of 3.2 g/cm3, which is near the values found for the density of many asteroids [23, 24]. For all the examples that will be presented in the case of a rotating cube, the orbit is circular. Figure 7 gives the PI (m/s) as a function of the inclination for different angular velocities about the z axis.Figure 7
Perturbation integral versus inclination for different rotation periods around thez axis for circular orbits with semimajor axis of 10 km.The letterR in the figures stands for the rotation period, where R=0 means that the cube is not spinning. The physical reason to have different values when the cube is rotating is that the spacecraft passes by different parts of the cube. If it spends more time near the faces, the total perturbation increases. If it spends more time near the equilibrium position, which is a position where the potentials of a sphere and a cube are the same, the total value of the perturbation decreases. From Figure 7, it is noticeable that the results are not symmetric around the polar orbit, since prograde and retrograde orbits have different rotations with respect to the cubic body and so the regions where the subsatellite point passes are different. The equatorial orbits are the ones where the rotation generates smaller effects. The polar orbits have stronger effects when compared to the equatorial orbits, showing a decrease in the total perturbation effects. The largest reduction of the perturbation occurs near the inclination of 1 rad. In general, the rotations decreased the value of the PI. Retrograde orbits are different from direct orbits. Since they are rotating in the opposite direction of the cube, the edges and vertices spend more times by the spacecraft and it increases the perturbation. It is important to note the large reduction of the PI, going from near 0.000042 m/s for an equatorial orbit to near 0.000026 m/s in the minimum situation (inclination near 1 rad), a reduction of about 38%. The period of the rotation of the cube also has some effects in the perturbation, with the largest one simulated (24 hours) being the one that generates the smallest perturbations. Anyway, those differences are small when compared to the reduction generated by the inclusion of the rotation.Figure8 shows the perturbation integral as a function of the inclination for a cube rotating around the y axis. Each line color represents a different angular velocity, and the blue line shows the case of a fixed cube. The same representation is used in all similar figures. The figure is now symmetric with respect to the polar orbit and the minimum is near 0.55 rad. The reduction of the PI is again very large when compared to the fixed case, reaching a value near 33% at the minimum perturbation point. For the polar orbits, the perturbation is almost the same with and without rotation of the cube. When the period of rotation is 10 hours, the rotating situation generates even more perturbation. All those effects are a result of the averaged distance between the spacecraft and the mass distribution of the cube.Figure 8
Perturbation integral versus inclination for different rotation periods iny axis for circular orbits with semimajor axis of 10 km.The last case of rotation is in thex axis, showed in Figure 9. It is noticed that, now, some periods of rotation decrease too much the PI, generating orbits with lower perturbation. The results are now more dependent of the period of the rotation. The results have symmetry around the polar orbit again. In general, polar and equatorial orbits are more perturbed in the nonrotating case. When the body is rotating, the inclination of the orbit plays an important role in the stability of the orbit. For example, there are studies that suggest that retrograde orbits are more stable against the effects of perturbation around nonuniform gravitational fields [25–27]. Also, depending on the rotation axis and the period of rotation, the orbit could be synchronous or resonant. It is clear that the minimum perturbation situation changes with the inclination of the orbit, so the decision of the best orbit to place a satellite needs to take into account the period of rotation of the body.Figure 9
Perturbation integral versus inclination for different rotation periods around thex axis.The acceleration that the spacecraft feels due to a rotating cube is very dependent on the inclination of the orbit. Figure10 shows the acceleration for an orbit with different values for the inclination, with a rotation period of 5 hours around the z axis.Figure 10
Acceleration for one revolution, for different inclinations, with rotation around thez axis and period of 5 hours.Still considering a rotation period of 5 hours, Figure11 shows the acceleration for different inclinations for a cube with rotation about the y axis. In this case, the results for retrograde orbits are the same as for the prograde orbits, because the rotation axis in this case influences the orbits on the equatorial plane in the same way.Figure 11
Acceleration for one revolution, for different rotation inclinations with rotation around they axis, with rotation period of 5 hours.The last graph showing the acceleration for a rotating cube is around thex axis and can be seen in Figure 12. For this direction of rotation, the inclination of 45° is affected differently compared to the rest of the initial inclinations used.Figure 12
Acceleration for one revolution, for different rotation inclinations with rotation around thex axis, with rotation period of 5 hours.


## 3.3. Lagrange’s Planetary Equations

Considering that there is a gravitational perturbation, the variation of each orbital element can be obtained using Lagrange’s planetary equations. The disturbing potential of the cube with respect to a sphere is used, and six differential equations are obtained to express the evolution of the Keplerian elements of the orbit of the spacecraft. These equations are given by(13)dadt=2na∂V∂M,dedt=1-e2na2e∂V∂M-1-e2na2e∂V∂o,didt=-1na21-e2sin⁡i∂V∂Ω+cot⁡ina21-e2sin⁡i∂V∂o,dωdt=1-e2na2e∂V∂e-cot⁡ina21-e2sin⁡i∂V∂i,dΩdt=1na21-e2sin⁡i∂V∂i,dMdt=n-2na∂V∂a-1-e2na2e∂V∂e.The variation rate of the semimajor axis, eccentricity, and inclination, as a function of the semimajor axis, is shown in Figures13, 14, and 15, respectively. The results are for one revolution considering an equatorial orbit, a polar orbit, and an orbit with inclination of 45°. For the integration of the Lagrange equations, the magnitude was used to avoid compensations that leads to near zero integrals, but that still disturbs the orbits. This view of the mapping of the orbits is different from the view that comes from the integral approach, and it complements the study. The integral approach shows the magnitude of the total perturbation in the orbit of the spacecraft, independent of which orbital element is affected. The approach using Lagrange’s planetary equations studies the perturbation of each orbital element individually, so the relation between the results has to be seen under this perspective. In general, the graphs showed that the evolution of the orbital elements as functions of the semimajor axis, using Lagrange’s planetary equations, have a good correlation with the results presented for the perturbation integral method. It shows that the further from the body the spacecraft is, the smaller the change in the orbital elements is. It is noticed that there are some cases where the rate of variation is zero, like the polar orbit for the variation of the eccentricity and the equatorial orbit for the variation of the inclination. This point is not against the integral approach, which shows that the orbit is perturbed, because the rates of variation of the other elements are not zero. This situation just shows that the existing perturbation is disturbing another orbital element.Figure 13
d
a
/
d
t as a function of the semimajor axis around a cube.Figure 14
d
e
/
d
t as a function of the semimajor axis around a cube.Figure 15
d
i
/
d
t as a function of the semimajor axis around a cube.A complete comparison between both approaches is not possible due to the singularities of Lagrange’s planetary equations in the eccentricity and inclination. A future work can extend this comparison in further details by using nonsingular variables. One more example is possible and is shown in Figure16 that illustrates another example of the connection of the PI with the Lagrange equations, showing the variation rate of the semimajor axis as a function of the eccentricity. It shows a near zero value at the beginning, with only a slight motion for all the curves, increasing drastically after eccentricity 0.5. Similarly to what happens for the PI, the fast increase of the time derivatives of the Keplerian elements starts near e=0.5.Figure 16
d
a
/
d
t as a function of the eccentricity around a cube.


## 4. Conclusions

This work presented an analysis of the gravitational perturbation due to a nonspherical body. For this study, the closed form equation for the potential of a cube was considered, allowing working analytically most of the time. Using a methodology called the Perturbation Integral, or PI, it was possible to quantify the perturbation for the specific initial conditions tested. This integral gives the change in the velocity over the time and, therefore, it is an efficient way to search for the least perturbed orbits. An alternative tool to study the variation of the orbital elements was used, which consisted in the analysis of the rate of change of the orbital elements separately by using Lagrange’s planetary equations, and it was shown that both methods have a good correlation. Comparing the results coming from a sphere with the cube it was proved that, close to the body, the orbit acts different from a Keplerian orbit, and therefore the gravitational perturbation should be considered properly. Since small celestial bodies are rotating, it was also considered the case of a spinning cube, which was tested for rotation in the three principal axes.Using the equations for a cube is helpful to test the method because its shape has a predictable behavior and it could also be used to model an irregular body with the sum of cubes. In space missions, it is important to know the perturbations acting on the spacecraft, because it could displace the vehicle from the desired trajectory or even damage it. In the case of missions to small bodies like asteroids, a perturbation that must be considered is due to the irregularity of the shape of the body, especially very close to the body. Therefore, the search for the least perturbed and more stable orbits is an important topic.The results showed that circular orbits are less perturbed than eccentric orbits and they quantified the expected decrease of the perturbation with the increase of the semimajor axis. Regarding the inclination of the orbits, the existence of a minimum for some values of the inclination is shown, which is an indication that some orbits are less perturbed than others. This fact can be used for the orbit selection during a preliminary mission design. The results also show and explain the inversion of the global maximum of the PI with respect to the inclination when the eccentricity of the orbit increases, as well as the shift to lower inclinations of the respective minimum.The addition of the rotation of the cube showed interesting results, usually decreasing the perturbation of the orbit when compared to the nonrotating case. The inclination of the orbit that generates minimum perturbation changes with the rotation, so this fact should not be neglected.



---

*Source: 493903-2015-03-02.xml*