International Journal of Aerospace Engineering
Volume 2017 (2017), Article ID 1235692, 17 pages
https://doi.org/10.1155/2017/1235692

Research Article

Satellite Constellation Orbit Design Optimization with Combined Genetic Algorithm and Semianalytical Approach

Tania Savitri,¹ Youngjoo Kim,¹ Sujang Jo,² and Hyochoong Bang¹

¹Department of Aerospace Engineering, Korea Advanced Institute of Science and Technology (KAIST), Daejeon 305-701, Republic of Korea
²Korea Aerospace Research Institute (KARI), Daejeon 305-701, Republic of Korea

Correspondence should be addressed to Hyochoong Bang; hcbang@ascl.kaist.ac.kr

Received 27 February 2017; Revised 26 April 2017; Accepted 27 April 2017; Published 31 May 2017

Academic Editor: Christian Circi

Copyright © 2017 Tania Savitri et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

This paper focuses on maximizing the percent coverage and minimizing the revisit time for a small satellite constellation with limited coverage. A target area represented by a polygon defined by grid points is chosen instead of using a target point only. The constellation consists of nonsymmetric and circular Low Earth Orbit (LEO) satellites. A global optimization method, Genetic Algorithm (GA), is chosen due to its ability to locate a global optimum solution for nonlinear multiobjective problems. From six orbital elements, five elements (semimajor axis, inclination, argument of perigee, longitude of ascending node, and mean anomaly) are varied as optimization design variables. A multiobjective optimization study is conducted in this study with percent coverage and revisit time as the two main parameters to analyze the performance of the constellation. Some efforts are made to improve the objective function and to minimize the computational load. A semianalytical approach is implemented to speed up the guessing of initial orbital elements. To determine the best parametric operator combinations, the fitness value and the computational time from each study cases are compared.

1. Introduction

Compared to a single satellite, a constellation may provide better coverage and higher reliability under failure of some satellites, ensuring a higher rate of survival and mission success. A constellation can offer unique capabilities that are often difficult to achieve through different means, for instance, enhanced temporal coverage [1]. A satellite orbit constellation design usually relies on the Walker approach [2] or streets of coverage method [3, 4]. Besides, a ground track-based approach has been developed [5], and recently the sliding ground track concept applied to constellations composed of one or more orbital planes has been introduced [6].

Related to the coverage mission, the idea of satellite constellation designs for complex coverage was presented by Ulybyshev [7]. Indrikis and Cleave, in their work about SPACE EGGS or the Satellite Coverage Model for LEO (Low Earth Orbit) Constellations [8], also tried to access the effectiveness of global, regional, and area coverage for proliferated small satellite constellations in low altitude orbits, stressing the capability of conventional analytical techniques. The idea of staged deployments of satellite constellations in LEO was also proposed by including the uncertainty feature of the expected number of users and activity level [9]. In the process, a satellite’s life cycle cost and its capacity is traded off to satisfy the minimum per-channel performance requirement. Another approach to optimize a satellite constellation design by using tiers of satellites with variations in the orbit’s altitude and inclination parameters was proposed Razoumny et al. in [10]. This approach tries to reduce the redundancy in the Walker approach by dividing the regional coverage into several latitude regions that can be addressed by different pairs of satellites’ altitude and inclination, allowing for a better revisit time or a reduced number of satellites.

Besides the analytic methods, the genetic algorithm (GA) is known for its robustness in obtaining a global optimum solution for nonlinear multivariable problems through its stochastic and heuristic search algorithms. In their thesis, Pegher and Parish [11] tried to compare the coverage optimization and the revisit time of sparse military satellite constellations using traditional approaches and GA. Another method of hybrid satellite constellation design called Genetic Satellite Constellation (GSC) was proposed in [12] by using single GA optimization. This approach demonstrates the performance of hybrid LEO/MEO, LEO/GEO, and MEO/GEO constellations. A study on Simulated Annealing and GA approaches to satellite constellation design for coverage of a limited latitude region was conducted by Crossley and Williams [13], where both methods outperformed the conventional Walker approach at low Earth central angles.

Multiobjective genetic algorithm (MOGA) is a variation of GA, which is known to be a robust technique to solve multiobjective optimization, resulting in a Pareto optimal set solution [14]. A class of fast and elite MOGA was introduced by Deb et al. in 2002 [15]. This new MOGA called NSGA-II stands for Nondominated Sorting GA, which uses a fast nondominated sorting and a crowding distance assignment in its algorithm. Several multiobjective optimization studies have been performed. For example, Ely used two-branch tournament GA for a constellation design with eccentric orbits, aiming to minimize the maximum constellation’s orbit altitude and the total number of satellites [16]. Several works from Ferringer are also focused on multiobjective optimization of satellite constellation designs including the trade-off analysis [17–20]. Related to MOGA, Mason et al. introduced a type of MOGA approach, called the Modified Illinois Nondominated Sorting Genetic Algorithm (MINSGA), combined with Satellite Tool Kit (STK) software to produce several constellation designs that provides a continuous global coverage [21]. Another work was also done by Confessore et al. to optimize a satellite constellation’s orbit design with GA aiming to minimize the number of satellites and maximize region coverage [22].

However, to the best of our knowledge, there is no specific unified design approach for a local continuous coverage or surveillance mission over a region. This paper addresses coverage and revisit time problems using a sparse satellite constellation design with limited coverage capability, as is often encountered in low-cost satellite missions. Some low-cost small satellites tend to use small cameras for image acquisition [23, 24]. Several recent developments made in this area are the BRITE mission, which is an explorer targeting bright stars, conducted by the Canadian Space Agency, partnering with the University of Vienna and Graz University of Technology from Austria and the Copernicus Astronomical Center from Poland [16], and the PRISM Earth-imaging validation mission developed by the University of Tokyo using a CMOS imager on a narrow-angle camera [25].

The contribution of this paper is twofold: the first is a method to address a complex coverage area by using sparse satellite constellation design with a limited coverage capability and the second is an efficient approach to lower the computational burden when performing the GA optimization. In this paper, a case study is conducted on a sparse satellite constellation using multiple circular LEO satellites equipped with imaging sensors. A GA-based optimization is chosen due to its ability to locate a global optimum solution for a wide class of nonlinear problems that may be applicable to the case of this study. Several constellations’ figures of merits [26], such as the area percent coverage and revisit time, are selected to evaluate the constellation performance using a multiobjective optimization algorithm. As much as GA seems to be a compatible optimization method to solve such a problem, it is often avoided due to its high computational load [27]. To minimize the computational load, a semianalytical approach is also proposed, which brings about a significant savings in computational time. The proposed approach, a combination of GA and an efficient semianalytical approach, will be a viable option for optimum satellite constellation designs.

2. Problem Statement

2.1. Mathematical Formulation of a Target Area

The area of interest is defined by inputting latitude-longitude coordinates of the main cities located in the country of interest. From these points, a convex hull polygon is generated using a computational geometry method, as depicted in Figure 1. The convex hull of a set of points is the smallest convex set that contains the points [28]. The area of interest is divided into several effective grid points, which are then evaluated to determine the target access time.

(a)

(b)

Figure 1: Example region (a) and the polygon made by convex hull mapping (b).

The grid points are defined in geodetic latitude and longitude , where and is the total number of grid points. Furthermore, zero ellipsoidal height of the target points is assumed.

2.2. Orbit Design

Initial orbit design consists of six standard orbital elements: semimajor axis , eccentricity , inclination , longitude of ascending node , argument of perigee , and mean anomaly .

The maximum Earth-centered half-angle and the maximum swath width of the satellite sensor occur when the satellite is at apogee . As shown in Figure 2, a satellite footprint projection on the Earth’s surface with elevation angle defines a coverage circle with a radius through the following equation [29]:

Figure 2: Satellite footprint projection.

3. Semianalytical Initial Guess

Using a semianalytical initial guess, the time step for propagation could be adjusted in order to reduce the overall computational load. Comparison between this method and the numerical propagation method in terms of the area coverage difference is presented as well in this paper. This method is based upon the correlation between nonsingular orbital elements and the target and . This approach is especially useful for observation over a repeating ground track orbit.

In comparison with numerical propagation, this approach allows for information about the satellite position at a specific time without having to propagate the satellite’s orbital elements or the satellite’s position and velocity vector over a time step. The accuracy of numerical propagation results depends on the size of the time step, which, on the other hand, conflicts with the size of available computing memory and computational load.

By the mean element theory [30], a satellite orbit is described using its mean orbital elements. The time rate of secular changes of mean , , and is given by [31]

where and the mean motion is given by

Drag perturbation is described using the following equations [32]:

3.1. Latitude Access

The correlation between the satellite mean orbital elements and the access to a target region’s latitude is defined in the following form [32]:

where represents the angular distance along the ground track measured from the ascending node at time = 0. The equation is only applied to an orbit with small eccentricity. On the assumption that the satellite in this study is equipped with a conical sensor, (6) becomes

For a noncircular orbit, the maximum swath width achieved on the perigee point is given by

It is recommended to use the time derivative of for a noncircular orbit with drag and perturbation included.

The latitude access time repeats for every period of mean argument of latitude or every . is defined as a set of latitude access times such that

3.2. Longitude Access

The correlation between satellite orbital elements and the area longitude is defined as

where is the mean argument of latitude and represents the entrance and exit of the access. A similar approach using a straightforward expansion to solve this equation is explained in [29] with a difference in the usage of instead of the usual , since the mean argument of latitude lies between and . The analytical solution approach is described briefly by first defining several variables. Let , , and be defined as follows:

Substituting (13)–(15) into (12) results in the following equation:

Equation (16) can be solved by using a solver in MATLAB®. is a quadrant-sensitive inverse of , which returns a value in the range . Using Taylor series expansion, can be defined as

Let us put a value of into the equivalent Kepler equation in terms of the mean element and derive :

is defined as a set of latitude access time for which

The time for intersection of a target point can be found from the time window intersection of each latitude and longitude access. The intersection time consists of access times, where

The flowcharts of the algorithm are depicted in Figures 3 and 4.

Figure 3: Latitude and longitude target access determination flowchart.

Figure 4: Target access time intersection flowchart.

4. Multiobjective Optimization

4.1. Pareto Approach

Multiobjective optimizations are usually expressed by introducing different objectives in a cost function, for which each objective importance is represented by a weight parameter for the -th objective, as in [33].

where

Using such a weighted sum approach, the objective functions should always be normalized or scaled so that their objective values are in similar magnitudes. However, according to [34], there is a certain drawback of minimizing the weighted sums of objectives for Pareto set generation in multicriteria optimization problems. The reason is that an evenly distributed set of weighting vectors cannot guarantee an evenly distributed representation of the Pareto optimal solutions.

4.2. Nondominated Sorting Genetic Algorithms Optimization (NSGA-II)

As described in Figure 5, NSGA-II, which stands for Nondominated Sorting GA, uses a fast nondominated sorting and a crowding distance assignment [8]. GA optimization itself is a heuristic, stochastic global search based on the Darwinian concepts of evolution and natural selection and was first introduced by John Holland in mid-1960s [33]. An elitist GA favors individuals with a better fitness measure (rank). A controlled elitist GA, which is used in this algorithm, facilitates the optimization process to converge into a set of diverse solutions, even if the solution might consist of several lower performing individuals. NSGA-II is selected for this research because it has been shown to exhibit a better spread and convergence, relatively close to the true Pareto optimal front compared to PAES and SPEA [34]. SPEA stands for Strength Pareto Evolutionary Algorithm, and PAES is Pareto Archived Evolution Strategy algorithm. PAES is a multiobjective optimizer approach, which uses an archive of previously found solutions in order to identify the dominance ranking of the current and candidate solution vectors. SPEA provides a form of elitism by using an archive of the nondominated set, which is maintained separately from the population of candidate solutions.

Figure 5: NSGA-II schema.

In the first generation, parent population of individuals creates an offspring population also in individuals through the common tournament Selection, Crossover, and Mutation operators. Next, population with a size of individuals goes through a nondominated sorting. After the first generation, crowded tournament selection is used to select the top individuals from to form the next generation . A solution that wins a tournament has either the highest rank or a better crowding distance parameter. The crowding distance parameter itself is defined as the largest cuboid surrounding a solution, for which no other solution is present [35].

Two parameters used in this algorithm to control the elitism are the Pareto fraction, which is the nondominated sorting using the Pareto set, and the distance function, which uses a crowding distance assignment. The Pareto fraction limits the number of individuals that can be put into the Pareto front (elite members). The distance function itself favors individuals who are located relatively far from the front since it can create a diverse, noncrowded set of solutions.

In this paper, the parametric operator is varied for the number of generations and population. The number of generations is usually set at the beginning of the optimization and is usually used as a stopping criterion for the optimization. Another useful parameter is the stall generation, where the optimization stops if the best fitness function value does not stop after hitting several generations until the stall generation is met. Population, which has a constant value during the optimization, is the number of individual sets for every generation.

5. Simulations

A target region for the case South Korea, with its 37 cities’ latitude and longitude positions, is used as an input for simulation. The area of interest itself is represented as a polygon of boundary points that includes all of the 37 cities. This area of the polygon is then divided into a number of grid points. denotes the number of fast nondominated target points which is the number of grid points located inside the polygon, for which grids dividing the example area result in a total of 50 grid points. The basic design criteria and constraints are as follows. A maximum satellite constellation of three satellites with a 10-degree half-angle of a conical sensor is used. This conical sensor is set as well so that only one satellite cannot cover the whole area of interest. Although from (1) it was mentioned that the maximum swath width occurs at the apogee, due to the circular orbit used in this simulation, this parameter is a matter of the satellite altitude only, assuming that the pointing error is neglected.

The constellation is placed into a circular orbit at an epoch date on August 1, 2020. Up to six hours of the satellite’s lifetime is studied. In comparison with the simulation done in [29], the satellite’s orbiting time period is much smaller; therefore, the argument of perigee and mean anomaly are included as optimization variables to compensate for the shorter observation time.

The 6-hour slot of the satellite lifetime was studied as an example of applying the proposed optimization algorithm to satellite constellation orbit design. The part of the satellite lifetime was selected to compare performance of the proposed semianalytical algorithm to that of typical genetic algorithms in a perspective of computational time. To get optimized results for whole time period of repeating ground track seemed unnecessary since the satellites in the constellation are on the same configuration every about 30 days and therefore much more time for simulation is needed.

Imaging sensors is assumed for this study. The illuminating conditions will change at every pass because the orbit is not Sun-synchronous. However, we believe it is feasible to use imaging sensors in such constellation as existing systems equipped with CMOS or thermal infrared camera.

The constraints for initial orbit elements used as optimization variables are defined in Table 1.

Table 1: Parameters variation for constellation visibility GA simulation.

Variation parameter
Value

Inclination
0–90 degrees

Semimajor axis
6538 to 6978 km, orbit height 160 to 600 km (LEO)

Longitude of ascending node
0–360 degrees

Argument of latitude
0–360 degrees

A comparison between the 4th order Runge-Kutta (4RK) numerical integration and the semianalytical approach is performed as well. The results match within an error of 0.15% to 4.33% over a 1.5- to approximately a 13-day orbit propagation time. In a single application, the semianalytical approach can save up to four times computational time compared to the numerical integration approach, using a computer with an Intel ® Core™2 Duo CPU E7500 2.93 GHz processor and only 2.0 GB RAM. This test was performed using a random target site at 37° latitude and 14°E longitude. A satellite located in orbit with 0.005 eccentricity,  km, , and 30° is propagated over five orbital periods. Figure 6 shows the error growth between the 4RK coverage results and the semianalytical approach, which can reach 4 minutes and 42 seconds during the longest 13-day propagation time. From the data, it is predicted that the difference will grow as the propagation time increases.

Figure 6: Coverage time comparison between 4th order Runge-Kutta propagation versus semianalytical approach.

In overall application, the optimization’s computational load can be reduced up to nine times using the semianalytical algorithm. If the satellite does not cover any target latitude, the algorithm automatically assigns a high value to the cost function and skips the access computation for longitude sections. Such bad performing individuals are automatically marked by the NSGA-II in the overall computation to exclude this individual in the next generation evaluation. This approach can speed up the entire computational time by skipping at least half of the calculation.

The fitness function of each individual is evaluated for every generation, and while the evaluation does not hit the stopping criteria, the selected parents continue to mate and reproduce new generations via Crossover and Mutation. The stopping criteria of the NSGA-II optimization algorithm used in this paper use MATLAB® default numbers, except for the TolFun or tolerance function being set at 10⁻³ (driven by the total number of grids and the definition of the fitness function formula). Optimum constellation design then covers as many grid points as possible and has a minimum average revisit time (ART), also with the lowest maximum revisit time.

Five out of six orbital elements, semimajor axis , inclination , longitude of ascending node , argument of perigee , and mean anomaly , are used as optimization design variables, while eccentricity is kept as a constant. As a note, and are included as design variables, despite the circular orbit assumption, to maximize the range of individual diversity. Using (4) and (5), the drag perturbation effects can be modeled and show increasing over propagation time because drag tends to lower the semimajor axis of a satellite. The data for drag is obtained from [35] in Table 2.

Table 2: Drag parameters for a sample satellite (330 km LEO alt., e = 0, i = 79°).

Satellite parameter
Value

Mass (kg)
2.0

Area (m²)
0.01

2.0

The effect of drag is presented in Figure 7.

Figure 7: Effect of drag perturbation on the angle.

There are four optimization objectives in this study, which is divided into two main case studies as follows: the coverage factor and the revisit time factor. The percent coverage for any point on the grid is described as the number of points covered, which is divided by the total number of points in the grid. The percent of coverage directly shows the effectiveness of satellite coverage. However, it cannot provide any information about the distribution of the gaps [32]. The coverage gap is defined here as the length of time where a point is not covered by any of the satellites in the constellation.

In the following, those two figures of merit of satellite constellations are represented in four objectives. Regarding the previous discussion about the convexity of the Pareto solution, the cost function in this study is represented in several objectives instead of using a weighted cost function. The area coverage is different from the time coverage since the time coverage is a summation of total satellite access time divided by the number of grids and the total simulation time. The following algorithm for GA-based optimization is depicted in Figure 8:

(1)Maximum area percent coverage

(2)Maximum average area time coverage

(3)Minimization of maximum coverage gap time

(4)Average coverage gap time

Figure 8: Schema for GA-based optimization algorithm.

The reward is defined as a negative value that is added to the cost function if 100% area coverage can be achieved by less than the maximum initial three satellites. Several case studies that were conducted are presented in Figure 9.

Figure 9: Parameters variations for percent coverage and ART by NSGA-II simulation.

Using the semianalytical technique, the time of access accuracy is independent of the time step since it does not use numerical integration such as the Runge-Kutta orbit propagation. In this simulation, test cases are analyzed by varying two types of parameters: the optimization algorithm parameters (or the NSGA-II parameters in this case) and the constellation orbit design parameters. The NSGA-II parameters that can be varied include the population, generation, and the stall criteria. Other adjustments that could be made for Selection, Crossover, and Mutation parameters for the population mating and generation are kept at their default values this time since the effect of changing those parameters has been discussed previously [33, 36].

In Case , a different number of populations and generations are tried. By increasing the number of populations and generations, there is a greater probability to achieve a global optimum for any nonlinear optimization problem as the number of search points increases. The analogy for this can be seen in Figure 10, where the increment in the number of populations looks like a mesh figure with an increasing number of mesh points. This results in a bigger search space throughout the overall optimization framework. In our case, it indicates that the Pareto front generated will have a better value, as can be seen in Case . However, Case suffers from longest computational time due to the increased number of cost function counting.

Figure 10: Mesh figure for nonlinear problem.

Case iterates the constellation orbit design parameters by increasing the number of satellites in the constellation from 3 to 6. The increment in the number of satellites can increase the area percent coverage and coverage time as well while also reducing both ART and maximum revisit time. However, the increment achieved is still lower than the results achieved in Case , even though the computational load is about two times lower. The role of the inclination angle in this satellite constellation design is to maximize the launch vehicle capability by choosing the inclination that is closest to the target latitude.

Cases , , and all employ and only as optimized variables. In comparison with Case , the number of populations and generations is increased in Case , and only two orbital planes are used in Case . It can be seen from Figures 14–16 that the largest percent coverage from all four cases (a comparison with Case included) is achieved by Case on 88% coverage. Both the lowest ART and maximum revisit time are achieved by Case at 1.614 and 0.896 hours, respectively. Case achieves the best overall performance in all four cases, although the CPU computational load significantly increases up to 4.5 times the average CPU computational load for all the three cases. To strengthen the results, Case is compared with Case , both with 750 population and 100 generations. From these two cases, the results are almost identical in terms of the area percent coverage with only 2% difference. In terms of the time percent coverage, Case exhibits slightly better coverage with 5 minutes difference, while in terms of ART, Case achieves 16.8 minutes shorter time compared to Case . Figures 11–16 show the results for various cases.

Figure 11: NSGA-II result for and objectives (Cases , , , and ).

Figure 12: NSGA-II result for and objectives (Cases , , , and ).

Figure 13: NSGA-II result for and objectives (Cases , , , and ).

Figure 14: NSGA-II result for and objectives (Cases , , , and ).

Figure 15: NSGA-II result for and objectives (Cases , , , and ).

Figure 16: NSGA-II result for and objectives (Cases , , , and ).

Case is conducted to seek an optimum design for a constellation of the same three satellites in two orbital planes in order to reduce the launch cost by using a combination of , , and . Seven variables are used in total for three satellites and two orbits, two slots for , two for Ω, and three for ω variation. However, the separation between two satellites in the same plane is considered too small with only 13° difference in argument of latitude. Therefore, an additional case is attempted using five optimized variables in total, which are , , , , and in Case . In Case , a constellation of three satellites distributed in two orbital planes is used. In comparison with Case , although the area coverage performance is similar, the time coverage achieved in Case is about 15 minutes longer. Also, it can be seen that Case has a greater number of solutions in the Pareto front with near optimum values. From the previous trends, regarding optimization results with respect to the number of optimized variables, one can see that the optimization using NSGA-II leads to a better set of solutions using a lower number of optimized variables (only and , compared to variations of , , , , and ) within a considerably lower number of populations and generations. Figures 17–19 show the results for Cases and .

Figure 17: NSGA-II result for and objectives (Cases and ).

Figure 18: NSGA-II result for and objectives (Cases and 12).

Figure 19: NSGA-II result for and objectives (Cases and ).

Cases and were conducted to see the results if the number of populations and generations is increased separately. In comparison with Cases and , Case yields the best performance in all of the four optimized objectives, while the second-best performance is by Case . The results of the simulation are given in Table 3. From Table 4, similar or better performance with Case is achieved by Cases , , and with approximately 15% less computational load. By increasing the number of populations, the NSGA-II optimization approach allows a greater chance of finding the global optimum prior to converging at a cost of increased computational load. This argument explains why the results from Case are similar to Case despite only half of the function evaluations in Case . Figures 20–22 show the results for Cases , , , and .

Table 3: Results for the NSGA-II simulation.

Case number
points in Pareto front
Computation time (hr)

1
105
4.104

2
199
16.134

3
125
8.308

4
175
4.279

5
16
4.116

6
263
18.864

7
175
4.098

8
175
3.716

9
175
12.584

10
263
5.386

11
175
3.913

12
175
3.942

Table 4: Best result for all objective functions (, , , ).

Case number
(%)
(min)
(hour)
(hour)

1
42
12.40
3.070
2.938

2
86
48.25
1.727
0.899

3
42
10.01
1.951
1.346

4
44
13.05
1.683
1.148

5
64
26.25
1.614
0.896

6
88
53.37
1.715
1.179

7
76
47.54
3.398
1.193

8
58
22.39
3.039
1.359

9
62
21.48
1.682
1.179

10
54
19.29
2.945
1.235

11
56
22.69
3.089
1.240

12
74
32.99
2.493
1.078

Figure 20: NSGA-II result for and objectives (Cases , , , and ).

Figure 21: NSGA-II result for and objectives (Cases , , , and ).

Figure 22: NSGA-II result for and objectives (Cases , , , and ).

Figures 23–25 present a comparison between Cases , , , and results with a different number of orbital elements used for optimization variables. The comparison indicates that Case yields the best results in terms of overall objectives (percent coverage and revisit time) followed by Case with the best area coverage at 56% and then Case at 44%. However, in terms of both the maximum revisit time and ART, the second rank is held by Case , while Case offers a longer maximum revisit time at 3.089 hours compared to 1.683 hours in Case . The least performing constellation is Case with a variation of , , , and as optimized variables, while the best performing constellation is Case with only and as optimized variables. Figures 26–29 show satellites’ full ground tracks for various conditions and Tables 5–8 summarize orbit elements for the solutions.

Table 5: Orbit elements for a solution in Case (72% area coverage and 1 hour and 19 minutes of ART).

Orbit element
Sat
Sat
Sat

(km)
6978
6978
6978

0
0
0

(deg)
42
42
44

(deg)
90
90
131

(deg)
324
326
284

(deg)
0
0
0

Table 6: Orbit elements for a solution in Case with best maximum revisit time at 3 hours and 24 minutes and 32% area coverage.

Orbit element
Sat
Sat
Sat

(km)
6978
6978
6978

0
0
0

(deg)
42
42
53

(deg)
154
154
136

(deg)
261
274
76

(deg)
0
0
0

Table 7: Orbit elements for a solution in Case (best 74% area coverage and 9 hours of ART).

Orbit element
Sat
Sat
Sat

(km)
6936
6936
6884

0
0
0

(deg)
38
38
75

(deg)
199
199
9

(deg)
251
159
158

(deg)
335
163
202

Table 8: Orbit elements for a solution in Case with best maximum revisit time at 2 hours and 30 minutes and 44% area coverage.

Orbit element
Sat
Sat
Sat

(km)
6938
6938
6871

0
0
0

(deg)
39
39
78

(deg)
198
198
118

(deg)
251
245
291

(deg)
333
342
23

Figure 23: NSGA-II result for and objectives (Cases , , , and ).

Figure 24: NSGA-II result for and objectives (Cases , , , and ).

Figure 25: NSGA-II result for and objectives (Cases , , , and ).

(a)

(b)

Figure 26: Satellites’ full ground track (a) and close-up (b) for a solution in Case with 72% area coverage and 1 hour and 19 minutes of ART. Red line: ground track of sat , pink line: ground track of sat , and blue line: ground track of sat . Green triangle: area covered by sat , pink triangle: area covered by sat , and blue star: area covered by sat .

(a)

(b)

Figure 27: Satellites’ full ground track (a) and close-up (b) for a solution in Case with best maximum revisit time at 3 hours and 24 minutes and 32% area coverage. Red line: ground track of sat #1, pink line: ground track of sat #2, and blue line: ground track of sat #3. Red star: area covered by sat #1 and blue triangle: area covered by sat #3. In this case, area covered by satellites #1 and #2 is similar, so the paths are stacked on top of each other.

(a)

(b)

Figure 28: Satellites’ full ground track (a) and close-up (b) for a solution in Case (best 74% area coverage and 9 hours of ART). Red line: ground track of sat , pink line: ground track of sat , and blue line: ground track of sat . Pink triangle: area covered by sat and red star: area covered by sat .

(a)

(b)

Figure 29: Satellites’ full ground track (a) and close-up (b) for a solution in Case with best maximum revisit time at 2 hours and 30 minutes and 44% area coverage. Purple line: ground track of sat and blue line: ground track of sat . Red star: area covered by sat .

In comparison with the best area coverage result in Case , it can be seen that, by introducing as an optimized variable in Case , the satellites’ ascending node in the same orbit plane (satellite 1 and 2) is separated up to 264°. Nonetheless, this separation does not create significant changes in terms of percent coverage nor the revisit time performance compared to the results of Case .

In summary, the overall objective function can be improved using several methods that are described in Table 9. This table was made by taking the average of performances from the cases that have either the same orbital elements as optimized variables or cases represented by the same number of populations and/or generations.

Table 9: Summary of parameter variation optimization results.

Variation
Percent coverage
Average revisit time
Computational load
Pareto front

Increasing number of populations
Better up to 29%
Better up to 54%
Worse up to 31%
Increasing number of solutions in Pareto front

Increasing number of generations
Better up to 48%
Better up to 60%
Worse up to 206%
N/A

Increasing number of populations & generations
Better up to 71%
Better up to 69%
Worse up to 293%
Increasing number of solutions in Pareto front

Using a variation of , , , , and
Better up to 33%
Better up to 82%
N/A
N/A

Using a variation of , , and
Better up to 5%
Better up to 61%
N/A
N/A

Using a variation of and
Better up to 29%
Better up to 70%
N/A
N/A

Fuel penalty used
Better up to 38%
Better up to 54%
N/A
N/A

Semisparse constellation (Cases & )
Better up to 25%
Better up to 61%
N/A
N/A

Increasing number of satellites in constellation
Same value
Better up to 54%
Worse up to 102%
N/A

The maximum coverage time is around 53 minutes on average per satellite orbit, which was achieved in Case . Also, 88% maximum area coverage was achieved in Case . In terms of the revisit time, the minimum ART is 54 minutes in Case , while the lowest maximum gap time is around 1 hour and 37 minutes for Case .

6. Conclusion and Future Work

By implementing a semianalytical approach, this study showed that computational load was reduced by over nine times within 0.5% error without having to numerically integrate satellite position or orbital elements over time to get the area coverage. Several figures of merits were chosen to analyze the performance of the optimized constellation orbit design, including maximum coverage, average coverage time, ART, and maximum revisit time over a target area. Circular LEO satellite constellation was used as a case study with a combination of , , , and ω as optimization design variables. The satellites in this constellation design case were assumed to have a small conical sensor angle, which means that the satellites were subjected to a limited coverage. Therefore, the study cases were defined as well to achieve maximum coverage efficiency from the available number of satellites. Several cases were analyzed to see how the NSGA-II parameter or the orbital parameter changes might impact the optimizations results. To determine the best parametric operator combination, NSGA-II employed the Pareto concept to give a set of points that are able to satisfy all objectives.

From the conducted study cases, it was shown that results could be improved by increasing the number of populations and generations in NSGA-II or by increasing the number of satellites in the constellation. The latter option was considered inefficient since the results of those two cases were not much different, and the option of adding more satellite in the constellation may cost more. Also, by minimizing the number of optimization variables, better results were obtained since the complexity of the problem was decreased with the decreased number of variables. The percent of coverage was also improved by introducing a fuel penalty. Future work, including analysis of sensor pinpointing accuracy, may yield more meaningful results and make this problem closer to a real case.

Nomenclature

: | Reference cross-sectional area of the satellite, m²
: | Orbit semimajor axis, km
: | Drag coefficient of the satellite
: | Eccentricity
: | -th cost function
: | Orbit altitude, km
: | Maximum number of satellites in the constellation
: | Orbit inclination, deg
: | Earth oblateness gravity harmonic coefficient,
: | Mean argument of latitude, deg
: | Mean anomaly, deg
: | Mass of the satellite, kg
: | Total number of target points
: | Orbit mean motion, deg/s
: | Number of grids
: | Number of individuals
naked: | Number of uncovered grids by satellite denoted with in percent ratio of all total available grids
: | Spherical radius of the Earth, 6378.1363 km
: | Satellite’s apogee height, km
: | Time, s
: | Sensor half-angle of conical center, deg
: | Elevation angle, deg
: | Target location longitude, deg, where
: | Gravitational parameter of the Earth,
: | Longitude of ascending node, deg
: | Atmospheric density, kg/m³
: | Target location latitude, deg, where
: | Earth-centered half-angle, deg
: | Argument of perigee, deg
: | The Earth’s rotation rate,  rad/s
: | True anomaly, deg.

Conflicts of Interest

The authors declare that there are no conflicts of interest.

References

A. da Silva Curiel, A. Cawthorne, and M. Sweeting, “Progress in small satellite technology for Earth observation missions,” in Small Satellites for Earth Observation: Selected Proceedings of the 5th International Symposium of the International Academy of Astronautics, Berlin, April 4–8 2005, p. 50, 2005.

J. G. Walker, “Satellite constellations,” Journal of the British Interplanetary Society, vol. 37, pp. 559–572, 1984.

R. D. Luders, “Satellite networks for continuous zonal coverage,” ARS Journal, vol. 31, no. 2, pp. 179–184, 1961.

T. J. Lang and W. S. Adams, “A comparison of satellite constellations for continuous global coverage,” in Mission Design & Implementation of Satellite Constellations, pp. 51–62, Springer, 1998.

J. C. King, Quantization and Symmetry in Periodic Coverage Patterns with Applications to Earth Observation, 1975.

C. Circi, E. Ortore, and F. Bunkheila, “Satellite constellations in sliding ground track orbits,” Aerospace Science and Technology, vol. 39, pp. 395–402, 2014.

Y. Ulybyshev, “Satellite constellation design for complex coverage,” Journal of Spacecraft and Rockets, vol. 45, no. 4, pp. 843–849, 2008.

J. Indrikis and R. Cleave, 'SPACE EGGS'-Satellite Coverage Model for Low Earth Orbit Constellations, 1991.

O. L. De Weck, R. D. Neufville, and M. Chaize, “Staged deployment of communications satellite constellations in low earth orbit,” Journal of Aerospace Computing, Information and Communication, vol. 1, pp. 119–136, 2004.

Y. N. Razoumny, P. G. Kozlov, V. Y. Razoumny, and A. A. Moshnin, “On optimization of Earth coverage characteristics for compound satellite constellations based on orbits with synchronized nodal regression,” in Proceedings of the 2nd IAA Conference on Dynamics and Control of Space Systems, Roma, Italy, March 2014.

D. J. Pegher and J. A. Parish, Optimizing Coverage and Revisit Time in Sparse Military Satellite Constellations: A Comparison of Traditional Approaches and Genetic Algorithms, DTIC Document, 2004.

M. Asvial, R. Tafazolli, and B. G. Evans, “Genetic hybrid satellite constellation design,” Parameters, vol. 1, p. 100, 2003.

W. A. Crossley and E. A. Williams, “Simulated annealing and genetic algorithm approaches for discontinuous coverage satellite constellation design,” Engineering Optimization, vol. 32, no. 3, pp. 353–371, 2000.

C. M. Fonseca and Fleming P. J., “Genetic algorithms for multiobjective optimization: formulation discussion and generalization,” in Proceedings of the 5th International Conference on Genetic Algorithms (ICGA '93), pp. 416–423, 1993.

K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan, “A fast and elitist multiobjective genetic algorithm: NSGA-II,” IEEE Transactions on Evolutionary Computation, vol. 6, no. 2, pp. 182–197, 2002.

T. Ely, “Satellite constellation design for zonal coverage using genetic algorithms,” in Proceedings of the 8th AAS/AIAA Space Flight Mechanics Meeting, pp. 443–460.

M. P. Ferringer and D. B. Spencer, “Satellite constellation design tradeoffs using multiple-objective evolutionary computation,” Journal of Spacecraft and Rockets, vol. 43, no. 6, pp. 1404–1411, 2006.

M. P. Ferringer, R. S. Clifton, and T. G. Thompson, “Efficient and accurate evolutionary multi-objective optimization paradigms for satellite constellation design,” Journal of Spacecraft and Rockets, vol. 44, no. 3, pp. 682–691, 2007.

M. P. Ferringer, R. S. Clifton, and T. G. Thompson, “Constellation design with parallel multi-objective evolutionary computation,” in Proceedings of AIAA/AAS Astrodynamics Specialist Conference and Exhibit, Keystone, Colo, USA, 2006.

M. P. Ferringer, D. B. Spencer, P. M. Reed, R. S. Clifton, and T. G. Thompson, “Pareto-hypervolumes for the reconfiguration of satellite constellations,” in Proceedings of the AIAA/AAS Astrodynamics Specialist Conference and Exhibit, pp. 1–31, August 2008.

W. J. Mason, V. Coverstone-Carroll, and J. Hartmann, Optimal Earth Orbiting Satellite Constellations via A Pareto Genetic Algorithm, University of Illinois at Urbana-Champaign, 2001.

G. Confessore, M. Di Gennaro, and S. Ricciardelli, “A Genetic Algorithm to Design Satellite Constellations for Regional Coverage,” Operations Research Proceedings, vol. 2001, pp. 35–41.

N. Orr, J. Eyer, B. Larouche, and R. Zee, “Precision formation flight: the CanX-4 and CanX-5 dual nanosatellite mission,” in Proceedings of the 21st Annual AIAA/USU Conference on Small Satellites, 2007.

H. Ashida, K. Fujihashi, S. Inagawa et al., “Design of Tokyo Tech nano-satellite Cute-1.7+APD II and its operation,” Acta Astronautica, vol. 66, no. 9-10, pp. 1412–1424, 2010.

M. Komatsu and S. Nakasuka, “University of Tokyo Nano Satellite Project 'PRISM',” Transactions of Space Technology Japan, vol. 7, 2009.

J. R. Wertz, Mission Geometry; Orbit and Constellation Design and Management: Spacecraft Orbit and Attitude Systems, vol. 13 of Space Technology Library, Microcosm, Boston, Mass, USA; Kluwer Academic Publishers, El Segundo, Calif, USA, 2001.

E. S. H. Hou, N. Ansari, and H. Ren, “Genetic algorithm for multiprocessor scheduling,” IEEE Transactions on Parallel and Distributed Systems, vol. 5, no. 2, pp. 113–120, 1994.

D. Sunday, The Convex Hull of a 2D Point Set or Polygon, vol. 109, 2006.

P. Sengupta, S. R. Vadali, and K. T. Alfriend, “Satellite orbit design and maintenance for terrestrial coverage,” Journal of Spacecraft and Rockets, vol. 47, no. 1, pp. 177–187, 2010.

Y. Kozai, “The motion of a close earth satellite,” The Astronomical Journal, vol. 64, pp. 367–377, 1959.

H. Schaub and J. L. Junkins, Analytical Mechanics of Space Systems, AIAA, 2003.

J. R. Wertz, Mission Geometry; Orbit and Constellation Design and Management, Space Technology Library, 2001.

H.-D. Kim, O.-C. Jung, and H. Bang, “A computational approach to reduce the revisit time using a genetic algorithm,” in Proceedings of the International Conference on Control, Automation and Systems (ICCAS '07), pp. 184–189, IEEE, Seoul, Republic of Korea, October 2007.

I. Das and J. E. Dennis, “A closer look at drawbacks of minimizing weighted sums of objectives for Pareto set generation in multicriteria optimization problems,” Structural Optimization, vol. 14, no. 1, pp. 63–69, 1997.

D. Ilyas, Orbital Propagation and Formation Flying of CubeSats within QB50 Constellation, Science and Technology (SpaceMaster), 2011.

T. Savitri, O. Jung, and H. Bang, “Optimization of grid points coverage by satellite constellation using genetic algorithm,” in Proceedings of the Asia-Pacific International Symposium, Takamatsu, Japan, 2013.