Journal of Applied MathematicsVolume 2011 (2011), Article ID 871304, 19 pagesdoi:10.1155/2011/871304

Research Article

Air Drag Effects on the Missile Trajectories

F. A. Abd El-Salam1,2

¹Department of Mathematics, Faculty of Science, Taibah University, Al-Madinah Al-Munawarah, Saudi Arabia²Department of Astronomy, Faculty of Science, Cairo University, Giza 12613, Egypt

Received 24 May 2011; Accepted 3 August 2011

Academic Editor: Michela Redivo-Zaglia

Copyright © 2011 F. A. Abd El-Salam. 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.

The equations of motion of a missile under the air drag effects are constructed. The modified TD88 is surveyed. Using Lagrange's planetary equations in Gauss form, the perturbations, due to the air drag in the orbital elements, are computed between the eccentric anomalies of the burn out and the reentry points [], respectively. The range equation is expressed as an infinite series in terms of the eccentricity e and the eccentric anomaly E. The different errors in the missile-free range due to the drag perturbations in the missile trajectory are obtained.

1. Introduction

The ballistic missile problem is typically represented in astrodynamics as a simplified planar two-body problem with a nonrotating Earth. The missile is only guided during the relatively brief initial and final powered phases of flight, while most of its course, the free flight range, is subsequently governed by the laws of orbital mechanics and ballistics. From this information we can usually determine the range of a projectile, the initial flight path angle, the orbit eccentricity and semimajor axis, and its true/eccentric anomaly angles at launch and impact points. It is of interest to expand upon this basic model in order to add realism to the problem.

In a previous work, Abd El-Salam and Abd El-Bar [1] computed the different errors in the free range of the ballistic missiles taking into account an oblate Earth model retaining the zonal harmonics of the geopotential up to . They obtained explicit expressions for the errors in the missile range due to the in-orbit plane and out-of-orbit plane changes in the missile range. In this work, the author aims to compute the effect of the atmospheric drag perturbations solely based on the usage of the TD88 density model on the missile trajectory. On one hand, the ballistic missile trajectories are very sensitive to any kind of perturbations, but on the other hand the missile trajectory takes approximately 30 minutes (nearly half an only one orbital period of a low Earth orbit, in brief LEO). Therefore, these perturbations can be safely added up linearly to the gravity perturbations on the missile trajectories obtained previously by the authors mentioned above.

The problem of the orbital motion of celestial bodies in a resistive medium can be traced back to Newton's Principle Newton, (bk. II, Sections I/IV) [2] and Laplace's analysis of a possibly finite velocity of gravitation Laplace [3]. However, with the advent of the space age, there arose the need of precisely predicting the orbital behavior of the newly launched bodies. Compared to the gravitational perturbations, the adequate allowance for atmospheric drag turned out to be a problem for unforeseen complexity Brouwer [4]. Besides the obstacles in the analytical treatment of the drag based on the idealized spherical nonrotating density model, the upper atmosphere revealed itself to behave unexpectedly. The irregular density variations derived from orbital analyses of the first launched satellites were attributed to the influence of the Sun, as they displayed periodical changes of approximately 27 days. Jacchia [5]; King-Hele and Walker [6] introduced the overall dependence of the drag on the angle between the perigee and the Sun. In the subsequent years, together with the advances in theory King-Hele [7], direct evidence has grown indicating that the solar and geomagnetic effects strongly influence the thermospheric density both periodically and randomly Jacchia [8–10], Jacchia and Slowey [11, 12].

The literature on the atmospheric drag effects on artificial satellites is very rich; see for example, King-Hele [7], Hoots and Roehrich [13] Gooding [14], Helali [15] Hoots and France [16, 17], Danielson et al. [18], Fonte et al. [19], Sehnal [20–22], Sehnal and Pospíšilová, [23, 24], Rossi et al. [25], and Abd El-Salam and Sehnal [26].

Although the literature on the perturbed missile trajectories is very limited, we can apply the concepts of atmospheric drag effects on artificial satellites on the perturbed missile trajectories with a modification on the integrating limits.

Isaacson and Vaughan [27] described a method of estimating and predicting ballistic missile trajectories using a Kalman Filter over a spherical, nonrotating Earth. They determined uncertainties in the missile launch point and missile position during flight. Vinh et al. [28] obtained a minimum-fuel interception of a satellite, or a ballistic missile, in elliptic trajectory in a Newtonian central force field, via Lawden's theory of primer vector. McFarland [29] treated ballistic missile problem by modeling spherical Earth, rotation, and addition of atmospheric drag using the state transition matrix. Akgül and Karasoy [30] developed a trajectory prediction program to predict the full trajectory of a tactical ballistic missile. Kamal [31] developed an algorithm that includes detection of cross-range error using Lambert scheme in free space, in the absence of atmospheric drag. Bhowmik and Chandrani [32] investigated the advantages and performance of extended Kalman Filter for the estimation of nonlinear system where linearization takes place about a trajectory that was continually updated with the state estimates resulting from the measurement. Forden [33] described the integration of the three degrees of freedom equations of motion, and approximations are made to the aerodynamic for simulating ballistic missiles. Harlin and Cicci [34] developed a method for the determination of the trajectory of a ballistic missile over a rotating, spherical Earth given only the launch position and impact point. The iterative solution presented uses a state transition matrix to correct the initial conditions of the ballistic missile state vector based upon deviations from a desired set of final conditions. Kamal [35] presented an innovative adaptive scheme which was called “the Multistage Lambert Scheme”. Liu and Chen [36] presented a novel tracking algorithm by integrating input estimation and modified probabilistic data association filter to identify warhead among objects separation from the reentry vehicle in a clear environment.

2. Air Drag Force Model

As an artificial Earth satellite moves through the Earth’s atmosphere, its orbital motion is perturbed by the resisting force. The main effects of the force are to bring out large secular decreases in the elements , the semimajor axis, and , the eccentricity, that cause the satellite plunges toward Earth. This perturbing influence is one of the most important perturbing forces in the altitude regime from 150 to 600 km, which is called the drag force. The drag force, , per unit mass of the satellite can be represented by where is a dimensionless drag coefficient, is the velocity of the satellite relative to the atmosphere (called the ambient velocity), represents the effective cross-sectional area of the satellite, which is to be found by averaging all possible projected areas of the satellite onto a plane perpendicular to , is the mass of the satellite, and is the density of the atmosphere. The density depends on position and time in a very complicated manner. Its variation with height depends on the form of the atmosphere, since it decreases rapidly with slight growing in the altitude, the effect of the flattening of the atmosphere is rather significant.

3. Rotation of the Atmosphere

This effect is so small, so a simple approximation is justified. The differential rotation of the atmosphere complicates the numerical calculation, and it gives no sense in this very small effect. So we will assume that the atmosphere rotates with the same angular velocity of the Earth , then where is the orbital velocity of the satellite and its radius vector. If are the base vectors of the rectangular geocentric equatorial frame associated with the spherical coordinates , then  rad./sec. Thus where is the orbital inclination, is the satellite's declination, is the integral of areas, and is given by The third term in may usually be neglected, being less than about , and the second term which is of order of may be evaluated at perigee, where the density is much greater than elsewhere on the orbit, thus we may effectively regard as a constant.

4. Drag Force Components

Let be unit vectors along the radius vector, normal to the radius vector in the orbital plane, and normal to the orbital plane, that is, along to the angular momentum axis, respectively; then the drag force given by (2.1) can be written in the components form as where

4.1. Modified TD88 Density Model

To integrate the equations of motion, an analytical expression is required for the density of the atmosphere. This density has many irregular and complex variation both in position and time. it is largely affected by solar activity and by the heating or cooling of the atmosphere. This time dependence is very difficult to be introduced in an analytical expression. Since the atmosphere is not actually spherically symmetric but tends to be oblate, we have to allow for this oblateness in any expression for the air density. The complications which come from the accuracy demands may cause the impossibility to use these classical aeronomic empirical models, for example, CIRA 86 Hedien [37], CIRA 72 Oliver [38], and DTM Barlier et al., [39] for an analytical work. In this work, we will use the modified version of the TD88 model first proposed and introduced by Sehnal [20, 22] and Sehnal and Pospíšilová [23] and then developed by Abd El-Salam and Sehnal [26] in terms of the eccentric anomaly as where , the main solar and geomagnetic activity parameters affecting the angular density, are given by with solar flux measured on 10.7 cm wavelength for a previous day, solar flux averaged over three solar rotations (81 days), geomagnetic index 3 hours before local time , the reciprocal of the density scale height , the equatorial radius of the Earth, the Earth's flattening, , , and numerical constants of the model given by Tables 1, 2, and 3.

Table 1: Constants .

i
1
2
3
4
5
6
7
8

0.007
0.2875
0.04762
0.0471
7.0
7.0
0.3333
15.0

Table 2: Constants .

n

1

2

3

4

5

6

7

Table 3: Phases .

n
3
4
5
6
7

263
−263
−29.41
8.0913
10.0813

The required nonvanishing coefficients for the model are where the nonvanishing coefficients , , , and were computed previously by the author, see Abd El-Salam and Sehnal [26].

5. Integration of the Equations of Motion

The Lagrange's planetary equations are a very well-known system of six first-order differential equations. They are used to evaluate the perturbations in the usual Keplerian orbital elements . Gauss modified this system so as to be applicable for all kinds of force especially nonconservative force, see Roy [40], with , , where are the mean motion in the orbit, the semi-latus rectum of the orbit, and the gravitational parameter, respectively, and are the drag force components given by (4.2).

To integrate the equations of motion (5), we first do the following simplifying steps:(1)introducing the density given by the modified TD88 model to the drag force components ;(2)using , to replace the explicit time with the eccentric anomaly as an independent variable;(3)using the binomial theorem and neglecting orders of and ;(4)since the series expansion of , the real numbers, is convergent, then we replace appeared in the density expression with its equivalent series expansion;(5)let us define.

5.1. The Increment in the Semimajor Axis

Using the first equation from the equations set (5) and following the above-mentioned simplifying steps, then after some lengthy algebraic manipulations we obtain the change in the semimajor axis as where the nonvanishing coefficients are using the trigonometric identity where is the usual Kronecker Delta symbol and is the greatest integer function, , and is the binomial coefficients. Now (5.4) can be written in the form Now (5.5) can be integrated between the eccentric anomalies of the burn out and the reentry points , respectively, then we have integrals of the form from which we can get Similarly, we can evaluate which yields where terms involving the integrals of the type within the expression (5.5) have no contribution to the problem. as usual to the drag theory of artificial satellites. Now substituting the relevant integrals from (5.7) into (5.5) yields the solution as Using the similar procedure to integrating the rest of (5) between the eccentric anomalies of the burn out and the re entry points , respectively, we obtained the following explicit expressions for the other orbital elements as follows.

5.2. The Increment in the Eccentricity

where the nonvanishing coefficients are:

5.3. The Increment in the Inclination

where the nonvanishing coefficients and are

5.4. The Increment in the Ascending Node

where the nonvanishing coefficients and are

5.5. The Increment in the Argument of the Perigee

where the nonvanishing coefficients are

5.6. The Increment in the Eccentric Anomaly

where the nonvanishing coefficients are

6. Drag Effect on the Free Flight Range

The usual ICBM problem concerns the determination of the free flight range angle. In order to improve the results, we will take into account the perturbations in the orbital elements coming from the atmospheric drag. Following the treatment given by Helali [15] the free flight range equation reads where defines the dimensionless parameter, is the magnitude of the position vector of the missile relative to the Earth's center, is the missile speed at any point in its orbit, is the corresponding missile circular speed at this point, and is the gravitational parameter. From the orbital mechanics of the two body problem and the symmetry of the free-flight portion (see [41, Figure and Figure , pages 283–285]) we can express the range equation in terms of the eccentricity and the eccentric anomaly as follows: Since the range equation is a function the eccentricity and the eccentric anomaly and and denote the variations in the free flight range angle due to the variation of the eccentricity and due to the variation of the eccentric anomaly, respectively, under the effect of the atmospheric drag, therefore, we can write The involved partial derivatives in (6.3) are given by

6.1. Computation of the Quantities and

Substituting from (5.11) and (6.4) into the first term of (6.3) yields

Substituting again from (5.11) and (6.4) into the second term of (6.3) yields

6.2. Errors in due to

For the sake of the simplicity, let us take the burnout point on the equator. If the reentry point, for some reason, was displaced by an amount, , and if all of the orbital elements were kept fixed except the longitude of the ascending node, then this displacement could be interpreted as a change in the longitude of the ascending node . Due to this change, a cross-range error, , at impact occurs, the value of which is obtained by applying the law of cosines in the spherical trigonometry. Hence, Since and are small angles, then we have Substituting (5.15) into (6.8) yields

6.3. Errors in due to

Again assume that the burnout point were on the equator and assume also the actual launch azimuth differs from the intended value by an amount, , and if all of the orbital elements were kept fixed except the orbital inclination, then this displacement could be interpreted as a change in the orbital inclination . Due to this change a cross-range error, , at impact occurs, the value of which is obtained in a similar manner as the previous subsection by applying the law of cosines in the spherical trigonometry. Hence, Since and are small angles, then we have Substituting (5.13) into (6.11) yields

7. Conclusion

In this work, the equations of motion of a missile are constructed. Using Lagrange's planetary equations in Gauss form, the perturbations due to the air drag in the orbital elements are computed between the eccentric anomalies of the burn out and the re entry points , respectively. The range equation is expressed in terms of the eccentricity and the eccentric anomaly . The different errors in the missile free range due to the drag perturbations in the missile trajectory are obtained. These obtained perturbations can be safely added up linearly to the gravity perturbations on the missile trajectories obtained previously by Abd El-Salam and Abd El-Bar [1] to constitute a more precise theory of the missile trajectories and the errors in the free flight range.

Acknowledgments

The author is deeply indebted to the team work at the deanship of the scientific research—Taibah University for their valuable help and critical guidance and for facilitating many administrative procedures. This research work was financed (supported) by Grant no. (627/1431) from the deanship of the scientific research at Taibah university, Al-Madinah Al- Munawwarah, Saudi Arabia.

References

F. A. Abd El-Salam and S. E. Abd El-Bar, “Computation of the different errors in the ballistic missiles range,” ISRN Applied Mathematics, vol. 2011, Article ID 349737.

I. Newton, Philosophiae Naturalis Principia Mathematica, 1987, English translation by F. Cajori, Newton's Principia, University of California Press, Berkeley, 1934.

P.S. Laplace, Traite' de Mécanique Céleste, vol. 4, part 2, Courcier, Paris, France, 1805.

D. Brouwer, “Review of celestial mechanics,” Annual Review of Astronomy and Astrophysics, vol. 1, pp. 219–234, 1963.

L. G. Jacchia, “Two atmospheric effects in the orbital acceleration of artificial satellites,” Nature, vol. 183, no. 4660, pp. 526–527, 1959.

D. G. King-Hele and D. M. C. Walker, “Irregularities in the density of the upper atmosphere: results from satellites,” Nature, vol. 183, no. 4660, pp. 527–529, 1959.

D. King-Hele, Theory of Satellite Orbits in an Atmosphere, Butterworth, London, Uk, 1964.

L. G. Jacchia, “Solar effects on the acceleration of artificial satellites,” SAO Special Report 29, 1959.

L.G. Jacchia, “A variable atmospheric-density model from satellite accelerations,” SAO Special Report 39, 1960.

L. G. Jacchia, “The atmospheric drag of artificial satellites during the October 1960 and November 1960 events,” SAO Special Report 62, 1961.

L. G. Jacchia and J. Slowey, “Preliminary analysis of the atmospheric drag of the twelve-foot balloon satellite (1961 d1),” SAO Special Report 84, 1962.

L. G. Jacchia and J. Slowey, “Accurate drag determinations for eight artificial satellites: atmospheric densities and temperatures,” SAO Special Report 100, 1962.

F. R. Hoots and R. L. Roehrich, “Models for propagation of NORAD element sets,” Project Spacecraft Report no. 3, Aerospace, 1980.

R. H. Gooding, “On the generation of satellite position (and velocity) by a mixed analytical-numerical procedure,” Advances in Space Research, vol. 1, no. 6, pp. 83–93, 1981.

Y. E. Helali, “Perturbations in the perigee distance due to atmospheric drag for artificial earth satellites,” Bulletin of the Astronomical Institutes of Czechoslovakia, vol. 38, pp. 329–334, 1987.

F. R. Hoots and R. G. France, “An analytic satellite theory using gravity and a dynamic atmosphere,” Celestial Mechanics, vol. 40, no. 1, pp. 1–18, 1987.

F. R. Hoots and R. G. France, “The future of artificial satellite theories: hybrid ephemeris compression model,” Celestial Mechanics and Dynamical Astronomy, vol. 66, no. 1, pp. 51–60, 1996.

D. A. Danielson, C. P. Sagovac, B. Neta, and L. W. Early, “Semi analytic satellite theory (SST): mathematical algorithms,” Tech. Rep. NPS-MA-95-002, Naval Postgraduate School, Monterey, Calif, USA, 1995.

D. J. Fonte, B. Neta, C. Sabol, D. A. Danielson, and W. R. Dyar, “Comparison of orbit propagators in the research and development goddard trajectory determination system (R&D GTDS). Part I: simulated data,” in Proceedings of the American Astronomical Society/American Institute of Aeronautics and Astronautics (AAS '95), vol. 431, 1995.

L. Sehnal, “Thermospheric total density model TD,” Bulletin of the Astronomical Institutes of Czechoslovakia, vol. 39, pp. 120–127, 1988.

L. Sehnal, “Theory of the motion of an artificial satellite in the earth atmosphere,” Advances in Space Research, vol. 10, no. 3-4, pp. 297–301, 1990.

L. Sehnal, “Comparison of the thermosphere total density model TD 88 with cira 86,” Advances in Space Research, vol. 10, no. 6, pp. 27–31, 1990.

L. Sehnal and L. Pospíšilová, “Thermospheric model TD88,” Preprint No. 67 of the Astronomical Institute of Czechoslovakian Academy of Sciences, 1988

L. Sehnal and L. Pospíšilová, “Lifetime of the ROHINI A satellite,” Bulletin of the Astronomical Institutes of Czechoslovakia, vol. 42, pp. 295–297, 1991.

A. Rossi, L. Anselmo, A. Cordelli, P. Farinella, and C. Pardini, “Modelling the evolution of the space debris population,” Planetary and Space Science, vol. 46, no. 11-12, pp. 1583–1596, 1998.

F. A. Abd El- Salam and L. Sehnal, “An analytical second order atmospheric drag theory based on usage TD 88 thermospheric density model,” Celestial Mechanics and Dynamical Astronomy, vol. 90, pp. 361–389, 2004.

J. A. Isaacson and D. R. Vaughan, Estimation and Prediction of Ballistic Missile Trajectories, RAND Corp, Santa Monica, Calif, USA, 1996.

N. X. Vinh, T. Kabamba, and T. Takehira, “Optimal interception of a maneuvering long-range missile,” Acta Astronautica, vol. 48, no. 1, pp. 1–19, 2001.

J. S. McFarland, “Modeling the ballistic missile problem with the state transition matrix: an analysis of trajectories including a rotating earth and atmospheric drag,” Semester Report, Virginia Tech, Blacks-burg, Va, USA, 2004.

A. Akgül and S. Karasoy, “Development of a tactical ballistic missile trajectory prediction tool,” Istanbul University Journal of Electrical and Electronics Engineering, vol. 5, no. 2, pp. 1463–1467, 2005.

S. A. Kamal, “Cross range error in the Lambert Scheme,” in Proceedings of the 10th National Aeronautical Conference, pp. 255–263, College of Aeronautical Engineering, PAF Academy, Risalpur, NWFP, Pakistan, April 2006.

S. Bhowmik and R. Chandrani, “Application of extended Kalman filter to tactical ballistic missile re-entry problem,” http://arxiv.org/abs/arXiv:0707.1918

G. Forden, “GUI missile flyout: a general program for simulating ballistic missiles,” Science and Global Security, vol. 15, pp. 133–146, 2007.

W. J. Harlin and D. A. Cicci, “Ballistic missile trajectory prediction using a state transition matrix,” Applied Mathematics and Computation, vol. 188, no. 2, pp. 1832–1847, 2007.

S. A. Kamal, “The multi-stage-lambert scheme for steering a satellite- launch vehicle (SLV),” in Proceedings of the 12th IEEE International Multitopic Conference, pp. 294–300, Bahria University, Karachi, Pakistan, December 2008.

C. Y. Liu and C. T. Chen, “Tracking the warhead among objects separation from the reentry vehicle in a clear environment,” Defence Science Journal, vol. 59, no. 2, pp. 113–125, 2009.

A. E. Hedin, “MSIS-86 Thermospheric Model,” Journal of Geophysical Research, vol. 92, no. A5, pp. 4649–4662, 1987.

CIRA 72, COSPAR Intl, Reference Atmosphere, Academic Verlag, Berlin, Germany, 1972.

F. Barlier, C. Berger, J. L. Falin, G. Koskarts, and G. Thuillier, “A thermospheric model based on satellite drag data,” Annales Geophysicae, vol. 34, p. 9, 1978.

A. E. Roy, Orbital Motion, Adam Hilger, Bristol, UK, 1978.

R. R. Bate, D. D. Mueller, and J. E. White, Fundamentals of Astrodynamics, Dover Publications Inc., New York, NY, USA, 1971.