Abstract: The solar activity is known to have a significant influence on the upper atmosphere and near Earth environment. When solar activity is high, ultraviolet radiation from the Sun heats and expands the Earth’s upper atmosphere which in turn, increases the atmospheric drag and orbital decay rate of the satellites. In the present study, the variation of the orbital elements of the satellites in Low Earth Orbit (LEO) environment due to the effect of the air drag force has been calculated. The variation of the rotation velocity of the atmosphere is considered. The Lagrange Planetary Equation in Guass form is used to calculate the variation of the orbital elements numerically. The United States Naval Research Laboratory Mass Spectrometer and Incoherent Scatter Radar Exosphere (NRLMSISE-00) model which maps out the atmosphere from sea level to 1000 km is used to estimate the density of the atmosphere. As a numerical application, the Hubble Space Telescope (HST), CORONAS I and PRIRODA satellites are taken to express the LEO environment. It is found that the difference in the variation of the orbital elements during the maximum and minimum solar activity reaches three orders in magnitude and that the air drag force has a larger impact on the satellites during the solar maximum epoch than during the minimum activity epoch.
INTRODUCTION
The drag force is one of the major perturbing forces on satellites in low earth orbit where the atmospheric density is the most variable element. A Low Earth Orbit (LEO) satellite is defined as a satellite that orbits below an altitude of approximately 1000 km (Samwel, 2014). Although the air density is much lower than that near the earth’s surface, the air resistance in those layers of the atmosphere, especially at altitudes ranging from 300-800 km above the earth’s surface, is still strong enough to produce drag and pull the satellites closer to the earth. Hence, disregarding the geopotential perturbing forces, the air drag force is considered to be the predominant among other perturbing forces as stated by Eshagh and Alamdari (2007).
Many efforts were made in the past to model the air drag using different mathematical solution techniques. Brouwer and Hori (1961) concerned the presence of the oblateness perturbations and the drag effect in a single solution. They used the canonical variables which have particular advantages in the treatment of this problem because of the ease with which the necessary transformation of variables can be made. They were the first to treat the joint effects of drag and earth oblateness. Thereafter, Delhaise (1991), developed an analytical solution for the motion of artificial earth satellites subjected to the combined effects of earth oblatness and atmospheric drag.
Khalil (2002) has proposed a theory of an Artificial Earth Satellite under the joint effects of earth oblateness and atmospheric drag. The Hamilton’s equations of motion are derived and the atmospheric model was taken as an oblate rotating model. Also, Saad et al. (2008) used KS regular variables to predict orbital lifetime numerically. They considered the perturbations due to the non-aspherecity of the earth and the atmospheric drag.
The impact of the solar activity variation on the satellite decay in LEO was addressed in many studies. Sehnal (1975) paid attention to the difficulty to obtain a precise analytic continuous solution due to the possible and unexpectable change of the density with the activity of the sun and due to the rapid increase of the drag constant above 400 km from the surface of the earth. Pardini et al. (2006) analyzed the orbital decay of eleven spherical satellites in the 150-1500 km altitude range over a full solar activity cycle to highlight the performances of the JR-71 and MSISE-90 density models. They concluded that the two models were not able to correctly represent the air density at each altitude and environmental condition. They found that below 400 km, both models overestimated the atmosphere density. However, MSISE-90 was proved to be the best model to compute the air density, in low solar activity conditions while JR-71 was generally more accurate over the maximum of the 23rd solar cycle. Lately Nwankwo and Chakrabarti (2013) studied the effects of the adverse space weather condition induced by solar events, especially the solar flares and Coronal Mass Ejections (CMEs) on the orbital decay of low earth satellite orbits. They showed that the result strongly depends on the phase of the solar cycle. They also showed that a major CME event can cause sufficient heating and expansion of the atmosphere so that the orbital radius may go down by a few km in a single event.
The unwanted orbital changes occur due to the drag force are known to be related to atmospheric density changes which are primarily induced by solar extreme ultraviolet irradiance variations (Tobiska, 2003). This ambient atmosphere is considered to be important environment at the LEO and polar earth orbit altitudes. For the other orbits, it is too low to be of significance (Hasting and Garrett, 1996).
Therefore, in the present study, the variations of the orbital elements which define the orbital trajectory of the satellite due to the effect of the air drag force at LEO environment during maximum and minimum solar activity are evaluated. The variation of the rotation velocity of the atmosphere is considered. The Lagrange planetary equation in Guass form is used to calculate the variation of the orbital elements numerically and the NRLMSISE-00 model which maps out the atmosphere from sea level to 1000 km is used to estimate the density of the atmosphere.
As a numerical application, the orbital elements of the Hubble Space Telescope (HST), CORONAS I and PRIRODA satellites, which express the Low Earth Orbit (LEO) satellites are used. The variation of the orbital elements of these satellites under the effect of the drag force for both maximum and minimum solar activity are calculated.
METHODOLOGY
Atmospheric Density model (NRLMSISE-00): The NRLMSISE is an empirical, global model that models temperatures and densities of atmospheric components from the ground to the exosphere. It is developed by the U.S. Naval Research Laboratory (NRL) and based on the data of mass spectrometer and Incoherent Scatter Radar (Hedin, 1991; Picone et al., 2002). The main advantage of the usage of these datasets is that they consist of independent observations of both temperature and number densities for the atmospheric constituents. The MSIS-86 model replaced Jacchia-71. An extension to the model was published as MSISE- 90 (Hedin, 1991), which is identical to MSIS-86 for the thermosphere region, but extends down to zero height. At the end of the 1990s, NRLMSISE-00 model was released (Picone et al., 2002) including additional mass spectrometer and incoherent scatter radar data, accelerometer data, as well as database of Drag Temperature Model (DTM) and Jacchia satellite orbit decay (Doornbos, 2012).
The NRLMSISE-00 comprises the main drivers of the upper atmosphere: The solar extreme ultraviolet (EUV) flux and geomagnetic heating. The 10.7 cm solar radio flux (F10.7) is an indicator of the overall solar activity level. It measures the radio flux at a wavelength of 10.7 cm in units of 10–22 W m–2 Hz–1 and is used as a proxy of the extreme ultraviolet radiation (Samwel, 2014) while the Ap daily geomagnetic index measures the geomagnetic component of space weather. The F10.7 index is obtained from solar data at the National Geophysical Data Center (NGDC) http://www.swpc.noaa.gov/ and The National Oceanic and Atmospheric Administration world data center for geomagnetism in Kyoto http://wdc.kugi.kyoto-u.ac.jp/.
The NRLMSISE-00 has also been compared to the previous MSIS models and Jacchia-70 and represented a noticeable improvement (Picone et al., 2002). In the present study, the NRLMSISE-00 model for density calculations is used. The variation of the total mass density (ρ) with the true anomaly (f), during the maximum and minimum solar activity, is obtained and fitted using Fourier series of degree 4 as represented in Eq. 1. The values of the fitting coefficients (a0, ai and bi) are outlined in Table 2:
| (1) |
Drag force: The aerodynamic drag force
| (2) |
where,
The Ballistic number (B) which is obtained from the TLE outlined in Space-Track website, is used (https://www.space-track.org/) where:
It is usually assumed that the atmosphere rotates with the same angular velocity
| (3) |
which realizes that as:
and as:
If the satellite is at
| (4) |
If
| (5) |
| (6) |
where:
| (7) |
and:
| (8) |
where,
In the case of Lagrange’s planetary equations, the perturbing force decomposed into radial, transverse and normal directions is needed. The unit vector
| (9) |
with:
| (10) |
where, μ is the earth’s gravitational parameter, a is the semi-major axis, e is the eccentricity of the satellite orbit and the transverse unit vector
| (11) |
| (12) |
and:
| (13) |
The drag force in the radial, transverse and normal components (R, T, N) are:
| (14) |
| (15) |
| (16) |
where:
| (17) |
| (18) |
| (19) |
| (20) |
| (21) |
| (22) |
where, i, ω, f are the inclination of the orbit on the equator, argument of the perigee and the true anomaly of the spacecraft orbit, respectively.
Lagrange planetary equations with air drag force: The Lagrange planetary equations are useful only for evaluating forces that can be presented in terms of a potential. However, there is a Gaussian form of these equations that can be used to study the perturbing forces (air drag, in our case) (Brouwer and Clemence, 1961).
The evaluation of such perturbing forces analytically can become quite complex. Therefore, we used the numerical integration to calculate the variation of the orbital elements due to air drag force. The Lagrange’s planetary equations in Gauss form are:
| (23) |
| (24) |
| (25) |
| (26) |
| (27) |
| (28) |
where, R, T, N are defined by Eq. 14-16, Ω is longitude of ascending node, M is mean anomaly and n is mean motion:
The above system of differential equations (Eq. 23-28) numerically is solved, using the fourth order Runge-Kutta method to get variations of orbital elements (Δa, Δe, Δi, ΔΩ, Δω, ΔM) as functions of true anomaly, due to the impact of the drag force. Introducing the true anomaly f as follows:
| (29) |
RESULTS AND DISCUSSION
In the present section, the total atmospheric density using the NRLMSIS-00 model during both maximum and minimum solar activity is calculated, taking three LEO satellites to express the LEO environment at different altitudes. In addition, the variations of the orbital elements due to the impact of air drag force are calculated.
Atmospheric density calculations: The NRLMSIS-00 is used to calculate the total atmospheric density for the maximum and minimum solar activity for the Hubble Space Telescope (HST), CORONAS I and PRIRODA satellites, as a representative of the LEO satellites.
The NRLMSIS-00 model needs a time input to reflect the current state of the atmosphere. It also accounts for the fluctuations in the atmosphere from latitude, longitude, F10.7 (Previous day), F10.7 (81 day average) and the geomagnetic index Ap.
The Two Line Elements (TLE) of the three satellites are obtained from Space-Track website https://www.space-track.org/ and their orbital elements are outlined in Table 1 during maximum and minimum solar activity.
The maximum level of solar activity employed for this simulation, is taken to start at the orbital epoch 2000, March 31. Thus, considering F10.7 (Previous day) = 205.1×10–22 Wm–2 Hz–1, F10.7 (81 day average) = 189.1593 and Ap = 23.
The Minimum level of solar activity employed for this simulation, is considered to start at the orbital epoch 1997, March 31. Thus, considering F10.7 (Previous day) = 73.8×10–22 Wm–2 Hz–1, F10.7 (81 day average) = 73.56049 and Ap = 5.
Figures 1 and 2 draw the total atmospheric density variation with the true anomaly through one orbital revolution for the three satellites during solar maximum and minimum respectively, using the NRLMSIS-00 model through the SPENVIS website https://www.spenvis.oma.be/. The figures represent that the total atmospheric density is higher during the epoch of maximum solar activity than that during minimum solar activity and that it is increasing with increasing altitude.
The data obtained is fitted using Eq. 1, The values of the fitting coefficients obtained in Eq. 1 for the three satellites are outlined in Table 2.
Orbital element variations due to air drag force impact: In the present sub-section, the variation of orbital elements (Δa, Δe, Δi, ΔΩ, Δω, ΔM), due to air drag force, as function of true anomaly are obtained by solving Eq. 23-28 numerically for the three satellites mentioned above. The orbital element variations are calculated during both maximum and minimum solar activity.
The variation in the semi major axis and eccentricity due to the air drag force during the maximum and minimum solar activity for the three satellites, HST, CORONAS-I and PRIRODA are represented in Fig. 3-5, respectively. Figure 3 represents the orbital variation of two orbital elements, namely, the semi major axis (a) and eccentricity (e) of the HST, for one revolution, during both maximum and minimum solar activity. In addition, Fig. 4 and 5 represent the variation in the same orbital elements for the satellites CORONAS-I and PRIRODA, respectively. The orders of magnitude of Δi, ΔΩ, Δω and ΔM are outlined in Table 3.
| Table 1: | Orbital elements of HST, CORONAS I and PRIRODA during maximum and minimum solar activity |
| Fig. 1(a-c): | Total atmospheric density variation with the true anomaly through one orbital revolution for (a) HST, (b) CORONAS-I and (c) PRIRODA satellites during the maximum solar activity |
| Table 2: | Coefficients of the fitting using Fourier series of the total mass density for the three satellite HST, CORONAS-I and PRIRODA |
Eshagh and Alamdari (2007), in their study, represented that the air drag force is reducing the semi major axis of the satellite CHAMP by 6 m in one revolution which is comparable to the result obtained in the present study for the satellite CORONAS-I at the epoch of maximum solar activity. The variation of the other orbital elements due to the air drag force obtained by Eshagh and Alamdari (2007), are nearly similar to the values obtained in the present study as represented in Table 3 and Fig. 4. They used the direct on-orbit measurements of the non-gravitational satellite acceleration instead of the air density models and used the high-low technique of satellite-to-satellite tracking to determine the precise orbit of CHAMP which was launched in summer 2000 into a circular polar orbit at an altitude of about 454 km. As it is represented in Table 1, satellite CORONAS-I is orbiting the earth in a nearly polar orbit with an altitude of about 450 km during the epoch of solar maximum.
From Table 3 and Fig. 3-5, it is found that the differences in the variations of the orbital elements at maximum and minimum solar activity are about three orders of magnitude.
| Fig. 2(a-c): | Variation of the total atmospheric density with the true anomaly through one orbital revolution for (a) HST, (b) CORONAS I and (c) PRIRODA satellites during the minimum solar activity |
| Fig. 3(a-d): | Variations of the semi major axis and eccentricity during (a-b) Maximum solar activity and (c-d) Minimum solar activity for HST |
| Fig. 4(a-d): | Variations of semi major axis and eccentricity during (a-b) Maximum solar activity and (c-d) Minimum solar activity for CORONAS I satellite |
| Fig. 5(a-d): | Variations of semi major axis and eccentricity during (a-b) Maximum solar activity and (c-d) Minimum solar activity for PRIRODA satellite |
| Table 3: | Order of magnitude of (Δi, ΔΩ, Δω and ΔM) due to the air drag force during maximum and minimum solar activity for the three satellites, HST, CORONAS I and PRIRODA |
This difference may be interpreted as when the sun is active, the sun adds extra energy to the atmosphere, hence the low density layers of the atmosphere at the LEO altitudes rise and are replaced by higher density layers that were previously at lower altitudes. As a result, the satellite flies through higher density layer and experiences a stronger drag force.
In addition, Fig. 3-5 and Table 3 clarify that the variation of the orbital elements are much larger in a case of low altitude satellites such as PRIRODA, incomparable to higher altitude satellites, such as CORONAS I and HST (Table 1). This can be interpreted that the density of the atmosphere decreases exponentially with altitude and consequently the air drag force has a larger impact on satellites orbiting in lower altitude that at higher altitudes.
CONCLUSION
In the present study, The variations of the orbital elements of the satellites in Low Earth Orbit (LEO) environment due to the effect of the air drag force during maximum and minimum solar activity are addressed. The variation of the rotation velocity of the atmosphere was considered. The Lagrange’s planetary equations in the Gauss form to analyze the time rates of classical orbital elements resulting from the total air drag force is employed. The (NRLMSISE-00) model was used to estimate the density of the atmosphere. As a numerical application, HST, CORONAS I, PRIRODA were taken to express the LEO satellites.
It is found that:
| • | The variations of the orbital element during the epoch of the maximum solar activity is much larger than those during the epoch of the minimum solar activity. These differences in the variations reach three orders of magnitude |
| • | Due to the air drag force, the variations of the orbital elements for a low altitude satellite is much larger than those of higher altitude satellites |
Therefore, the present study shows that the air drag force has a larger impact on the satellites during the maximum solar activity and on lower altitude satellites.
The present study, is a preliminary investigation of the orbital element variations for one revolution and will be extended for more than one revolution for a comprehensive study for a future work.