# Station Keeping of Constellations Using Multiobjective Strategies

Authors: Evandro Marconi Rocco; Marcelo Lopes de Oliveira e Souza; Antonio Fernando Bertachini de Almeida Prado
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

The goal of the present paper is to study the problem of station keeping maneuvers of satellites that belong to a constellation. The objective is to perform those maneuvers with low fuel consumption and but also including a time constraint. This type of problem has several aspects to be considered, like fuel consumption, duration of the maneuvers. These aspects are applied to all the satellites of the constellation, so a strategy has to be found to consider the global optimization of those variables that can control the geometry of the constellation. So, this is a multiobjective problem. To find a solution, strategies are formulated that minimize the fuel consumption considering all the satellites and also considering the duration of the maneuver. Numerical examples are shown.


---

## Body

## 1. Introduction

The majority of the spacecrafts that have been placed in orbit around Earth use the basic concepts of orbital transfers. During the launch, the spacecraft is placed in a parking orbit that can be different from the final orbit for which it was designed. Therefore, to reach the desired final orbit, the spacecraft must perform orbital maneuvers. Besides that, the spacecraft orbit must be corrected periodically because there are perturbations acting on it. This becomes more complex when the spacecraft is part of a constellation because in this case, the spacecraft position and velocity with respect to the other members are more important than the inertial position and velocity. The correction maneuvers can place the spacecraft in its nominal position or in a different one, in order to save fuel or time. Examples of applications of those orbital maneuvers to individual satellites can be the placement of a satellite in geosynchronous orbit, the maintenance of this orbit, the rendezvous between two spacecrafts, and so forth.Examples of applications of these maneuvers to constellations can be the maintenance made in the satellites of the constellations INMARSAT, GLOBALSTAR, and NAVSTAR GPS, which are constellations that require maneuvers to maintain the system operating. Besides the maintenance maneuvers, new satellites are sometimes injected and, at the same time, others are removed from those constellations, because the satellites have a limited lifetime. These placement and removal maneuver are periodic activities. Currently, in Brazil, there are studies to place a constellation of small satellites to help to observe and distribute Internet services in remote areas, as the interior of the Amazonian region, so those studies are of great interest. It is also possible to mention the cases of similar satellites in simultaneous operation. Although they do not form a strict constellation, the orbit maintenance of these satellites is very similar to the case of constellations, due to the similarities of the missions carried out by them, as it is the case of the Chinese-Brazilian CBERS series. These satellites should be maintained as if they formed a constellation, because they have similar orbital elements and they carry out the same type of mission. Therefore, it would be advantageous to maintain the satellites well distributed with respect to each other in order to avoid that two satellites cover the same point of the terrestrial surface at the same time, while other points remain uncovered.In this way, it is clear that there are two possible basic maneuvers: the transfer and the correction maneuvers. Both maneuvers are usually calculated with the goal of minimizing:(1) the fuel consumption, because it is not renewable in space; (2) the time spent in the maneuver, because it minimizes the energy consumption of the batteries, or because it is necessary to perform the maneuver in visibility. There are many studies in the literature that consider both goals individually, either assuming a low thrust propulsion system [1–12] or using an impulsive propulsion system [13–28].However, few of them consider both conditions of time and fuel consumption restrictions simultaneously. For constellations, they can correct the positions and/or relative or absolute velocities. They can make the corrections to place the spacecraft in its nominal values or in another value that provides a better strategy in terms of saving time and/or fuel, number of maneuvers required by the satellite, and the total duration and cost of the use of the control centers, ground tracking stations, and so forth. Examples of those studies are shown in the next parts of this paper.


## 2. Objective

The multiobjective optimization process consists of the minimization (or maximization) of a vector of objectivesZ(x) that may be subjected to constraints or bounds. This type of problem is very common in station keeping of satellite constellations, because during the maneuvers it is necessary to optimize many parameters at the same time. So, the objective of this work is to study the multiobjective optimization problem applied to the orbital maneuvers (transfers and corrections) of symmetrical constellations of satellites, with minimum fuel consumption, using impulsive maneuvers, and subjected to time constraints.To perform the station keeping of a constellation ofn satellites is a complex problem, because it is necessary to optimize the maneuvers of all the satellites at the same time, taking into account the fuel consumption as well as the time constraints. Therefore, the goal of this work is to formulate and study strategies that, in some way, make possible to obtain solutions with small fuel consumption considering all the satellites that form the constellation. When it is considered that all the satellites of the constellation should stay as close as possible to their nominal positions, it is necessary to maneuver them periodically to eliminate the drift caused by the disturbances. But the relative positions of the satellites are more important than the absolute positions, because, in most constellations, the main goal of the maneuvers is to keep the ground coverage, so the decision of when to perform the maneuvers depends on the analysis of the relative positions. In some cases it is considered, as a first approximation, that all the satellites (in the case of symmetrical constellations) suffer the same disturbances and, therefore, they present the same drift; so the relative positions among the satellites stay constant. Thus, the terrestrial coverage stays practically constant and, for the case of small drifts, it is not necessary to maneuver the satellites. On the other hand, if one of the satellites is maneuvered, the relative positions of this satellite with respect to the others change, so all the satellites have to be maneuvered. Besides that, when the maneuvers are performed, it may be more economical to place the satellites a little higher than their nominal orbits because, in this way, the interval of time between maneuvers is increased. Therefore, the decisions of when to maneuver the satellites and how to perform those maneuvers are very important and very complex when involving constellations.


## 3. Definition of the Problem

This work considers symmetrical constellations of satellites of the Rosette type for a first study. This type of constellation was developed independently by Walker [29] in Great Britain and by Mozhaev [30] in the former Soviet Union. It consists of a constellation using several orbital planes, all of them with the same inclination, distributed symmetrically along the Equator.This problem can be defined as the problem of finding a method that, when applied to a satellite constellation, is capable of keeping the satellites in their nominal positions, executing the maneuvers in a certain time, such that the fuel consumption of all the satellites in the constellation is minimized.Considering that the spacecraft propulsion system is able to apply an impulsive thrust and that the maneuver for each satellite comprises of two impulsive velocity increments with magnitudesΔv1 and Δv2, the total velocity increment for each satellite is given by
(1)V=Δv1+Δv2=F(X),
where X is an arbitrary variable for the transfer.The time spent in the maneuver for each satellite is given by(2)T=G(X).Therefore, the problem is to minimizeV for a prescribed value of  T. If this time is T0, the constraint relation is
(3)T-T0=0.In this way, the performance index for each satellite is given by(4)J=V+λ(T-T0).Considering then satellites of the constellation we have
(5)J1=V1+λ1(T-T1),J2=V2+λ2(T-T2),⋮Jn=Vn+λn(T-Tn).Therefore the performance index for the entire constellation is given by(6)J=∑i=1n(Vi+λi(T-Ti))=∑i=1nJi.This problem can be treated in an analytical and/or numerical approach, making use of numerical routines of integration and the solution of the Two Point Boundary Value Problem. However, due to the complexity of the performance index when considering all the satellites of the constellation, the optimal solution becomes quite difficult to obtain. The difficulty increases if constraints in the positions are considered, what certainly should be the case for constellations because the relative positions of the satellites should be considered. Those position constraints depend on the type of the constellation studied (Rosette [31], Polar [32], Ellipso [33], Polyhedron [34], etc.). Therefore, obtaining a suboptimal solution becomes more practical in many cases.


## 4. Presentation of the Method

A possible approach to solve this type of problem is the use of algorithms that search for suboptimal solutions considering the entire constellation. In the literature there are some works that use algorithms of this type: Shah et al. [32], Agnse and Brousse [33], Brousse et al. [34], Brochet et al. [35], Folta et al. [36], Kechichian [37], Carroll et al. [38], Krasilshikov and Sypalo [39], and so forth. However, in this work, the algorithm shown in Figure 1 is used. The functions carried out by the several boxes shown in Figure 1 are concisely described in Table 1.Table 1
Description of the flow diagram of Figure1.
(
1
)
Establishes contact between the satellite and the ground tracking station
(
2
)
Makes the position and velocity measurements, soX, Y, Z, X˙, Y˙, and Z˙ are determined in theinertial reference system
(
3
)
Verify if the satellite was already maneuvered during the same iteration of the control program,through the verification of the actual time of the satellite
(
4
)
If the actual time is smaller than the final time, the orbit of the satellite is determined by themeasurements of X, Y, Z, X˙, Y˙, and Z˙
(
5
)
It is verified which should be the orbital nominal elements of the satellite
(
6
)
The differences between the current elements and the nominal elements are calculated
(
7
)
The decision is made if the maneuver should be performed
(
8
)
It is defined which are the orbital elements that should be reached after the maneuver
(
9
)
It calculates the optimal transfer maneuver with time constraint
(
10
)
The maneuver is performed
(
11
)
The current orbital elements of the satellite are considered and we do not perform anycorrection maneuver
(
12
)
The position and velocity are found, which means thatX, Y, Z, X˙, Y˙, and Z˙ after themaneuver are determined
(
13
)
Propagation of the orbital elements for the next iteration
(
14
)
The processing time is transmitted to the satellite
(
15
)
The processing time is updated.Figure 1
Flow diagram.The problem can be divided in three main parts:(1) the choice of when and which satellites should be maneuvered; (2) the calculation of the optimal maneuver with time constraint for each satellite; (3) the propagation of the orbit of each satellite.The calculation box developed in Rocco [40] and Rocco et al. [41] is then used. The equations presented in Eckel and Vinh [42] are used, which provides the transfer orbit between noncoplanar elliptical orbits with minimum fuel and fixed time transfer or the minimum time transfer for prescribed fuel consumption. The case of minimum fuel and fixed time was studied in Rocco et al. [41] and the case of minimum time for prescribed fuel consumption was studied in Rocco et al. [43]. So, the calculation box is able to consider both cases. Besides that, still using the equations presented by Eckel and Vinh [42], the case of the coplanar orbital transfers were developed and implemented.The main contribution of the present work is the development of the decision box. In this decision box the satellites that will be maneuvered must be selected. To do this, it should be verified if the differences between the nominal elements and the actual elements are inside the tolerance allowed for the mission and if these differences are increasing or decreasing along time. Besides that, it should be verified that the fuel consumption and the time spent in the maneuver for each satellite, as well as for all the satellites will be maneuvered. After that, it is necessary to know how these maneuvers will affect the mission (coverage area, interruption of the operation of the satellites, etc.). To determine the fuel consumption, calculations have to be made to determine the optimal transfer maneuver for each satellite. In this way, a feedback is necessary from the decision box after the optimal transfer maneuvers have been obtained in the calculation box. This feedback is necessary because there are situations where, even punishing the fuel consumption or the time spent with a certain satellite, it is possible to obtain an economy when considering the whole constellation. Therefore, the decision box should have the capacity to compare all possible combinations of maneuvers and to choose one which is feasible and that minimizes the fuel consumption of the entire constellation. So, it is a multiobjective problem.According to Cohon [44], a static optimization problem with one objective can be written as
(7)MaximizeZ(x)Subjecttogi(x)≤0i=1,2,…,mSubjecttoax≥0GivenZ(·).
or
(8)MaximizeZ(x)Subjecttox∈FdGivenZ(·),Fd,
where Fd is the practicable area of the decision space that is defined by
(9)Fd={x∣gi(x)≤0,i=1,2,…,m;x≥0}.A multiobjective optimization problem can be written as(10)MaximizeZ(x)=[Z1(x),Z2(x),…,Zp(x)]Subjecttox∈Fd.
Therefore, the objective function, in this case, is a vector with dimension p.In problems of one-dimensional optimization (when there is only one objective), the possible solutions (x∈Fd) can be compared by means of the objective function, that is, given two solutions x1 and x2 it is possible to compare Z(x1) with Z(x2) and then determine the optimal solution x* such that it does not exist x∈Fd that implies in Z(x)>Z(x*).In problems of multidimensional optimization (multiobjective problems), in general, it is not possible to compare all the possible solutions, because the comparison on the basis of one objective can conflict with the comparison based on another objective. Namely, supposing that(11)Z(x1)=[Z1(x1),Z2(x1)],Z(x2)=[Z1(x2),Z2(x2)],x1 is better than x2 if and only if
(12)Z1(x1)>Z1(x2),Z2(x1)≥Z2(x2),
or
(13)Z1(x1)≥Z1(x2),Z2(x1)>Z2(x2).IfZ1(x1)>Z1(x2) and Z2(x1)<Z2(x2) it is not possible to conclude anything regarding x1 and x2.A possible solution for the multiobjective problem would be to combine all the objectives to obtain only one final objective that is formed by the average of the original objectives multiplied by an influence factor. This eliminates the necessity of the use of a more complex optimization algorithm, but introduces new parameters, that are the influence factors. Thus, the solution depends on the correct determination of these factors. Therefore, this determination process becomes an optimization by itself.Another possibility would be the use of the algorithm developed by Pareto (Carroll et al. [38]). This algorithm can be used in problems where the objectives compete with each other. The optimal solution is an element of the group of solutions that is considered equally good with respect to all the components of the objective vector. So, the algorithm should select, starting from the practicable area of the decision space, a group of solutions that supports one “a priori” defined criterion. This can be made through a methodology that systematically makes a comparison among candidates. A solution x can only be considered optimal for a certain group of objectives, if a better solution y, considering all the objectives, does not exist. A solution x is called “nondominated” if a possible other solution y does not exist such that:
(14)Z(y)≥Z(x),
or
(15)Zk(y)≥Zk(x),k=1,2,…,p.Ify exists, then x is a dominated solution (or inferior) and so it is not considered a Pareto's optimal solution, which means that it is not the best solution in none of the objectives.According to Kuhn and Tucker [45], if x is a nondominated solution, then multipliers ui≥0,i=1,2,…,m  and  wk≥0,k=1,2,…,p should exist such that
(16)x∈Fd,uigi(x)=0,i=1,2,…,m,∑k=1pwk∇Zk(x)-∑i=1mui∇gi(x)=0.The conditions of Kuhn-Tucker, (16) are the necessary conditions for the solution x to be nondominated. They are also sufficient if Zk(x) are concave for k=1,2,…,p,Fd is convex, and wk>0 for all k.However, it is noticed that there are different domain levels. A dominated solution is always inferior with respect to some nondominated solution, but a dominated solution can be dominated by other dominated solution. Those domain levels furnish a characterization of the practicable area of the decision space, sorting the solutions in categories with different levels of optimality. Using this concept, it is clear that the Pareto's optimization process can be described basically as being a search for nondominated solutions. This search consists of sorting candidate solutions in groups ofdominatedand nondominated. So, the nondominated solutions form a set of solutions for the multiobjective optimization problem.To illustrate this methodology, in Tables2 and 3, a hypothetical example is presented. It is desired to choose a maneuver that minimizes the Δθ (rad), the position constraint; Δv (km/s), the velocity increment necessary to perform the maneuver; the time spent T in seconds.Table 2
Optimal maneuvers.
Δ
θ
Δ
v
T
(
1
)
0.565637
1.07333
1660
(
2
)
0.445267
1.02667
1700
(
3
)
0.314329
0.99600
1800
(
4
)
0.118357
0.88733
2110
(
5
)
0.217551
0.85493
2250
(
6
)
0.241902
0.84267
2305
(
7
)
0.275673
0.80400
2550
(
8
)
0.278017
0.79867
2600
(
9
)
0.291305
0.78667
2705
(
10
)
0.289064
0.77867
2800
(
11
)
0.314594
0.76400
2910
(
12
)
0.337307
0.76212
2990Table 3
Optimal maneuvers using nondimensional units.
Δ
θ
/
Δ
θ
max
Δ
v
/
Δ
v
max
T
/
T
max
(
1
)
0.808053
0.894442
0.533333
(
2
)
0.636096
0.855558
0.566667
(
3
)
0.449041
0.830000
0.60000
(
4
)
0.169081
0.739442
0.703333
(
5
)
0.310787
0.712442
0.75000
(
6
)
0.345574
0.702225
0.768333
(
7
)
0.393819
0.670000
0.850000
(
8
)
0.397167
0.665558
0.866667
(
9
)
0.416150
0.655558
0.901667
(
10
)
0.412949
0.648892
0.933333
(
11
)
0.449420
0.636667
0.970000
(
12
)
0.481867
0.635100
0.996667Examining Table2 it is possible to select maneuvers 4, 12, and 1 as nondominated solutions for the entire group of candidates. However, it is possible to continue the sorting of the solutions to obtain other groups with different levels of Pareto's optimality: level 1, maneuvers 4, 12, and 1; level 2, maneuvers 5, 11, and 2; level 3, maneuvers 6, 10, and 3; level 4: maneuvers 7 and 9; level 5: maneuver 8. With this classification it is possible to choose the best maneuver that satisfies the constraints. According to Carroll et al. [38], it is possible to choose any of the maneuvers of the optimality level 1 that would obtain the highest level of Pareto's optimality in the satisfaction of the multiobjective problem. But, which maneuver should be chosen? All of them are good solutions for a certain objective. Choosing any of those solutions would also imply in choosing a certain objective as a priority and, in this way, the problem of the choice of influence factors for the objectives, which is an optimization process by itself, would be back. In this situation the Pareto's optimization algorithm becomes unnecessary.However, if points in a three-dimensional graph that represents the maneuvers selected in each optimality levels are made, a graph is obtained, for three objectives, showing triangles whose vertices represent the optimal solutions for a certain objective. Figures2 and 3 show this fact.Figure 2
Solutions of levels 1 and 2.Figure 3
Solutions of levels 1, 2, and 3.These solutions are the extreme ones, because each one of them only considers one objective as a priority. In this way, to obtain a less extreme solution that does not assume any objective as a priority, it would be necessary to choose an intermediary solution.An attempt to find this intermediate solution would be to find the barycenter of the triangle. It seems that the optimal solution for the multiobjective problem would be the central point of the figure that has the optimal solutions for each objective as vertices. Therefore, for problems with three objectives, the solution would be the center of a triangle and, forn objectives, the solution would be the center of a n-dimensional figure. This is known as the “Smallest Loss Criterion”. Since it is not possible to satisfy all the objectives simultaneously, the best multiobjective solution is the one that implies in the smallest loss for the objectives when considered independently (Rocco et al. [46]).Therefore, for the example shown here, the solution, considering the triangle of level 1, isΔθ = 0.340434 rad; Δv = 0.907593 km/s; t = 2233.333 s. However, this maneuver does not exist in Table 2. So, the best maneuver seems to be that whose parameters are the one, that is, closest to the coordinates of the barycenter. Calculating the distances, using nondimensional units, as shown in Table 3, between the barycenter and the points of the maneuvers, it is concluded that the best maneuver to be applied is maneuver 6. Another possibility would be to use the program that calculates the optimal maneuvers again with time 2233.333 s as input, and then calculate a maneuver that approaches the barycenter coordinates.To exemplify this methodology applied to a constellation of satellites, an example of a constellation composed by three satellites with circular and equatorial nominal orbits is used. Therefore, the orbital elements are given by (distances in km, velocities in km/s and angles in radians)(17)e=0.00000000,a=7100.00000000,l=7100.00000000,i=0.00000000,ω=0.00000000,Ω=0.00000000,
where e is the eccentricity, l is the semilatus rectum, ω is the argument of perigee, a is the semimajor axis, i is the inclination, and Ω is the right ascension of the ascending node.It is considered that, to attend the specifications of the mission, the satellites should be positioned such that the differences among the true longitudes of the satellites (θ1, θ2, and θ3) should be equal to 120°.At the initial instant, the satellite 1 is entering in the cone of visibility of the ground tracking station and the measurements ofX, Y, Z, X˙, Y˙, Z˙ in the inertial reference system are shown below:
(18)X=7000.00000000,X˙=0.00000000,Y=0.00000000,Y˙=7.54605517,Z=0.00000000,Z˙=0.00000000,Thus, it is possible to determine the actual elements of the satellite 1 and calculate the differences between the nominal orbit and the actual orbit:(19)a1=7000.00000000km,Ω1=0.00000000,u1=0.00000000,e1=0.00000000,ω1=0.00000000,f1=0.00000000,i1=0.00000000,M1=0.00000000,θ1=0.00000000,whereM is the mean anomaly, u is the eccentric anomaly, f is the true anomaly, and θ is the true longitude (θ=ω+f).With the orbit propagation it is possible to obtain the actual elements for the other satellites, using the orbital elements determined in the last passage of each satellite by the ground tracking station.Satellite 2:(20)a2=7050.00000000km,Ω2=0.00000000,u2=4.10152228rad,e2=0.00000000,ω2=0.00000000,f2=4.10152228rad,i2=0.00000000,M2=4.10152228rad,θ2=4.10152228rad.Satellite 3:(21)a3=7100.00000000km,Ω3=0.00000000,u3=6.19591663rad,e3=0.00000000,ω3=0.00000000,f3=6.19591663rad,i3=0.00000000,M3=6.19591663rad,θ3=6.19591663rad.Figure4 shows the nominal positions of the satellites. With the actual true longitudes ofθ1, θ2, and θ3 the position constraints Δθ1, Δθ2, Δθ3, and Δθ are calculated:
(22)θ1≤θ2:Δθ1=(θ2-θ1)-2π3,θ2<θ1:Δθ1=(2π-θ1+θ2)-2π3,θ2≤θ3:Δθ2=(θ3-θ2)-2π3,θ3<θ2:Δθ2=(2π-θ2+θ3)-2π3,θ3≤θ1:Δθ3=(θ1-θ3)-2π3,θ1<θ3:Δθ3=(2π-θ3+θ1)-2π3,Δθ=|Δθ1|+|Δθ2|+|Δθ3|3.Figure 4
Nominal positions of the satellites.If the difference between the nominal and the actual elements and the position constraintΔθ do not satisfy a tolerance previously specified, at least a correction maneuver becomes necessary; otherwise it is considered that the actual elements and propagation of the orbit of the satellite 1 is made to predict which will be the elements of this satellite in the next passage by the ground tracking station.Considering that it is necessary to execute the maneuver, several possible maneuvers are calculated, each one of them with different values of the semimajor axis of the final orbit and different values of the time spend by the maneuver. Thus, the values of the orbital elements of the transfer orbit are obtained, whereα1 is the true anomaly of the location of the first impulse, α2 is the true anomaly of the location of the second impulse, Δv1 and Δv2 are the velocity increments generated by the first and second impulses, Δv is the total velocity increment, and T is the time spent in the maneuver.Maneuver 1:afinal=anominal= 7100 km,
(23)a=7058.47757275,α1=3.15689637,Δv2=0.24778731,e=0.03350302,α2=3.58355988,Δv=0.49527448,ω=1.80337175,Δv1=0.24748417,T=400.00000000.Maneuver 2:afinal=1.01anominal= 7171 km,
(24)a=3665.74352360,α1=3.39673299,Δv2=6.40939972,e=0.97703009,α2=3.79335480,Δv=13.08244661,ω=0.29589943,Δv1=6.67304689,T=400.00364957.Maneuver 3:afinal=1.02anominal= 7242 km,
(25)a=7186.41901289,α1=2.75652874,Δv2=0.68732022,e=0.09296305,α2=3.12445310,Δv=1.37312307,ω=1.37794679,Δv1=0.68580285,T=350.00000000.Maneuver 4:afinal=1.03anominal= 7313 km,
(26)a=7213.51705205,α1=2.25536227,Δv2=0.61884640,e=0.08487354,α2=2.77724241,Δv=1.23430207,ω=0.95995823,Δv1=0.61545567,T=500.0000000.Maneuver 5:afinal=1.04anominal= 7384 km,
(27)a=7362.14584995,α1=2.28558873,Δv2=1.09147041,e=0.14827650,α2=2.64836455,Δv=2.17926114,ω=0.90897257,Δv1=1.08779073,T=350.00000000.Maneuver 6:afinal=1.05anominal= 7455 km,
(28)a=7377.96016234,α1=2.24060920,Δv2=1.00299979,e=0.13749403,α2=2.70383841,Δv=1.99966828,ω=0.92020617,Δv1=0.99666849,T=450.00000000.Maneuver 7:afinal=1.06anominal= 7526 km,
(29)a=7467.28224595,α1=2.25598085,Δv2=1.16017752,e=0.15936377,α2=2.71609123,Δv=2.31204450,ω=0.93670814,Δv1=1.15186698,T=450.00000000.Maneuver 8:afinal=1.07anominal= 7597 km,
(30)a=7566.14143619,e=0.18135158,ω=0.87937375,α1=2.19754033,Δv1=1.30699181,Δv=2.62449979,α2=2.65460419,Δv2=1.31750798,T=450.00000000.Therefore, Tables4 and 5 are obtained, considering  Δθmax = 0.7 rad, Δvmax = 3.0 km/s, and Tmax = 600 s.Table 4
Optimal maneuvers.
Δ
θ
Δ
v
Time
(
1
)
0.76883425
0.49527148
400
(
2
)
0.90869753
13.08244661
400
(
3
)
0.49794018
1.37312307
350
(
4
)
0.16074759
1.23430207
500
(
5
)
0.18054781
2.17926114
350
(
6
)
0.14717616
1.99966828
450
(
7
)
0.15534470
2.31204450
450
(
8
)
0.11435334
2.62449979
450Table 5
Optimal maneuvers using nondimensional units.
Δ
θ
/
Δ
θ
max
Δ
v
/
Δ
v
max
T
/
T
max
(
1
)
1.09833464
0.16509049
0.6666667
(
2
)
1.29813933
4.36081554
0.6666667
(
3
)
0.71134311
0.45770769
0.5833333
(
4
)
0.22963941
0.41143402
0.8333333
(
5
)
0.25792544
0.72642038
0.5833333
(
6
)
0.21025166
0.66655609
0.7500000
(
7
)
0.22192100
0.77068150
0.7500000
(
8
)
0.16336191
0.87483326
0.7500000Examining Table4 it is possible to select maneuvers 8, 4, 3, and 5 as nondominated solutions for the entire group of candidates (level 1). The maneuvers of level 2 are numbers 6 and 7. Maneuvers 1 and 2 were not considered because they exceeded the values of Δθmax and Δvmax.The barycenter coordinates of level 1, using the nondimensional units shown in Table5, are 0.34056747; 0.61759884; 0.6875000 which are equivalent to Δθ = 0.23839723 rad; Δv = 1.85279652 km/s; Time = 412.5 s.Calculating the distances between the barycenter and the points determined by the maneuvers of level 1, the results are(31)d(bar,m8)=0.31855550,d(bar,m3)=0.41700174,d(bar,m4)=0.27581938,d(bar,m5)=0.17182121.Using the program to calculate the optimal maneuver with time 412.5 s as input, the following maneuver is obtained:(32)a=7614.92854399,α1=2.25471058,Δv2=1.43892611,e=0.19756377,α2=2.67368856,Δv=2.86809986,ω=0.91571645,Δv1=1.42917375,T=412.50000000,
which results in a position constraint Δθ/Δθmax = 0.21922727, a velocity increment Δv/Δvmax = 0.95603329, and a time spent in the maneuver of T/Tmax = 0.6875. Thus, the distance among the barycenter and the point determined by the maneuver is d(bar,m) = 0.35952930.Comparing all the distances calculated, it is noticed that the maneuver that better approaches the coordinates of the barycenter is maneuver 5. In this case, even calculating the transfer maneuver again, using as input the time specified by the barycenter, it was not possible to obtain a better maneuver. In this way, maneuver 5 is the best maneuver.Therefore, in the diagram presented in Figure1, maneuver 5 is selected by the decision box. Then, a command is sent to the maneuvering box in order to implement the solution found. The process is repeated for each satellite as they come into visibility by the ground tracking station.Instead of using a multiobjective methodology, each satellite could be maneuvered to its nominal orbit, or to another specified orbit, whenever the position constraints necessary to maintain the constellation were not respected. This is what normally occurs in station keeping of constellations. A comparison of these two methods can be found in Rocco et al. [47–49].Another example to illustrate the advantage of the multiobjective approach is the orbital maintenance of a constellation in Low Earth Orbit. In this situation, the effect of the atmospheric drag is important. Thus, the variation of the semimajor axis is accentuated and the software developed is submitted to a situation of worst case with respect to the orbital maintenance of the satellites. Three satellites in circular and equatorial nominal orbits compose the constellation studied. In this simulation all satellites will be maneuvered. Therefore, the nominal elements of the satellites are given by (a and l in km and angles in radians)
(33)e=0.00000000,a=6603.13900000,l=6603.13900000,i=0.00000000,ω=0.00000000,Ω=0.00000000.It may be assumed that, at the initial instant, satellite 1 is entering in the visibility cone of the ground tracking station and that the actual orbital elements of this satellite are given by the following.Satellite 1:(34)a1=6601.93044509,Ω1=0.00000000,u1=0.00000000,e1=0.00000000,ω1=0.00000000,f1=0.00000000,i1=0.00000000,M1=0.00000000,θ1=0.00000000,
where M is the mean anomaly, u is the eccentric anomaly, f it is the true anomaly, and θ is the true longitude. The actual orbital elements of the others satellites are given by the following.Satellite 2:(35)a2=6601.93044509,Ω2=0.00000000,u2=2.09439522,e2=0.00000000,ω2=0.00000000,f2=2.09439522,i2=0.00000000,M2=2.09439522,θ2=2.09439522.Satellite 3:(36)a3=6601.93044509,Ω3=0.00000000,u3=4.18879042,e3=0.00000000,ω3=0.00000000,f3=4.18879042,i3=0.00000000,M3=4.18879042,θ3=4.18879042.To assist the specifications of the mission, the satellites should be positioned in such a way that the difference between the true longitudes (Δθ1, Δθ2, and Δθ3) should be 2.09439435 rad (120 degrees). With the actual true longitudes θ1, θ2, and θ3 it is calculated the position constraints Δθ1, Δθ2, Δθ3, and Δθ, which represent the position errors of the satellites, in the same way done in the previous example.If the differences between the nominal and the actual elements and the position constraintΔθ do not satisfy the tolerance previously specified, at least one correction maneuver becomes necessary. The actual position errors for the satellites are Δθ1 = Δθ2 = Δθ3 = Δθ = 0, but the error in the semimajor axis is −1.20855491 km. Considering that it is necessary to execute the maneuver, several possible maneuvers are calculated, each one of them with different values of the semimajor axis for the final orbit and different values of the time spent by the maneuver. The semimajor axis and the time vary from predefined values that belong to an operation range for the satellite. Thus, the orbital elements of the transfer orbit are obtained, where Δv is the total velocity increment and T is the time spent in the maneuver. Note that the semimajor axis of the transfer orbit is much smaller than the final orbit because its eccentricity is larger.Maneuver 1:afinal=anominal= 6603.139 km,
(37)a=3813.50753965,e=0.84297654,ω=3.29684376,Δv=10.43563617,T=701.00770076,Maneuver 2:afinal=1.01anominal= 6603.364 km,
(38)a=3838.94341245,e=0.82943908,ω=3.30375052,Δv=10.16880281,T=702.00882124,Maneuver 3:afinal=1.02anominal= 6602.914 km
(39)a=3832.39271231,e=0.82972598,ω=3.30147425,Δv=10.15092314,T=692.00916737,Maneuver 4:afinal=1.03anominal= 6603.589 km,
(40)a=3826.99885076,e=0.83567624,ω=3.30053543,Δv=10.28923606,T=701.10831091,Maneuver 5:afinal=1.04anominal= 6602.689 km,
(41)a=3832.00577843,e=0.83006536,ω=3.30139725,Δv=10.15886021,T=692.50914418,Maneuver 6:afinal=1.05anominal= 6603.814 km,
(42)a=3836.77670075,e=0.83146479,ω=3.30331350,Δv=10.21398945,T=704.45390971,Maneuver 7:afinal=1.06anominal= 6602.464 km,
(43)a=3836.67153417,e=0.82958744,ω=3.30024999,Δv=10.16508475,T=699.00894038,Maneuver 8:afinal=1.07anominal= 6604.039 km,
(44)a=3834.28236173,e=0.83205939,ω=3.30253352,Δv=10.21971060,T=702.00859016,Maneuver 9:afinal=0.996anominal= 6602.239 km,
(45)a=3834.83954066,e=0.82894987,ω=3.30220747,Δv=10.14092168,T=694.00916232,Applying the “Smallest Loss Criterion” to select the best maneuver, using normalized measurements withΔθmax=0.7 rad, Δvmax=11km/s, and Tmax=1000s to calculate the barycenter, Table 6 can be obtained.Table 6
Choise of the best maneuver.
Δ
θ
Δ
v
Time
Classification
(
1
)
0.34398092
9°
10.43563617
9°
701.00770076
5°
9°
(
2
)
0.33559710
5°
10.16880281
5°
702.00882124
8°
5°
(
3
)
0.33072576
1°
10.15092314
2°
692.00916737
1°
1°
(
4
)
0.33920262
8°
10.28923606
8°
701.10831091
6°
8°
(
5
)
0.33118940
2°
10.15886021
3°
692.50914418
2°
2°
(
6
)
0.33815711
7°
10.21398945
6°
704.45390971
9°
7°
(
7
)
0.33415824
4°
10.16508475
4°
699.00894038
4°
4°
(
8
)
0.33730895
6°
10.21971060
7°
702.00859016
7°
6°
(
9
)
0.33122408
3°
10.14092168
1°
694.00916232
3°
3°
B
0.33122408
10.14092168
694.00916232Therefore, maneuver 3 is selected, because it is the maneuver closest to the barycenter. Now, it is considered that the satellite 3 is in visibility (θ3=0). Therefore, the satellites present the following orbital elements.Satellite 1:(46)a1=6602.68166455,Ω1=0.00000000,u1=1.59913823,e1=0.00000000,ω1=0.00000000,f1=1.59913823,i1=0.00000000,M1=1.59913823,θ1=1.59913823.Satellite 2:(47)a2=6601.36490371,Ω2=0.00000000,u2=4.18996899,e2=0.00000000,ω2=0.00000000,f2=4.18996899,i2=0.00000000,M2=4.18996899,θ2=4.18996899.Satellite 3:(48)a3=6601.36490371,Ω3=0.00000000,u3=0.00000000,e3=0.00000000,ω3=0.00000000,f3=0.00000000,i3=0.00000000,M3=0.00000000,θ3=0.00000000.
So, the actual position errors for the satellites are Δθ1 = −0.49643566; Δθ2 = 0.00117879; Δθ3 = 0.49525687;Δθ = 0.33095711. Calculating the maneuvers for the satellite 3, the following results are obtained.Maneuver 1:afinal=anominal= 6603.139 km,
(49)a=3820.10647821,e=0.83967098,ω=3.29880702,Δv=10.37217642,T=702.20795943,Maneuver 2:afinal=1.01anominal= 6603.364 km,
(50)a=3838.33164219,e=0.82962949,ω=3.30368228,Δv=10.17237866,T=701.80880615,Maneuver 3:afinal=1.02anominal= 6602.914 km,
(51)a=3820.02868970,e=0.83647339,ω=3.29834695,Δv=10.28484948,T=692.32301320,Maneuver 4:afinal=1.03anominal= 6603.589 km,
(52)a=3831.09118064,e=0.83285319,ω=3.30165542v,Δv=10.22995622,T=699.50860888,Maneuver 5:afinal=1.04anominal= 6602.689 km,(53)a=3829.89510389,e=0.83097526,ω=3.30092382,Δv=10.17578723T=692.00905341,Maneuver 6:afinal=1.05anominal= 6603.814 km,
(54)a=3801.74263575,e=0.85048278,ω=3.29372503,Δv=10.59662285,T=704.00701406,Maneuver 7:afinal=1.06anominal= 6602.464 km,
(55)a=3821.94279967,e=0.83753216,ω=3.29916105,Δv=10.32212647,T=699.00822263,Maneuver 8:afinal=1.07anominal= 6604.039 km,
(56)a=3827.14559359,e=0.83613392,ω=3.30076197,Δv=10.30244922,T=702.80820995,Maneuver 9:afinal=0.996anominal= 6602.239 km,
(57)a=3821.52887093,e=0.83498276,ω=3.29863995,Δv=10.25135796,T=690.50874504,Applying the “Smallest Loss Criterion” to select the best maneuver, using normalized measurements withΔθmax=0.7rad, Δvmax=11km/s, and Tmax=1000s to calculate the barycenter, Table 7 can be obtained.Table 7
Choise of the best maneuver.
Δ
θ
Δ
v
Time
Classification
(
1
)
0.34323784
8°
10.37217642
8°
702.20795943
7°
8°
(
2
)
0.33646255
4°
10.17237866
1°
701.80880615
6°
4°
(
3
)
0.33606968
3°
10.28484948
5°
692.32301320
3°
5°
(
4
)
0.33739004
5°
10.22995622
3°
699.50860888
5°
3°
(
5
)
0.33236134
1°
10.17578723
2°
692.00905341
2°
1°
(
6
)
0.35149893
9°
10.59662285
9°
704.00701406
9°
9°
(
7
)
0.34016808
6°
10.32212647
7°
699.00822263
4°
6°
(
8
)
0.34122538
7°
10.30244922
6°
702.80820995
8°
7°
(
9
)
0.33416546
2°
10.25135796
4°
690.50874504
1°
2°
B
0.33432978
10.19984128
694.77553487Therefore, maneuver 5 is selected. Now, it is considered that the satellite 2 is in visibility (θ2=0). The satellites have the following orbital elements.Satellite 1:(58)a1=6602.25313837,Ω1=0.00000000,u1=3.69287951,e1=0.00000000,ω1=0.00000000,f1=3.69287951,i1=0.00000000,M1=3.69287951,θ1=3.69287951.Satellite 2:(59)a2=6600.67050929,Ω2=0.00000000,u2=0.00000000,e2=0.00000000,ω2=0.00000000,f2=0.00000000,i2=0.00000000,M2=0.00000000,θ2=0.00000000.Satellite 3:(60)a3=6602.43722227,Ω3=0.00000000,u3=1.59657315,e3=0.00000000,ω3=0.00000000,f3=1.59657315v,i3=0.00000000,M3=1.59657315,θ3=1.59657315.Thus, the actual position errors for the satellites areΔθ1 = 0.32893532; Δθ2 = 0.32240998;Δθ3 = 0.00652534; Δθ = 0.21929022. Calculating the maneuvers for the satellite 2, the results are as follows.Maneuver 1:afinal=anominal= 6603.139 km,
(61)a=3772.37771504,e=0.86553757,ω=3.28511293,Δv=10.90875480,T=700.50599655,Maneuver 2:afinal=1.01anominal= 6603.364 km,
(62)a=3823.66473495,e=0.83557977,ω=3.29961906,Δv=10.27556023,T=695.80850117,Maneuver 3:afinal=1.02anominal= 6602.914 km,
(63)a=3798.81533139,e=0.84792405,ω=3.29267829,Δv=10.51590970,T=692.00758398,Maneuver 4:afinal=1.03anominal= 6603.589 km,
(64)a=3813.45483938,e=0.84224231,ω=3.29699006,Δv=10.41646132,T=699.00782878,Maneuver 5:afinal=1.04anominal= 6602.689 km,
(65)a=3815.40470459,e=0.84007344,ω=3.29739101,Δv=10.36610269,T=696.00810762,Maneuver 6:afinal=1.05anominal= 6603.814 km,
(66)a=3822.87235742,e=0.83877470,ω=3.29977568,Δv=10.35876729,T=704.10794359,Maneuver 7:afinal=1.06anominal= 6602.464 km,
(67)a=3834.51473078,e=0.82946812,ω=3.30244565,Δv=10.15517129,T=695.50905470,Maneuver 8:afinal=1.07anominal= 6604.039 km,
(68)a=3832.27494773,e=0.83282349,ω=3.30220482,Δv=10.23427641,T=701.50853489,Maneuver 9:afinal=0.996anominal= 6602.239 km,
(69)a=3831.52178764,e=0.83010164,ω=3.30149199,Δv=10.16054631,T=692.50911115,Applying the “Smallest Loss Criterion” to select the best maneuver, using normalized measurements withΔθmax=0.7rad, Δvmax=11km/s, and Tmax=1000s to calculate the barycenter, Table 8 can be obtained.Table 8
Choise of the best maneuver.
Δ
θ
Δ
v
Time
Classification
(
1
)
0.02871318
9°
10.90875480
9°
700.50599655
7°
9°
(
2
)
0.00575679
3°
10.27556023
4°
695.80850117
4°
1°
(
3
)
0.01195674
7°
10.51590970
8°
692.00758398
1°
8°
(
4
)
0.01179272
6°
10.41646132
7°
699.00782878
6°
6°
(
5
)
0.00879068
5°
10.36610269
6°
696.00810762
5°
3°
(
6
)
0.01514568
8°
10.35876729
5°
704.10794359
9°
7°
(
7
)
0.00164066
2°
10.15517129
1°
695.50905470
3°
5°
(
8
)
0.00687667
4°
10.23427641
3°
701.50853489
8°
2°
(
9
)
0.0012963
1°
10.16054631
2°
692.50911115
2°
4°
B
0.00496478
10.27720910
693.34191661Therefore, maneuver 2 is selected.The comparison among the results obtained with the application of the multiobjective methodology developed in this work (“Criterion of the Smallest Loss” and the results obtained with the application of nominal maneuvers (each satellite is maneuvered to reach the nominal orbit), can be seen in Figures5 to 7. Figure 5 compares the position error after maneuvering the three satellites of the constellation. Figure 6 compares the sum of the velocity increment necessary to maneuver the satellites of the constellation. Figure 7 compares the sums of the time spent to maneuver all the satellites.Figure 5
Δ
θ
MULTIOBJECTIVE and ΔθNOMINAL.Figure 6
Δ
v
MULTIOBJECTIVE and ΔvNOMINAL.Figure 7
Δ
t
MULTIOBJECTIVE and ΔtNOMINAL.Figures5 to 7 show that the multiobjective methodology presented better results. The position errors obtained by the multiobjective methodology were considerably smaller than the error obtained by the nominal methodology. The necessary velocity increment and the time spent to maneuver the satellites were also smaller when the multiobjective methodology was applied. Therefore the multiobjective methodology applied to the orbital maintenance of the satellite constellation proved to be more advantageous than a methodology where each satellite is maneuvered into its nominal orbit. With the multiobjective approach and with the use of the “Smallest Loss Criterion,” it was possible to obtain better results for all individual objectives considered.


## 5. Conclusions

The present paper studied the problem of orbital maneuvers of satellites that forms a constellation. The main objective was to find maneuvers that considered the time required by the transfer as well as the fuel consumption. So, a strategy using the Pareto's idea of multiobjective optimization was used to solve the problem for a global optimization of the constellation.As an example, a simple constellation composed by three satellites was considered to test the algorithm developed. The tests have shown that the algorithm reaches the objectives, been able to consider all the constraints at the same time.



---

*Source: 476451-2013-02-19.xml*