Re-entry aerodynamics derived from space debris trajectory analysis

Re-entry aerodynamics derived from space debris trajectory analysis

00324633/92 %5.00+0.00 Pergamon Press Ltd Phet. Sprrre &I., Vol. 40. No. 5. PP. 641446. I992 Prinled in Great Britain. RE-ENTRY AERODYNAMICS DERIVED...

535KB Sizes 5 Downloads 54 Views

00324633/92 %5.00+0.00 Pergamon Press Ltd

Phet. Sprrre &I., Vol. 40. No. 5. PP. 641446. I992 Prinled in Great Britain.

RE-ENTRY AERODYNAMICS DERIVED FROM SPACE DEBRIS TRAJECTORY ANALYSIS DRA Space Department,

R. CROWTHER Royal Aerospace Establishment, Farnborough,

U.K.

(Received in~~~lfo~m 20 November 1991) Abstract-This paper considers the technique of orbital analysis as a means of determining the ill-defined gas-surface interaction between spacecraft and atmospheric molecules in low Earth orbit. The interaction is a major uncertainty in trajectory predictions for a body moving within an atmosphere. The rate of change of the orbital period of a debris object, the uncontrolled Salyut ‘i/Kosmos 1686 space station, is analysed in order to determine the free molecular drag coefficient. The results are compared with theoretical values for the drag coefficient calculated using a complex representation of the vehicle configuration and motion and applying the Monte Carlo Test Particle method. Results suggest a nature of re-emission very close to the classical diffuse, totally accommodated case was occurring at the surface of the debris object as it approached re-entry. However the determined drag coefficient and therefore the derived interaction are found to be very sensitive to the neutral density and therefore the atmospheric model used in the analysis.

1. INTRODUfftON

The large number of uncontrolled debris objects in orbit about the Earth pose a collision hazard to operational vehicles in intersecting orbits. Additionally larger objects in low Earth orbit (LEO) such as Skylab (1973-2714) present a threat resulting from impact on the surface of the Earth because of the appreciable mass of the vehicle which can survive re-entry. One such vehicle was the S&yut 7 (1982-33~)~~o~~o~ 1686 (198%86A) combination which re-entered the Earth’s atmosphere in February 1991 (Crowther, 1991). When it became clear that the orbit of the vehicle would soon decay and that it would re-enter the atmosphere, a trajectory prediction campaign was initiated by the relevant agencies (Klinkrad, 1991). The orbit decay and final re-entry time were difficult to predict because the major perturbing influence on the orbit, aerodynamic drag, is difficult to characterize. This is primarily due to the density variation around the orbit and the gas-surface interaction (between atmospheric molecules and the spacecraft surface) being ill defined. However it was realized that because the motion of the complex was not controlled, it would be possible to measure the perturbation due to the aerodynamic forces acting on the vehicle directly and therefore gain some insight into the nature of the perturbing processes occurring. The approach considered in this paper observes the motion of the centre of mass of the vehicie via orbital analysis. This technique assumes that given a series of accurateiy determined satellite orbits, it is possible to relate the changes in the elements representing the satellite orbit to the

aerodynamic forces acting on the vehicle. This type of analysis has been applied to the orbits of a number of vehicles, in particular the three-axis stabilized satellite ANS-1 (1974-70A) (Sehnal, 1983; Crowther and Stark, 1991). Analysis of the motion of ANS-1 used indirect measures of the satellite drag and lift coefficients derived from the orbital semi-major axis and inclination evolution. This course of action was not possible with the Sulyut l/Kosmos 1686 complex because the motion of the station about its centre of gravity resulted in no residual lift force. Therefore in this paper the orbital period, or more exactly the rate of change of orbital period, is analysed to determine the drag coefficient of the vehicle from which the gassurface interaction can be derived.

2. DRSCRIF’TIONOF THE VEHICLE Sulyut 7 (1982-334) was launched on 19 April 1982. It was joined with Cosmos 1686 (1985-864) on 2 October 1985. The complex had an overall mass of 40,150 kg and a combined length of 27 m. Salyut 7 had three solar arrays attached to it while the Kosmos vehicle had two, each 20m2 in area (see Fig. 1). During the period of observation the attitude of the vehicle was gravity gradient stabilized with Kosmos towards the centre of the Earth. There was a nominal rotation about the longitudinal axis of the vehicle and a precession around the stable radial direction at an estimated half cone angle of 20”. Radar observations (MehrhoIz and Magura, 1991) suggested that the vehicle remained in this orientation until shortly after 641

642

R. CROWTHEK 3. FREE MOLECULAR

AERODYNAMICS

The flow regime encountered by the Salyut complex during the period of the analyses is that of free molecular flow. The implications of this regime are : (i) the gas molecules have a Maxwellian distribution of thermal velocity superimposed upon their uniform mass velocity ; (ii) no boundary layer or shock regions are formed; (iii) the aerodynamic forces acting on a body are a direct result of the momentum exchange between the atmospheric species and the body; (iv) the incoming and reflected molecules and atoms do not influence each other.

I

Earth centre

The first point allows us to define the velocity distribution of particles impinging on the vehicle, the latter three points allow us to use the Monte Carlo Test Particle method to calculate the aerodynamic coefficients. The main steps in this method are to represent the spacecraft configuration as a combination of elementary surfaces such as cones, cylinders and flat plates which can be represented mathematically by an expression of the form :

FIG. 1. THESalyui f/Kosmos I686 COMBINATION.

fyz+gx+hytiz+k

= 0,

(1)

i.e.

the observation period at which time a slow complex tumble became the dominant motion. The orbital inclination remained at a nominal 51.6’, the orbital mean motion was of the order of 15.9 rev day-’ and the orbital eccentricity -0.0003. The rotation about the longitudinal body axis and the precession around the radial direction cause several problems when trying to calculate the aerodynamic coefficients of such a complex configuration as the Sulyut complex. In order to relate the change in the orbital period to the aerodynamic force acting on the vehicle one must integrate the aerodynamic coefficients about the spin and precession cycles. This is a difficult process as different surfaces of the vehicle will be exposed or shielded from the flow by other surfaces at different points around the spin and precession cycles. In addition concave surface elements such as the corner represented by the solar arrays and body of the Salyut station will experience multiple reflections where atmospheric molecules striking one surface of the spacecraft may upon re-emission strike another surface of the spacecraft. The most effective way of representing the concave flow phenomena is to adopt a numerical procedure such as the Monte Carlo Test Particle method (Bird, 1976).

a cylinder is represented by an equation of the form x2 +y2 = r2 where Yis the radius of the cylinder. The representation for the Safyut complex is shown in Fig. 2. A control volume is constructed around this mathematical model of the configuration and then particles with velocities characteristic of the flow conditions are generated serially at the surfaces of the control volume. Their paths are followed to determine intersection with the vehicle body and upon re-emission the momentum transfer is calculated at each point of intersection. This is then summed for the whole body to calculate the overall aerodynamic coefficients. Four re-emission cases were considered in the test: diffuse, total accommodation (single reflection and multiple reflection} and specular re-emission (single reflection and multiple reflection). The re-emission model used in the Monte Carlo procedure is based on the widely used Nocilla model (Nocilla, 1963). This represents the incident and re-emitted particles as drifting Maxwellian fluxes with incident molecular speed ratio sir represented by :

velocity of spacecraft relative to the atmosphere random thermal velocity of atmospheric molecules ’ (2)

Re-entry aerodynamics derived from space debris trajectory analysis

643

FIG. 2. REPRESENTATION OF VEHICLEUSING TEST PARTICLEMETHOD.

an incident angle Oi of the incident flow vector relative to the spacecraft surface, a re-emission molecular speed ratio s,, and a temperature ratio TJT, relating the static temperature of the re-emitted flux to that of the incident flux. Therefore for the diffuse case considered above, S, = 9,

s, = 0, 0, = 90”,

T,/T, = 1

and for the specular case s, = 9,

S, = 9,

0, = Oi,

T,/T, = 1.

Due to the rotation of the Salyut complex about its spin axis and the precession of this axis at an angle of 20” about the normal to the velocity vector it was necessary to carry out a large number of runs to calculate the coefficients at different rotation and precession angles and then integrate these values to get the spin and precession averaged values. The aggregate values are given in Table 1. It can be seen that there is very little difference between the single reflection and multiple reflection values for diffuse re-emission. This is to be expected because the re-emitted molecules will have a very small proportion of their initial momentum available after TABLE 1. SPIN AND PRECESSION-AVERAGED VALUESOF DRAG Re-emission

Reflection

Drag coefficient x area (m’)

Diffuse Diffuse Specular Specular

Single Multiple Single Multiple

277.3 216.9 253.9 305.6

total accommodation at the surface. The differences between the single reflection and multiple reflection results for the specular case are, however, quite different. Again we might expect this as the re-emitted molecules in this case will retain their initial momentum even after numerous impacts with the body. Perhaps more importantly there is a marked difference between the results for diffuse re-emission and the specular cases. This suggests that if we can accurately derive the drag coefficient from orbital analysis, it should be possible to determine the nature of re-emission occurring at the surface of the Sulyut complex during orbital decay. 4. ORBITAL THEORY It can be shown (King-Hele, 1987) that the change in the size of a satellite orbit denoted by the semimajor axis a, due to the influence of aerodynamic forces is given by :

da dE=

(1 +ecosE)3/2 -~2~F6(,-ecos~)“7’

(3)

where E is the eccentric anomaly around the orbit, e is the eccentricity of the orbit, p is the atmospheric density, Fis an atmospheric rotation factor (0.92) and 6 is the satellite ballistic coefficient, given by :

where S is the normalizing area of the body (square metres), m is the mass of the body (kilograms) and C, is the drag coefficient of the body.

644

R.

Using a simple exponential sphere such that : P=

CROWTHER

mode1 of the atmo-

ppexp{(rp-r)/W,

(5)

where H is the density scale height (kilometres), the subscript p denotes the perigee value and r is the altitude given by : r = a(1 -ecosE),

(6)

Hand pp are derived from the CIRA 1972 atmosphere given a particular altitude and level of solar activity prevailing at the time of analysis. This influence is shown in Figs 3 and 4. The values of neutral density derived from the model can be expected to differ from the actual values by less than 10%. However this is still very important and must be accounted for in the determination process. This uncertainty in the derived value of density and the fact that we are dealing with a vehicle moving through an atmosphere at low altitude during a period of very high solar activity means the influence of the diurnal variation in atmospheric density can be neglected for the Sulyut 7 case (King-Hele, 1987). Integrating E between 0 and 27r, the change in semimajor axis around one revolution for a near circular orbit (in the case of the S&W station e = 0.0003) is given by (King-Hele, 1987) :

0’

1

100

150

I

I

200

250 altitude

I

I

300

350

_!

400

(km)

FIG. 4. VARIATION OF SCALE HEIGHT WITH F,O., FROM CIRA

12.

Aa = -2na2F~

CDp, exp (ae/H)[I,(ae/H) + 2e1,

(4H)l,

(7)

where I,, and I, are Bessel functions. Since orbital period T cc usi and the time to complete one orbital revolution At = T, then the rate of change of orbital period with time, T is given by :

so that the rate of change of period is given by : ? = -3~:

C,,Fap, exp(ae/H)[I,,(ae/H) +2eZ, (ne/H)l.

(9)

Given that all the parameters apart from SCn can be derived from the orbital data or are already known from models, we can use a least-squares fitting procedure to give the value of SC, which gives the best fit to the observed data for change in orbital period.

I

l.OE-13’ 100

150

200

250 altitude

FIG. 3.

VARIATION

OF

NEUTRAL

300

350

400

(km) DENSITY

CIR4 72.

WITH &,,

FROM

5.

APPLICATION

OF THE THEORY

The source of the observation data is USSPACECOMMAND via the NASA two line element sets determined for target objects. The orbital period for

Re-entry aerodynamics

derived

from space debris trajectory

2 50

Solar Activity ( 104 Jy) 2 00

5451

:

1 50

542: Orbital Period(s)

I

,

ri

I

I

1345350355360365370375360385 Days FIG.

5.

ORBITAL

from t/l/90

PERIOD VARIATION FOR S&u

7

COMPLEX.

the observation interval is plotted against time in Fig. 5. Also shown are the varying levels of solar activity during this period which would influence the atmospheric density value. The rate of change in orbital period was calculated and is plotted in Fig. 6. It is apparent that this is very sensitive to the level of solar activity and a slight lag can be observed between rises in solar activity and increases in the rate of change of

,250

Solar Activity (IO’ Jy)

-: period change (s/day) _’

IL

( 340 I 345 I

I

350

I

355 Days

FIG. 6.

I

380

,

,

385

from

370

,

375

,

380

385

111190

RATE OF CHANGE OF ORBITAL PERIOD WITH TIME.

645

orbital period. The least-squares fitting procedure was applied to the orbital data and a value of SCo of 275.8 m* was determined. This is very close to the value of SCu predicted for a diffuse re-emission derived from the Monte Carlo simulation and given in Table 1. The orbital determination technique is a very robust method as errors in semi-major axis measurement are only of the order of 40 m and this translates to a negligible error in the determined value of SCu. There is an uncertainty of the order of 1% in the values of SC,, derived from the Monte Carlo Test Particle calculations. This is due to statistical fluctuations inherent in the numerical approximation. It was stated in Section 3 that the value of neutral density derived from CIRA could differ from the actual value by the order of 10%. Therefore it is necessary to perturb the values of density used in the determination process in order to evaluate the influence on this uncertainty on the derived SC, value. When this is carried out, assuming a 10% underestimate of true density value and a 10% overestimate, it is found that the bounds of uncertainty relating to the derived value of SCo extend to 250.62 and 306.30 m2. Comparison with Table I then shows that the whole range of possible SCu values calculated for the Salyut 7 complex can be accommodated within this range. This might suggest that the results are inconclusive. It can be argued however that the somewhat artificial biasing of the true values of neutral density relative to the CIRA results, i.e. CIRA values = 90% true density values or CIRA values = 110% true density are not realistic. In practice we would expect the CIRA values to be evenly distributed about these true values of neutral density. What we can confidently conclude from this analysis is that our ability to use this technique to determine drag coefficients is limited by our knowledge of the atmosphere.

6.

orbital

analysis

CONCLUSIONS

The free molecular drag coefficient for the debris object Sulyut 7 has been derived using both numerical Monte Carlo techniques and orbit determination. Through comparison of the results of the two methods the nature of re-emission at the surface of the vehicle appears to be very close to the totally accommodated, diffuse case. For this natue of re-emission the influence of multiple reflection on the free molecular drag coefficient for this vehicle is negligible. However it should be noted that the results are sensitive to the assumed atmospheric density values and if we assume that the CIRA atmospheric model consistently underestimates or consistently overestimates the actual density by lo%, then the specular single reflection and

646

R. CROWTHER

multiple reflection models can also be accommodated within the results.

International Symposium on Spacecraft Flight Dynamics, ESOC, Germany, October 1991. ESA Spec. Publ. 326.

REFERENCES Bird, G. A. (1976) Molecular Gas Dynamics. Oxford University Press, Oxford. CIRA 72 (1972) COSPAR International Reference sphere. Academic Verlag, Berlin.

applied to ANS-1 (1975704). Planet. Space Sci. 39,129. King-Hele, D. G. (1987) Satellite Orbits in an Atmosphere. Blackie. Glaseow. Klinkrad,‘H. (1991) Review of ESOC re-entry prediction results for Salyut l/Kosmos 1686, in Proceedings of 3rd

Atmo-

Crowther, R. (1991) Debris re-entry prediction: Salyut 7IKosmos 1686, in Proceedings of ESA Workshop, ESOC, Germany, April 1991. ESA Spec. Publ. 345.

Crowther, R. and Stark, J. (1991) The determination of the gas-surface interaction from satellite orbit analysis as

Mehrholz, D. and Magura, K. (1991) Radar observations of Salyut 7/Kosmos 1686, in Proceedings of ESA Workshop, ESOC, Germanv, Auril 1991. ESA Snec. Publ. 345. Nocilla, S. (1963)Thesurface re-emission law in free molecular flow. Proceedings of the 3rd Symposium on Rarefied Gas Dynamics, pp. 327-346. Academic Press, New York. Sehnal, L. (1983) Determination of the basic constants of satellite atmospheric interaction from the analysis of the motion of 1974-704. Ado. Space Res. 3,91.