Determination of thermospheric quantities from simple ionospheric observations using numerical simulation

_lmmal of Amosphcric and

Tenewial

Physics,

Vol. 39,

pp. 531 to 537. Pergamon Press, 1977. Printed in Northern Ireland

Determiuation of thermospheric quantities ionospheric observations using numerical

from simple simulation

D. A. ANTONIADIS Radioscience Laboratory, Stanford University Stanford, CA 94305, U.S.A. (Receiued 2 February 1976; in reoised form 29 September 1976) Abstract-Measured ionospheric electron content and peak electron concentration data are introduced into a numerical simulation of the ionosphere to yield values of induced plasma drifts and

exospheric neutral temperatures consistent with the observations. Data collected on 23-24 March 1970 on the East Coast of the U.S.A. are analyzed and the results are in agreement with incoherent radar measurements at Millstone Hill, Massachusetts. Neutral winds and meridional exospheric temperature pradients that aive rise to the computed plasma drifts are calculated through the use of a dynamic model of the thermosphere.

of the present

1. INTRODUCTION

The electron content, Z,, and the peak electron concentration, N,,,,,, of the ionosphere are two important quantities observable by means of relatively simple techniques. For this reason and also because of their practical importance in communication systems, both quantities are under continuous observation from several locations on the Earth. The most common method for observing Z, of Faraday is the monitoring rotation of radiowaves transmitted from synchronous satellites, using polarimeters, while N,,,, is obtained from vertical incidence sounding of the ionosphere, using ionosondes. Both Z, and N,,, data have been used extensively in empirical studies of the ionosphere. They have also been used by several workers for testing the validity of theoretical models of upper atmospheric processes. However, only few attempts have been made to develop methods using Z, and/or N,,,,, to determine thermospheric quantities. One example of such methods is the work of TITHERIDGE (1973), in which he attempts to determine the exospheric neutral temperature, T, , from the ionospheric slab thickness, T, defined as T= ZJN,,, In the present paper we outline the development of a novel method that allows the determination of the vertical plasma drift velocity W, with good temporal resolution using measurements of Z,.‘Drift’ in this work is defined as the bulk motion of plasma induced by neutral winds and/or electric fields (c.f. RISHBETH, 1972). The method consists of a computer simulation of the ionosphere based on a theoretical dynamic model. By means of a numerical technique, developed during the course

‘forced’

work, the simulated ionosphere is to follow the observed Z,. This is accom-

plished by letting the vertical plasma drift become a dependent variable of the simulation. We also present an extension of the method just described, that allows the calculation of the exospheric neutral temperature, as well as the calculation of drift velocity, by using observed N,,,, in addition to Z,. In this case both exospheric temperature, T,, and W are dependent variables. Finally, we present a technique for calculating the vertical profiles of the neutral wind, by coupling the ionospheric model to the thermospheric dynamic model presented by ANTONIADIS (1976a). The methods just described have been applied to a geomagnetically quiet day with data collected at the East Coast of the United States. The results obtained are in agreement with incoherent scatter radar measurements at Millstone Hill, Massachusetts, during the same day. 2. THE MODEL

A block diagram of the basic mode (mode 1) of the ionospheric model used in this work, is shown in Fig. 1. Many of the details of this model are given by ANTONIADIS (1967b) and thus only a brief description will be given here. As can be seen in the figure, the ionospheric model is composed of a neutral atmospheric model, the positive ion continuity equations, the plasma heat balance equations, and a photoionization and solar extreme ultraviolet (EUV) model. These building blocks of the ionospheric model are interrelated as shown in the same figure. The various inputs and outputs, as

531

D. A. ANTONIADIS

532

“k%”

SPECTRUM

Fig. 1. Block diagram of the ionospheric model. Also shown are the input parameters required for the basic mode (mode 1) of simulation. This mode is used to calculate the time dependent, induced

vertical plasma drift, W, from electron content data. well as the blocks themselves, are briefly described in the following paragraphs. The neutral atmosphere is assumed to be in diffusive equilibrium and the vertical distributions of 0, O?, N2 and ‘Z’,,are given by the analytic expression of BATES (1959). The necessary inputs to this atmospheric model are the concentrations at the bottom of the thermosphere, at 120 km, and the exospheric neutral temperature, T,. The ionospheric constituents considered are, O’(4S), O’(‘D), NZ+, 02+. The concentration of O’(“S), N1, and the vertical ion drift velocity, W, are derived from the continuity equation of 0’(4S), in the form:

aN,

at=

Q,-PN,-F $+A$'+BN, 1

(1)

where Q1 is the net production rate of 0+(‘S) resulting from photoionization processes and calculated using the solar EUV spectrum of HINTEREGGER (1970) multiplied by an adjustable factor F,,,; /3 is the O’(4S) loss rate coefficient; D, the effective diffusion coefficient for vertical transport related to the ambipolar diffusion coefficient by D, = Da sin’ Z, where Z is the geomagnetic inclination or dip angle. The coefficients A and B are functions of altitude and of neutral and plasma temperatures. The exact expression for the vertical plasma drift velocity is: W(z) = cos Z sin D + U, sin’ Z (2)

where x, y, z, are our Cartesian coordinates correspondingly eastward, northward and upward, E the electric field, U the velocity of neutral wind and D the magnetic declination. Since wind velocity is generally height dependent, it can be seen from (2) that W(z) would also, under most circumstances, be altitude dependent. However, in order to keep the solution of (1) in terms of W and N1, as simple as possible, we have made the assumption of an altitude independent ‘effective’ W. It turns out, as will be seen later on, that the value of this effective W is very close to the value of an altitude dependent drift at the peak of the F2-layer. The boundary conditions, b.c., for (1) are: chemical equilibrium at the lower boundary (z, = 120 km) and the vertical plasma flux, &, at the upper (zu = 800 km). & is the flux coupling the ionosphere with the protonosphere and is determined as described later on. In addition, the vertical O’(4S) content, I, is used as a b.c. so that a simultaneous solution for the vertical profile of N, and the vertical plasma drift, W, is obtained. Details of the numerical technique are given in ANTONIADIS (1976b). In order to derive I,, from the vertical electron content, Z, (which is obtained from measurements), it is necessary to determine the contribution of the other ions, primarily 0,’ and NO’ and also N,+ and O’(*D) to I,. This is accomplished by solving the continuity equations of those ions to yield their concentration profiles. These equations however can be and are simplified substantially by assuming negligible transport, and, for N,’ and O’(‘D) further assuming chemical equilibrium. The electron and ion temperatures are obtained from the heat balance equations (HERMAN and

533

Determination of thermospheric quantities CHANDRA, 1969):

winter conditions), the sensitivity of W to T, is minimal and therefore if both W and T, are dependent variables in the model they can be calculated separately. The second fact is that N,,, is very sensitive to T, and thus the latter can become a dependent variable in the model if N,,, is externally furnished. The block diagram of mode 2 is illustrated in Fig. 2. Note the replacement of the externally supplied Tk, by N,,,. The function of the T, control block is to calculate T, at every time step of the numerical simulation by determining the correction AT, that must be applied to the T, of the previous step. Denoting time steps by the index j the following empirical equation is used:

(3)

where Tb, T, are the electron and ion temperatures; (2H the electron heat input given by SWARTZ et al. electron(1972); L,,, L,,, L, the electron-ion, neutral, ion-neutral heat loss rates; K,, K, the electron and ion thermal conductivities; k Boltzmann’s constant and N,, Ni the electron and ion concentrations. Equations (3) are by no means comptete since mechanisms such as convection and non-local heating are neglected. However, they are adequate for our purposes. Boundary conditions for (3) are: T, = T, = T, at the lower boundary and the electron temperature at the upper boundary, TeZ,S. The differential equations of the model are solved numerically as discussed in ANTONIADIE (1967b). The resulting computer program is very efficient; for example, a complete solution for a single time step with 70 altitude cells, requires approximately 2 set on a SDS E-5 computer. A second mode of the ionospheric model (mode 2) has also been developed by slightly modifying the basic mode just described. This second mode allows the determination of W and also of T,, as explained further on by using observed values of peak electron concentration, N,,,, as an additional boundary condition to the ionospheric model. Two facts have permitted the development of mode 2. The first is that under most circumstances studied so far (medium to high solar activity, equinox and

where r = I,j/N,,,,,i is often referred to as ‘slab thickness.’ The initial values of T, as well as of all other variables of the ionospheric simulation are obtained by initially assuming the ionosphere (usually around local noon) to be in steady state and solving the model equations until all variables converge to constant values.

3. RESULTS The ionospheric model, in the configuration displayed in Fig. 1 (mode l), was used to derive the behavior of vertical plasma drifts for a 24 hr period during days 23-24 March 1970. The columnar contentwas measured at Sagmore Hill, Massachusetts using the Faraday rotation technique with beacon signals of the geostationary satellite ATS-3. The

t

MODE2

Fig. 2. Block diagram of the ionospheric model (mode 2) used to calculate time dependent, induced vertical plasma drift, W, and exospheric neutrai temperature, T,, from electron content and N,,, data.

D. A. ANTONIADIS

534

420 km ionospheric point was at 39”N, 72”W. Vertical plasma drifts for the same period of time had been derived by SALAH and HOLT (1974), from radar data obtained at Millstone (42.6”N, 71S”W). This proximity provides an excellent opportunity for evaluation of our technique.

s

LON:72OW LAT = 39ON DAY=82-83170 I

LON : 72’W LA, = 39. N DAY :82-83/70

Fig. 4 Exospheric neutral temperature observed at Millstone Hill and calculated from electron content and N,,, by means of mode 2 shown in Fig. 2.

EST 201-I

MEASUREMENTS CALCULATIONS

”

--004 EST 0

500

-

MEASUREMENTS CALCULATIONS

Fig. 3. Day 23-24 March 1970. From top to bottom: Observed columnar electron content, I,; calculated and observed induced vertical plasma drift; calculated and observed N,,,,,; calculated and observed h,.,.

The time dependent values of the upper boundary electron temperature, TeC,R,plasma flux, &, and exospheric neutral temperature, T,, were taken from published incoherent radar measurements (EVANS, 1971a, b; SALAH and EVANS, 1973). A comparison of our results with the observations is presented in Fig. 3. Shown from top to bottom are: the observed (vertical) columnar electron content; the vertical ion drift deduced from our calculation together with the drift at 300 km from radar data (SALAH and HOLT, 1974); the calculated peak ion concentration, N,,,, together with the N,,,,, measured by the Wallops Island ionosonde (38ON, 75”W); the calculated peak height, h,,,, together with the observed h,,, at Millstone Hill. As can be seen there is good agreement between calculations and observations. In Fig. 4 we present the exospheric temperature T, obtained by using our ionospheric model in mode 2, indicated in Fig. 2, in which N,,,,, is an additional input. Both I, and N,,, were smoothed (mostly to remove TIDs from N,,, which are less apparent in I,) by a 2 hr running mean filter and the N,,, time basis was shifted by +12 min to make it coincide with the time basis of the ATS-3 ionospheric point. The vertical plasma drift, yielded by this mode is essentially the same as that shown in Fig. 3, and therefore is not repeated here. Up to this point we have derived vertical plasma drifts by means of modes 1 and 2 of the ionospheric simulation without specifically considering their origin. These drifts have been assumed altitude independent as described in Section 2 where we defined such drifts as ‘effective’, meaning that they have the same effect on the plasma content as altitude dependent drifts under the same circumstances. If the calculated drifts were entirely the result of electric fields then, they could be immediately related to the eastward component of the electric field (assuming D ==O”), as can be seen

Determination

of thermospheric

~_--------__--

535

quantities

_-____

----

“E”

SPECTRUM t

t

W N,n

Ni

L-------.-__---___-_--_____

J v

MODE

v

3

EQUATIONS OF NEUTRAL AIR MOTION u, k!, uy (21 +I---+

ah a+

Fig. 5. Block diagram of the combined ionospheric-thermospheric model (mode 3). This mode allows the calculation of vertical profiles of neutral winds as well as exospheric temperature and its gradients, from electron content and IV,,,, data. from equation (2). However, since at midlatitudes dent) W(z), which at h,,, had the same value as the drift in (a). All other model parameters were (under geomagnetically quiet conditions), vertical kept the same. The results of the two runs are plasma drifts are primarily driven by neutral winds, shown in Fig. 6. It can be seen that during most of it is reasonable to relate the calculated drifts to the day the calculated contents agree well. The these neutral winds. larger differences at night are due to the sharp Neutral winds and the meridional gradient of increase of vertical drift with decreasing altitude exospheric temperature can indeed be calculated below h,,,, associated with the low nighttime ion from the results of the ionospheric simulation drag at these heights. under mode 2 through a dynamic thermospheric Neutral winds and meridional gradients of exosmodel in a way very similar to that presented by pheric temperatures were calculated using mode 3 ANTONIADIS (1976a). In the present case, instead of using incoherent radar data, the data produced LON = 7Z’W by the ionospheric simulation are coupled directly LAT = 4E’N DAY = 82.83/70 to the thermospheric model. We refer to this mode -___ 5 of coupled ionospheric-thermospheric simulation - ALT. INDEP W lo 1 as ‘mode 3.’ A functional block diagram of this --:ALT. DEP W lb) : E 4mode is given in Fig. 5. ” As discussed in ANTONIADIS (1976a), the al- R titude dependent horizontal components of the wind are specified through the thermospheric model from the value of the vertical plasma drift at a given height. At present the assumption is made that the value of this drift is equal to the value of 4 the effective drift, W, and that the appropriate height is the peak of electron concentration, h,.,, calculated through the ionospheric model. oc-l---i_l_i__1 I I I 1 We have tested the above assumption by com00 04 Of3 12 I6 20 00 LT paring the plasma content I,, calculated through the model in two ways (note that I. in this case was not an input): (a) the model was used with an altitude Fig. 6. Electron content calculated with: (a) altitude indeindependent plasma drift W, and (b) the model was pendent vertical plasma drift, and (b) neutral wind induced (altitude dependent) drift. used with a wind compatible (i.e. altitude depen-

D. A. ANTONIADIS

536

zl TO d Y

2

--

ROBLE

(1974)

L-II I;!

16

04

20

08

12

ES

Fig. 7. Day 23-24 March 1970. Calculations by means of model shown in Fig. 6. Top panel: The wind induced vertical plasma drift. Contours annotated in m/s. Positive velocities are upward. Bottom panel: The meridional neutral temperature gradient calculated here and the same as calculated by ROBLEet al. (1974) directly from

the radar data.

of the model and the I, and N,,,,, data already discussed. The top panel of Fig. 7 illustrates the calculated wind-induced plasma drift and the bottom one, the calculated meridional temperature gradient together with the gradient derived by ROBLE et al. (1974) for the same period. The solar EUV spectrum multiplier, Feuy, used in all modes of simulation was equal to 2. This is consistent with the fact that the analyzed time period occurred during relatively high solar activity (c.f. ROBLE, 1976). The neutral composition at 120 km was kept constant with 0, N2 and 0, concentrations of 1.35, 3.0 and 0.5 (X 10” cma3), respectively. These concentrations are in general agreement with those used by ROBLE (1975). However, it was not found necessary to vary the 0 concentration between day and night, as he did. Under the assumption that F,,, and n(0)iZO are constant with respect to time, they are easily adjusted by comparing the calculated and observed h,,,. A general sensitivity analysis of the calculated plasma drifts with respect to T,, and &is is

difficult. However, a series of numerical experiments for geomagnetically quiet, mid-latitude, equinox and winter conditions has led to the following generalizations: (a) The sensitivity to T,,, is very small throughout the day and a crude model of its behavior is sufficient. For instance, in the calculations presented, setting T,, constant and equal to 2500 K, 3000 K or 3500 K affects our results by a negligible amount. (b) The sensitivity to & is also small provided that (PUBis significantly smaller than the integrated rate of production or loss of ionization. This indicates that under most circumstances, & can be safely set equal to zero or any reasonable value when the solar zenithal angle is smaller than about 95”. However, during nighttime a much closer examination of individual cases is necessary. In the particular period under investigation the rate of nighttime ionization decay is much larger than 2 x 10’ cm-’ set-‘. Since the average downward flux observed at Millstone Hill was less than 5 x 10’ cm-’ set- ’ the sensitivity to a & of that order is very small. Indeed, setting &,l3 = 0 for the complete diurnal run barely affects our results. 4. CONCLUSIONS We have presented a new method of time dependent ionospheric simulation, The basic concept underlying this method is the concept of ‘forcing’ parameters of the simulated ionosphere to follow the behavior of the corresponding observed ionospheric parameters. In the processes of doing so, quantities in the simulation, that normally must be furnished externally, become dependent variables and are thus determined so as to be consistent with the observed data. Three variations or modes of the ionospheric simulation have been identified in relation to respective input and output quantities. The basic mode (mode 1) allows the calculation of the ‘effective’ vertical plasma drift velocity compatible with observed electron content. The second mode (mode 2) allows the calculation of both ‘effective’ drift and exospheric neutral temperature compatible, respectively, with observed electron content and N,,,,,. The third mode (mode 3) is an extension of mode 2 that permits the calculation of vertical profiles of neutral wind velocity compatible with the calculated ‘effective’ vertical drift and exospheric temperature. For this latter mode a thermospheric neutral model, presented in a previous paper, is coupled to the ionospheric model. Results produced by the introduced ionospheric simulation method under all three modes are in satisfactory agreement with corresponding results

Determination of thermospheric quantities of incoherent scatter radar measurements. It can thus be concluded that both the employed ionospheric model and external parameters (such as atmospheric composition and solar EUV flux) are reasonably realistic representations of natural processes, at least under the conditions prevailing during the time period of the simulation in this work. Furthermore in view of the simplicity in acquisition of the required ionospheric input parameters, it is

531

believed that the simulation method presented here, will be a useful new source of neutral wind and exospheric temperature data. Acknowledgements-Many valuable discussions with A. V DA ROSA and 0. K. GARRIOTTof Stanford University are

greatly appreciated. The author wishes to thank J. KLOBUCHARof AFGL for providing the ATS-3 data. The present work was supported by NASA contract NCR OS020-001.

REFERENCES

ANTONIADISD. A. BATES D. R. EVANS J. V. EVANSJ. V. HERMANJ. R. and CHANDRAS. HINTEREGGERH. E. RISHBETHH. ROBLE R. G. ROBLE R. G. ROBLE R. G., EMERY B. A., SALAHJ. E. and HAYS P. B. SALAH J. E. and EVANS J. V.

1976a 1959 1971a 1971b 1969 1970 1972 1975 1976 1974

.I. atmos. ten. Phys. 38, 187. Proc. R. Sot. A253, 451. Radio Sci. 6, 609. Radio Sci. 6, 843. Planet. Space Sci. 17, 815. Annls. Ghophys. 26,547. .I. afmos. rev. Phys. 34, 1. Planet. Space Sci. 23, 1017. J. geophys. Res. 81, 265. J. geophys. Res. 79, 2868.

1973

Space Research XIII, p. 267. Akademie-Verlag,

SALAHJ. E. and HOLT J. M. SWART~W. E. and NISBETJ. S. TITHERIDGEJ. E.

1974 1972 1973

Berlin. Radio Sci. 9, 301. J. geophys. Res. 77, 6259. Planer. Space Sci. 21, 1775.

Reference is also made to rhe following unpublished maferial: ANTONIADISD. A.

1976b

Determination of thermospheric quantities from ionospheric radio observations using numerical simulation, Tech. Rept. No. 18, SU-SEL-76013, Stanford Electronics Lab., Stanford University.