Abstract
The combination of atmospheric drag and lunar and solar perturbations in addition to Earth’s oblateness influences the orbital lifetime of an upper stage in geostationary transfer orbit (GTO). These high eccentric orbits undergo fluctuations in both perturbations and velocity and are very sensitive to the initial conditions. The main objective of this paper is to predict the reentry time of the upper stage of the Indian geosynchronous satellite launch vehicle, GSLV-D5, which inserted the satellite GSAT-14 into a GTO on January 05, 2014, with mean perigee and apogee altitudes of 170 km and 35975 km. Four intervals of near linear variation of the mean apogee altitude observed were used in predicting the orbital lifetime. For these four intervals, optimal values of the initial osculating eccentricity and ballistic coefficient for matching the mean apogee altitudes were estimated with the response surface methodology using a genetic algorithm. It was found that the orbital lifetime from these four time spans was between 144 and 148 days.
1. Introduction
The estimation of the orbital lifetime of decaying upper stages is important, due to the impact risks associated. In fact, the geostationary transfer orbit (GTO) is a high eccentric orbit which traverses low Earth orbit (LEO) and geostationary orbit (GSO), where the threat of collision with an operational satellite is high. Hence, minimizing the lifetime of an upper stage is a significant consideration in space debris mitigation endeavors. In order to reduce the collision probability, GSO satellites are placed in graveyard orbits at the end of their operational life.
Lunar and solar perturbations on geostationary transfer orbits account for the oscillation of the eccentricity keeping the semimajor axis constant [1]. The effect of lunar-solar perturbations in reducing the orbital lifetime of high eccentricity orbits is also the topic of interest of earlier studies [1–3]. An accurate simulation of the process requires a good estimate of the initial state [4]. Since the semimajor axis is directly related to the period of an orbit, it can be easily measured. Hence, the eccentricity and the ballistic coefficient which depends on the drag coefficient , the mass of the object , and the effective area are considered as uncertain parameters [5].
In this paper, a method to perform the lifetime estimation of an upper stage in GTO is presented as an optimal estimation problem [5, 6]. Using this method we predict the reentry time of the cryogenic stage of the Indian geosynchronous satellite launch vehicle GSLV-D5, employing the orbital data of the rocket body (R/B) in the form of two-line element sets (TLEs), which were downloaded from the Space Track Organization website (http://www.spacetrack.org/) from February 11, 2014, to April 6, 2014. Based on near linear variation of the mean apogee altitude, four intervals were identified for the reentry study [6]. The response surface method with genetic algorithm (GA) is effectively used to determine the optimal initial estimates of eccentricity and ballistic coefficient using the TLEs of each time interval [7]. The reentry time estimation in each case is computed using the Numerical Prediction of Orbital Events (NPOE) software by C. David Eagle (http://www.cdeagle.com/html/npoe.html). The influences of lunar-solar perturbations, Earth’s oblateness, and atmospheric drag are considered to predict the reentry time more accurately [1, 8, 9].
2. Method of Prediction
The orbital elements data given as two-line element sets (TLEs) were downloaded from the Space Track website from February 11, 2014, to April 4, 2014. Since the orbital period was greater than 225 minutes, the SDP4 orbit propagator model [10, 11] was used in the satellite tracking SatSpy software (http://www.satspy.com/) to get the state vectors, consisting of position and velocity, at the particular epochs chosen for the reentry time predictions of the GSLV R/B. SDP4 extends the SGP4 model to lunar-solar gravity perturbations, solar radiation pressure effects, and Earth resonance terms. The state vectors obtained are then given in input to the NPOE software to obtain the mean and osculating orbital elements at the initial state for the orbit simulation. The initial values of the uncertain parameters, ballistic coefficient and eccentricity, are estimated with the response surface technique using a genetic algorithm for the four time intervals where a near linear variation of the mean apogee altitude is observed. Using these initial values for each time interval, the reentry time of the GSLV R/B is predicted. The observed and predicted values of solar flux () and geomagnetic index () data utilized for density computations are obtained from the Solar Terrestrial Activity Report by Jan Alvestad (http://www.solen.info/solar/).
2.1. Response Surface Methodology Using Genetic Algorithm
The response surface methodology (RSM) is a combination of mathematical and statistical techniques used widely in various fields for the purpose of development of models and optimization [12]. The objective of using RSM is to optimize a response (output variable) which is influenced by a number of independent variables (input variables) [13, 14]. Consider the response influenced by two independent variables and , which can be represented mathematically as follows: where is the random error observed in the response . The surface characterized by is called the response surface.
The actual relationship between the response and the independent variables is usually unknown. Therefore, the first step in RSM is to find the approximate model function which shows the true relationship. Generally, the approximation is started with a low-order polynomial function. If the response is characterized by a linear function, then a first-order model is used as approximation function. If the response surface has a curvature, then a higher-order polynomial function is used, such as a second-order model. The objective of RSM is not only to understand the contour of the response surface, but also to estimate the values of independent variables called “design variables,” for which the optimal response is obtained.
In order to obtain the optimal response, a genetic algorithm (GA) which belongs to the group of evolutionary algorithms (EA) is used [15, 16]. In the genetic algorithm there are six steps to be taken.(1)Creation of an initial population: an initial set of solutions or chromosomes is generated randomly or by a seeding procedure.(2)Fitness evaluation: the quality of each solution of the population is evaluated using the fitness function.(3)Selection: with this process promising solutions are chosen to pass to the next generation at the expense of other solutions which are considered ill-equipped for the objective.(4)Crossover: this process consists of taking two strings or parents from the population and performing a random exchange of portions between them to form a new solution. This new chromosome has information from both parents. The crossover does not apply to the entire population string, but it is limited by the crossover rate.(5)Mutation: this involves making changes in individual values of variables in a solution. Mutations serve to maintain the diversity of the population, reducing the probability of finding a local minimum or local maximum rather than the global optimal solution.(6)Checking if the stopping criterion is satisfied: if the stopping criterion is not satisfied, the process returns to step number 3. If the criterion is satisfied, the algorithm finishes.
For the reentry time prediction, we consider the initial osculating eccentricity and the ballistic coefficient as the design variables to find the initial state. The mean apogee altitude generated for the considered initial osculating eccentricity and ballistic coefficient, referred to as mean apogee surface [4], is the response surface. The topography of the mean apogee surface reveals the dynamics of epochal motion. Four time intervals shown in Figure 2 are chosen based on near linear variation of the mean apogee altitude. Each time interval is referred to as a “zone.” Then, we select three initial values of the osculating eccentricity () and the ballistic coefficient () such that the observed mean apogee altitudes fall well within the lower and the upper bounds of the mean apogee surfaces. With the selected initial values of and , we generate nine mean apogee surfaces. Each surface is generated by propagating the trajectory with the selected initial values of and using NPOE. These mean apogee surfaces are used to compute the predicted mean apogee altitude at a specified epoch by interpolation [17]. The genetic algorithm is used to find the optimal solution of the initial values of osculating eccentricity and ballistic coefficient which match the mean apogee altitude. The population size is considered to be 20. Because all optimization parameters have a specified range, a binary coded GA is utilized and all parameters are coded in 60 bits. A single point crossover probability of 0.8 and mutation probability of 0.01 are selected. The cycle of the genetic algorithm is given in Figure 1.
3. Reentry Time Prediction for GSLV-D5
The TLEs downloaded from the Space Track Organization website are converted into state vectors using the SatSpy software and the mean and osculating orbital elements are computed using the NPOE software. The observed mean apogee and perigee altitudes computed from the TLEs are plotted in Figures 2 and 3, respectively. The identified four zones between February 11, 2014, and April 6, 2014, are labelled as A, B, C, and D, respectively. Then, the initial osculating eccentricity and ballistic coefficient are computed for each of the considered zones using the above methodology. The optimal values of and are used in the NPOE software to predict the reentry time, keeping the other orbital elements unchanged. For the orbit propagation using the NPOE software, the terms up to of the Earth gravity model based on GEM10B [18] are considered. For the atmospheric drag perturbations, the MSIS90 density model [19], which includes the observed and predicted values of solar flux and geomagnetic index, is used. The lunar and solar perturbation forces were also included. The mean perigee altitude has gone below 140 km in Figure 2 due to the solar perturbations which has dominated the forces acting on the object.
3.1. Case Studies
3.1.1. Zone A
The osculating orbital elements of GSLV R/B (NORAD No. 39499) as obtained from the TLE on February 13, 2014, 15 : 24 (UTC), are as follows: semimajor axis (km) = 24163.152291, eccentricity = 0.7279715192, inclination (°) = 19.33662306, argument of perigee (°) = 207.14804219, right ascension of the ascending node (°) = 196.09051738, true anomaly (°) = 152.63951389.
To generate a set of mean apogee surfaces for zone A, three values of initial osculating eccentricity (0.7278715192, 0.7279715192, and 0.7280715192) and three values of ballistic coefficient (60, 80, and 100) are selected to obtain nine grid points as plotted in Figure 4. The mean apogee altitudes for the osculating eccentricity and the ballistic coefficient are predicted using the NPOE software. Similarly, the mean apogee altitudes for (), (), (), (), (), (), (), and () are predicted using the NPOE software and plotted using MATLAB. Using the genetic algorithm, the values of and are obtained as 0.727963865 and 84.4065933, respectively. With the initial estimates of and , the reentry epoch of GSLV-D5 is found to be on June 1, 2014 (147 days from January 5, 2014). The results are given in Table 1.
The number of chromosomes is 24 and the mating pool is Mate =
A number of crossovers and mutations were performed on the chromosomes to get the best solution. The fitness or the quality of the solution is determined by where is the average dispersions in apogee, is observed apogee altitude, is the osculating eccentricity, is time instants for observations, is apogee surface, where , and is the number of observations.
3.1.2. Zone B
The osculating orbital elements of GSLV R/B as obtained from the TLE on March 2, 2014, 10 : 45 (UTC), are as follows: semimajor axis (km) = 24087.380496, eccentricity = 0.7279210409, inclination (°) = 19.33807553, argument of perigee (°) = 219.63775080, right ascension of the ascending node (°) = 189.22925587, true anomaly (°) = 140.11376942.
To generate a set of mean apogee surfaces for zone B, three values of initial osculating eccentricity (0.7278210409, 0.7279210409, and 0.7280210409) and three values of ballistic coefficient (50, 80, and 110) are selected to obtain nine grid points as plotted in Figure 5. Using the genetic algorithm, the values of and are obtained as 0.727896214 and 85.8965302, respectively.
With the initial estimates of and , the reentry epoch of GSLV-D5 is found to be on June 1, 2014 (147 days from January 5, 2014). The results are given in Table 2.
3.1.3. Zone C
The osculating orbital elements of GSLV R/B as obtained from the TLE on March 10, 2014, 14 : 32 (UTC), are as follows: semimajor axis (km) = 24031.141029, eccentricity = 0.7277057104, inclination (°) = 19.33698821, argument of perigee (°) = 225.83231865, right ascension of the ascending node (°) = 185.82625378, true anomaly (°) = 133.79218860.
To generate a set of mean apogee surfaces for zone C, three values of initial osculating eccentricity (0.7276057104, 0.7277057104, and 0.7278057104) and three values of ballistic coefficient (50, 80, and 110) are selected to obtain nine grid points as plotted in Figure 6. Using the genetic algorithm, the values of and are obtained as 0.727617443 and 81.6301346, respectively.
With the initial estimates of and , the reentry epoch of GSLV-D5 is found to be on June 2, 2014 (148 days from January 5, 2014). The results are given in Table 3.
3.1.4. Zone D
The osculating orbital elements of GSLV R/B as obtained from the TLE on March 24, 2014, 06 : 31 (UTC), are as follows: semimajor axis (km) = 23888.894645, eccentricity = 0.7268822734, inclination (°) = 19.29960868, argument of perigee (°) = 236.26578322, right ascension of the ascending node (°) = 180.05835079, true anomaly (°) = 121.02982099,
To generate a set of mean apogee surfaces for zone D, three values of initial osculating eccentricity (0.7267822734, 0.7268822734, and 0.7269822734) and three values of ballistic coefficient (50, 80, and 110) are selected to obtain nine grid points as plotted in Figure 7. Using the genetic algorithm, the values of and are obtained as 0.726806283 and 76.4813202 , respectively. With the initial estimates of and , the reentry epoch of GSLV-D5 is found to be on May 29, 2014 (144 days from January 5, 2014). The results are given in Table 4.
From the above four zones, the orbital lifetime of GSLV-D5 is found to be between 144 and 148 days, which is a small variation between the four lifetime values from four different epochs. The near linear variation of mean apogee altitude has shown the reentry time more accurately from the TLEs considered. Hence, it is proven once again that the method based on near linear variation of mean apogee altitude utilized in [6] is providing reasonably good estimates of the orbital lifetime of the rocket body. The genetic algorithm parameters are given in Table 5. The computed values of and , the predicted reentry time, and time interval are given in Table 6.
The justification for using RSM to predict mean apogee altitude, as opposed to using NPOE software to propagate the trajectory forward (starting with the initial values of and possessed by each population member in the GA) to the observed epoch time of the TLE, is the computational savings of interpolation over numerical propagation. While the present methodology is not guaranteed to improve reentry estimates in every scenario, these results indicate that the RSM technique combined with a GA is a promising approach. Moreover, the objects from GTO usually reenter after a long time but this object reenters between 144 and 148 days, making this paper very important for space mitigation efforts. The idea of identifying near linear variation of mean apogee altitude time intervals for reentry time prediction has helped to get better results with less time consumption.
4. Conclusions
The response surface technique with the genetic algorithm is utilized to obtain the optimal values of the initial osculating eccentricity and the ballistic coefficient of each of the selected time intervals based on the near linear variation of mean apogee altitude. Using these optimal values, the orbital lifetime of a GSLV-D5 rocket body is found to be between 144 and 148 days from its injection into the orbit on January 5, 2014.
Nomenclature
| GTO: | Geostationary transfer orbit |
| GSLV: | Geosynchronous satellite launch vehicle |
| LEO: | Low Earth orbit |
| GSO: | Geostationary orbit |
| : | Drag coefficient |
| : | Mass of the object |
| : | Effective area |
| : | Rocket body |
| TLE: | Two-line element |
| GA: | Genetic algorithm |
| NPOE: | Numerical prediction of orbital events |
| SGP: | Simplified general perturbations |
| SDP: | Simplified deep space perturbations |
| : | Solar flux |
| : | Geomagnetic index |
| RSM: | Response surface methodology |
| : | Independent variables |
| : | Random error |
| : | Response |
| : | Evolutionary algorithms |
| : | Osculating eccentricity |
| : | Ballistic coefficient |
| A, B, C, D: | Four zones |
| : | Zonal and tesseral harmonic terms |
| GEM10B: | Goddard Earth model 10B |
| MSIS90: | Mass-Spectrometer-Incoherent-Scatter-1990 atmosphere model |
| NORAD: | North American Aerospace Defense Command |
| UTC: | Coordinated universal time |
| (): | Three values of osculating eccentricity |
| (): | Three values of ballistic coefficient |
| : | First value of osculating eccentricity and ballistic coefficient |
| : | First value of osculating eccentricity and second value of ballistic coefficient |
| : | Average dispersions in apogee |
| : | Observed apogee altitude |
| : | Apogee surface |
| : | Time instants for observations |
| : | Number of observations |
| : | . |
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors are thankful to Dr. T. V. Christy, Director of School of Mechanical Engineering, and Dr. Pradeep Kumar, Head of Aerospace Department, for their constant support and motivation.