T H E O R E T I C A L MODELS OF N O N T H E R M A L PROCESSES IN T H E A T M O S P H E R E S OF V E N U S , E A R T H A N D M A R S By Gregory G. Arkos B.Sc. Hons. (Geophysics), University of Manitoba, 1990. A T H E S I S S U B M I T T E D IN P A R T I A L F U L F I L L M E N T O F T H E R E Q U I R E M E N T S F O R T H E D E G R E E O F D O C T O R O F P H I L O S O P H Y in T H E F A C U L T Y O F G R A D U A T E S T U D I E S D E P A R T M E N T O F E A R T H A N D O C E A N S C I E N C E S ( G E O P H Y S I C S ) We accept this thesis as conforming to the required standard T H E U N I V E R S I T Y O F B R I T I S H C O L U M B I A August 1997 © Gregory G. Arkos, 1997 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of Earth and Ocean Sciences (Geophysics) The University of British Columbia 2219 Main Mal l Vancouver, B C , Canada V6T 1Z4 Date: August, 1997. Abstract The present state of planetary exospheres is determined largely by satellite and ground based observations which are predominantly measurements of the emissions of exospheric constituents. Such observations are responsible for the growing recognition of the impor-tance of nonthermal collisional processes in determining the distribution and escape of species in planetary exospheres. Nonthermal processes provide an important enhanced escape mechanism for lighter species such as hydrogen. They may also make escape pos-sible for heavier species, such as oxygen, nitrogen and carbon, for which thermal escape is very small. Nonthermal processes have been employed in order to understand discrepan-cies in the terrestrial helium budget. They have also been used to reconcile discrepancies between observed and calculated escape fluxes for hydrogen on Earth. Nonthermal pro-cesses have also been utilized to explain features observed in the exospheres of other planets. These include the measured deuterium-to-hydrogen ratio on Venus and the ex-tended hot oxygen corona on Mars. Given the importance of nonthermal processes it is clear that exospheric conditions are determined to a large extent by collisional processes and that the collisionless model must be reconsidered. The formation of translationally energetic (or hot) oxygen coronae via the nonthermal process of dissociative recombination of 0^ in the atmospheres of Venus and Mars is examined using both hydrodynamic and kinetic theory approaches. Of interest is the distribution of hot oxygen at altitude resulting from production and transport from lower altitudes. It is found that an extended hot oxygen corona can be predicted from either approach, although the magnitude and extent of the predicted coronae vary significantly. Product velocity distribution functions describing the rate of production of hot atoms for the atomic systems H - H + , D-H - f , O - H , and O - D are calculated for a variety of nonthermal processes, including direct-elastic and charge-exchange collisions. The cal-culations are carried out using a kinetic theory approach, and utilize direct numerical integration techniques. The calculations incorporate realistic, quantum-mechanical col-lision cross sections for each system so as to accurately describe the kinematics of the collision process. Energy exchange rate coefficients for each of the atomic systems are calculated and compared with results obtained using a more complicated Monte Carlo approach. The product velocity distribution functions are also used to estimate the es-caping fractions of H and D as a result of nonthermal direct elastic energization by hot oxygen atoms. These kinetic theory calculations are compared to work done by other workers using Monte Carlo methods incorporating approximate and quantum mechanical cross sections. The calculations show that the fraction of hot deuterium produced via direct energization by hot oxygen, while less than the fraction of hot hydrogen, is not negligible as previously believed. A n altitude dependent, kinetic theory approach is used to calculate the rate of escape of atmospheric constituents, in the context of escape resulting from energization of neutral atmospheric species via nonthermal processes. The reduction of the escape rate by the ambient atmosphere is included through an altitude dependent parameter describing the probability of escape, although the effect of thermalization via collisions with the background is neglected. Temperature and density profiles used in the calculations are taken from available atmospheric data and from atmospheric models, and escape fluxes of hydrogen and deuterium are estimated for Venus and Earth. 111 Table of Contents Abstract ii List of Tables vii List of Figures ix Acknowledgements xii 1 Introduction and Basic Theory 1 1.1 Introduction 1 1.2 A n Overview of the Exospheric Problem: Basic Theory 3 1.3 Nonthermal Processes 10 1.4 Kinetic Theory and the Boltzmann Equation 16 1.5 A n Overview of the Thesis 17 2 Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 22 2.1 Introduction 22 2.2 Constant Temperature Model 31 2.3 Non-constant Temperature Model 41 2.4 The Linear Boltzmann Equation Model 54 2.4.1 Discretization Procedure and Solution 58 2.5 Collisional and Diffusional Timescales 77 2.6 Summary 80 iv 3 Nonthermal Production of Energetic Hydrogen and Deuterium 83 3.1 Introduction 83 3.2 Collision Cross Sections 92 3.2.1 Quantum Mechanical Scattering 93 3.2.2 Calculation of the Phase Shifts 95 3.2.3 Cross Sections for H+-H and D+-H 96 3.2.4 Cross Sections for O-H and 0-D 101 3.3 Energy Exchange Rate Coefficients 107 3.3.1 Theory 109 3.3.2 Time Evolution of the Average Test Particle Energy 114 3.3.3 Results and Discussion 116 3.4 Product Velocity Distribution Functions 137. 3.4.1 Theory 137 3.4.2 Total Collision Rates 144 3.4.3 Rate of Production of Escaping Atoms 145 3.4.4 Results and Discussion 146 3.5 Summary 177 4 Kinetic Theory Calculations of Nonthermal Escape Fluxes 179 4.1 Theory 182 4.2 Summary 189 5 Summary and Suggestions for Future Research 195 A Density Profiles in a Collisionless Exosphere 199 B Details from Chapter 2 206 B . l Modified Fourier's Law 206 v B.2 Symmetrization of the Collision Kernel 208 B.3 Symmetrization of the Discretized Collision Operator 210 B.4 Solution of the Eigenvalue Problem 212 C The Differential Col l i s ion Cross Section 215 D Scatter ing Theory 219 D . l Quantum Mechanical Scattering 219 D.2 Calculation of the Phase Shifts 223 D. 3 Semi-Classical (WKB) Phase Shift Approximation 228 E Detai ls of the P roduc t ion of Ho t Atoms 230 E. l Derivation of £ 2 - # 230 E.2 Derivation of g + g 231 E.3 Transformation from d£2 to d3d£ 231 E.4 Transformation from x to the scattering angle 9 232 E.5 Derivation of the hard sphere product velocity distribution function . . . 233 E.6 Interpolation scheme for the calculation of integrals over reduced energy . 238 F Der iva t ion of the H a r d Sphere Col l i s ion Frequency 240 Bib l iography 244 vi List of Tables 1.1 Planetary escape velocities and energies 8 1.2 Planetary exospheric values for hydrogen escape 8 2.3 Diffusional versus collisional timescales for Mars 78 2.4 Diffusional versus collisional timescales for Venus 79 3.5 Energy exchange rate coefficients for H + - H 117 3.6 Energy exchange rate coefficients for H + - H from Hodges and Breig . . . 117 3.7 Energy exchange rate coefficients for H + - D 124 3.8 Energy exchange rate coefficients for 0 -H and 0-D 129 3.9 Timescales for the time evolution of the average energy for H 133 3.10 Timescales for the time evolution of the average energy for H + 134 3.11 Comparison of timescales for the time evolution of the average energy . . 136 3.12 Energy exchange rate coefficients for H + - H 159 3.13 Energy exchange rate coefficients for H + - D . . 160 3.14 Energy exchange rate coefficients for O-H and O-D 161 3.15 Comparison of energy exchange rate coefficients for hard spheres 161 3.16 Rate of production of escaping H on Venus 164 3.17 Rate of production of escaping H on Earth 165 3.18 Rate of production of escaping H on Mars 166 3.19 Rate of production of escaping H on Earth from Hodges & Breig . . . . 169 3.20 Rate of production of escaping D on Venus 170 3.21 Rate of production of escaping D on Earth 170 vii 3.22 Rate of production of escaping D on Mars 171 3.23 Rate of production of escaping H on Venus from 0 - H and O D interactions 172 3.24 Rate of production of escaping H on Earth from 0 - H and O D interactions 172 3.25 Rate of production of escaping H on Mars from 0 - H and O D interactions 173 3.26 Escape fraction for hot H from 0 - H elastic collisions on Venus 176 3.27 Escape fraction for hot H and hot D from elastic collisions with 0 . . . . 176 A.28 Exospheric particle classes 203 viii List of Figures 1.1 Stratification of the terrestrial atmosphere 4 2.2 Observed Martian density profiles 27 2.3 Observed Venusian density profiles 28 2.4 Stratification of the atmosphere by collisional property 29 2.5 Martian hot oxygen density (diffusion equation) 35 2.6 Venusian hot oxygen density (diffusion equation) 36 2.7 Martian hot oxygen density temperature dependence (diffusion equation) 38 2.8 Venusian hot oxygen density temperature dependence (diffusion equation) 39 2.9 Hot oxygen densities at large altitudes for Venus and Mars 40 2.10 Solutions of the collisionless, sourceless momentum conservation equation 48 2.11 Time dependence of hot oxygen density and temperature on Mars . . . . 63 2.12 Time dependence of hot oxygen density and temperature on Venus . . . 64 2.13 Time evolution of the Martian exospheric distribution function 66 2.14 Time evolution of the Venusian exospheric distribution function 67 2.15 Energy distribution function of hot oxygen on Mars 68 2.16 Energy distribution function of hot oxygen on Venus 69 2.17 Maxwellian fit of the energy density distribution 71 2.18 Maxwellian fit of the Venusian energy density distribution 72 2.19 Hot oxygen density profile for Mars derived from ballistic component . . 73 2.20 Hot oxygen density profile for Venus derived from ballistic component . . 74 3.21 Schematic of the charge-exchange and LTA 91 ix 3.22 Interaction potentials for H + - H 97 3.23 Differential cross sections for H + - H 98 3.24 Differential cross sections for D + - H 99 3.25 Total elastic collision cross sections 100 3.26 Interaction potentials for O-H 102 3.27 Total elastic collision cross sections 104 3.28 Total elastic collision cross sections 105 3.29 Differential cross sections for O-H/O-D 106 3.30 Fits to the LTA momentum transfer cross section 120 3.31 Temperature dependence of k% for the LTA 121 3.32 Dependence of energy exchange rate coefficients on bath temperature for H+-H 122 3.33 Dependence of energy exchange rate coefficients on bath temperature for D+-H 126 3.34 Details of the integrand in calculation of the energy exchange rate coefficient 127 3.35 Convergence of the energy exchange rate coefficient 128 3.36 Time evolution of the average energy in a relaxing H + - H system 131 3.37 Time evolution of the average energy in a heating H + - H system 132 3.38 Comparison of the time evolution of the average proton energy 135 3.39 Schematic of an elastic binary collision 138 3.40 Coordinate system for the calculation of collisional production 141 3.41 Product velocity distributions for the D P E cross section for H + - H . . . . 148 3.42 Product velocity distributions for the C E cross section for H + - H 149 3.43 Product velocity distributions for the LTA cross section for H + - H . . . . 150 3.44 Product velocity distributions for H + - H 152 3.45 Product velocity distributions for the DIR cross section for H + - D . . . . 154 3.46 Product velocity distributions for the C E cross section for D + - H 155 3.47 Product velocity distributions for the DIR cross section for O - H 156 3.48 Product velocity distributions for the D I R cross section for O - D 157 3.49 Dependence of charge exchange rate coefficients on bath temperature for D+-H 163 3.50 Product velocity distribution functions for O - H and O-D escape producing collisions 175 4.51 Geometry for a plane-parallel atmosphere 190 4.52 Temperature and density profiles used as input for flux calculations on Venusl91 4.53 Density profiles used as input for flux calculations on Earth 192 4.54 Density profiles used as input for flux calculations on Earth 193 4.55 Temperature profile used as input for flux calculations on Earth 194 A.56 Graphical illustration of exospheric particle classes 202 C . 57 Typical scattering configuration 216 D. 58 Wave function scattering by a central potential 222 E . 59 Geometry of the dynamics of an elastic collision 233 xi Acknowledgments I would like to acknowledge the support, patient advice, and many helpful suggestions of Dr. B . D. Shizgal throughout the course of my graduate studies. His input and guidance have contributed greatly to this work. I would like to acknowledge the members of my supervisory committee, Dr. D. Oldenburg, Dr. G. Fahlman, and Dr. R. F. Snider, for their constructive comments and advice during the writing of this work. I also thank the many graduate students and post-doctoral researchers who took time from their own work to discuss various aspects of this work, to help out with the specifics of computer code, or just commiserate when things weren't working. K i , Andrew, Duncan, Heli, Francois, and Gianpierro - your support was much appreciated. I would also like to thank my sister, Julia, and my entire family for their unconditional support and encouragement throughout my graduate work. Their unwavering support and love helped me to focus on completion of this work through some difficult and tumultuous times. Without them, I am certain that I would not have been able to successfully complete this work. I wish to acknowledge the Chemistry Department and the computing resources which were made available to me throughout the course of this work. Finally, I would like to gratefully acknowledge the sources of financial support through-out my graduate studies. To the Natural Sciences and Engineering Research Council, the Canadian Society of Exploration Geophysicists Trust Fund, Canadian National Scholar-ship Program, the Department of Geophysics and Astronomy, and Dr. B . D. Shizgal I extend my many thanks for making my graduate studies possible. xn Chapter 1 Introduction and Basic Theory 1.1 Introduction The understanding of the composition, dynamics, and evolution of planetary atmo-spheres is important in many fields, including geophysics, astronomy, and atmospheric science. A critical part in the determination of the overall behaviour of planetary atmo-spheres is the behaviour of the exosphere. The exosphere is the high altitude region of a planetary atmosphere characterized by vanishingly small gas densities.1 As densities are so low, the mean free path of particles in the exosphere are exceedingly long, and gas kinetic energies may be sufficient to permit atoms to escape the planetary gravitational potential. The altitude where the mean free path is equal to the scale height, defined later, is the bottom of the exosphere, known as the exobase. Below the exosphere, the atmosphere is predominantly turbulent, well mixed, and dense enough to be collision dominated. In this region of the atmosphere, hydrodynamic or fluid theory is often used to characterize gas behaviour. Above the exobase, in the exosphere, densities are low, collisions are rare, and gases are distributed primarily according to diffusive equilibrium rather than turbulent mixing. Under such conditions, the kinetic theory of gases is more appropriate to describe the behaviour and character of exospheric species. In the real atmosphere, the transition from one regime to the other is gradual and continuous. However, early exospheric models considered the exobase as a discontin-uous division between collisional and collisionless regions of the atmosphere. Despite 1 Chapter 1. Introduction and Basic Theory 2 this unrealistic treatment of the atmosphere, the discontinuous exobase model has been extensively employed in the interpretation of observations of Lyman-o: and Balmer-a emissions of atomic hydrogen in the exosphere.2-8 Atoms in a range of altitudes above and below the exobase can attain speeds in excess of the planetary escape speed and escape from the planetary atmosphere. The escape process may thus be viewed as an 'evaporative' process where high energy particles are preferentially removed from the atmosphere. This thermal picture of escape was formulated by Jeans,9 and the escape flux resulting from this evaporative escape process is known as Jeans flux. Comparison of escape fluxes inferred from satellite measurements to those predicted by Jeans based on neutral exospheric temperatures clearly indicated that the thermal escape process alone was insufficient. It was suggested by Cole 1 0 that excitation of 'cold' atmospheric species could take place via nonthermal processes such as charge exchange. Such translationally energetic (or 'hot') nonthermal product atoms would then have sufficient energy to escape from the planetary atmosphere. Such hot atom populations have been observed by satellite and ground based techniques in the Earth's exosphere, 1 1 ' 1 2 and inferred from measurements made of Venus 1 3 , 4 and Mars . 6 ' 1 4 The interest in hot atoms in the atmospheres of Venus and Mars is motivated by sev-eral factors. A thorough understanding of the hot atom populations is necessary to un-derstand observations of emissions in planetary exospheres made by instruments aboard spacecraft and from the ground. The long term evolution of planetary atmospheres, including the loss of hydrogen and oxygen related to the removal of water, requires an understanding of the transport and distribution of atomic species in the upper atmo-spheres of the planets. Accurate calculation of escape fluxes requires an understanding of the magnitude and extent of planetary hot atom populations. Such a description is also important in understanding and describing the interaction of planetary atmospheres Chapter 1. Introduction and Basic Theory 3 with the solar wind. For example, the effect of photoionization of energetic atoms in plan-etary exospheres and the subsequent interaction with the solar wind (pick-up, knock-on, and precipitation of photo-ions). There is also the 'direct' interaction of the solar wind 'striking' the upper atmospheres of Venus and Mars, since these planets do not have (very strong) magnetic fields to deflect or slow the solar wind. A good overview of basic exospheric physics which emphasizes nonthermal escape was given in a recent paper by Shizgal and Arkos. 1 5 The review by Chamberlain 1 6 provides a very detailed theoretical description of the collisionless exosphere. Additional reviews dealing with the problem of the description of exospheric processes have been written by Hunten and Donahue,1 7 Tinsley, 1 8 Fahr and Shizgal1 and more recently by Hunten, 1 9 ' 2 0 and Mahajan and Kar . 2 1 1.2 An Overview of the Exospheric Problem: Basic Theory The basic concepts of planetary atmospheric escape were given by Jeans9 and ex-tended to the conventional collisionless thermal escape model by Chamberlain. 1 6 A pic-ture of the stratification of the atmosphere is illustrated in Figure 1.1. The lower region of the atmosphere where turbulent mixing of gases leads to a homo-geneous composition is called the homosphere. Above this region, turbulence ceases and the vertical distribution of individual atmospheric gases is determined by their respec-tive masses. In this region, called the heterosphere, the density profile n ?(r) of the ith constituent is determined by the balance of gravity and gas pressures. When these two opposing forces are equal we have the condition of hydrostatic equilibrium, given by 2 2 dpi dr = nirriig (1.2.1) Chapter 1. Introduction and Basic Theory 4 Figure 1.1: Atmospheric divisions as a function of altitude for the terrestrial atmosphere. Indicated are the various divisions of the atmosphere, including the homosphere, hetero-sphere, and exosphere. The curves represent average temperature profiles for the bulk neutral atmosphere for quiet and active solar periods. From Banks & Kockarts . 2 2 Chapter 1. Introduction and Basic Theory 5 where g = GM/r2 is the gravitational acceleration, m; and T8 are the mass and (con-stant) temperature of the ith constituent, k is the Boltzmann constant, G is the universal gravitational constant, M is the planetary mass, and pi = ra,-fcT,- is the (partial) hydro-static pressure as given by the ideal gas law. If we substitute for g in equation (1.2.1) we have .dm GM kTi— = -ra,-m;—r-or rz Integrating both sides yields rni drii _ GMrrii r dr •j Tii0 Tl{ kTi Jr0 In (n t-/n 0) = ^^ 7*' (l/r ~ l/ro) m{r) = n0exp[r/Hi(r)-r0/Hi(r0)] (1.2.2) where r0 is some reference level, n0 is the density of species i at r0, and where we have defined the scale height for species i as JcT-r2 hT *w = ^ = ^ ) (1-2-3) Equation (1.2.2) is the barometric or hydrostatic density distribution. It is important to note that while we have included the correct radial dependence of the gravitational acceleration in deriving equation (1.2.2) the result is nonphysical in that the density n,-(r) is finite as r —> oo. This is a result of the fact that our assumption of hydrostatic equi-librium is no longer valid. As the density becomes vanishingly small, collisions become extremely rare and hydrodynamics breaks down. If we write z = r — r0, we have -GMrrii m(r) — n0exp = n0exp kTi GMrrii i kTir0 (l/(r0 + z)-l/r0) (1/(1 + * / r 0 ) - 1) Chapter 1. Introduction and Basic Theory 6 Assuming that z/r„ < 1, we have that g « g(r0), constant, we may expand the first term of the exponential, yielding GMrrn rii(r) — n 0exp = n0exp kTirl z (1.2.4) The barometric density distribution given by equation 1.2.4 now vanishes as r —> oo, as we would physically expect, although we no longer correctly account for the radial dependence of the gravitational acceleration. The exobase or critical level, r c , is defined as the altitude for which the mean free path, the average distance between particle collisions, given by 1 = - ^ - (1.2.5) is equal to the barometric scale height of the heaviest constituent, that is, £(rc) — H(rc). In equation (1.2.5), a is the energy independent total elastic collision cross section for the diffusing species and the background and N = i s the total density. i In the standard Chamberlain model, 9 ' 1 6 the atmosphere is considered collisionless above the exobase and collision dominated below it. Particles that reach the exobase from below and move upwards with speeds in excess of the escape speed, given by 2GM c e s c = y — (1.2.6) will escape from the gravitational field of the planet. This model assumes that above the exobase the particles move on collision free trajectories determined by the plane-tary gravitational field. The classification of exospheric species into classes of particles such as ballistic, satellite, and escaping based on this model is discussed at length by Chamberlain, 1 6 and Fahr and Shizgal,1 and is summarized in Appendix A . In these Chapter 1. Introduction and Basic Theory 7 models a Maxwellian distribution of particle velocities mc2 2*rj ( 1 ' 2 ' 7 ) is assumed to exist at the exobase. The escape flux of particles moving radially outward with speeds in excess of the escape velocity (at the critical level) is determined by aver-aging over the outward directed velocity for speeds greater than the escape speed, that is e=o /•oo f6=0 Fj = 2n / fmax cos 0 c3 d(cos 6) dc (1.2.8) so that the thermal or Jeans escape flux is given by 9 Fj = y J ^ ( l + Ac)e'A< (1.2.9) I V rmr In equation (1.2.9), Tc and nc are the temperature and density of the escaping species at the exobase respectively, and the escape parameter, A c , is defined by A c = | ^ (1.2.10) kTc v ; where the total energy for escape is Eesc = \mc2esc. The escape speeds and escape energies for hydrogen, deuterium, and oxygen from the terrestrial planets are shown in Table 1.1. Mars, the least massive of the three planets, has the lowest escape energy. The escape flux is determined by the ratio of the escape energy relative to the thermal energy, that is, by A c in equation (1.2.10). Table 1.2 compares the values of A c for the terrestrial planets; the very large value for Venus, arising from the low exospheric temperature, is evident. For Venus, thermal escape is insignificant at present, and nonthermal processes play a dominant role. However, Donahue2 3 has discussed the escape of primordial atmospheres from planetesimals for which A c is small and Jeans flux is significant. It is important to note that in equation (1.2.8) the three dimensional velocity in-tegral is carried out only over the upper half of the velocity space (6 > 0). Since the Chapter 1. Introduction and Basic Theory 8 Planet cf° (km/sec) EL (eV) (eV) EZc (eV) Earth 10.8 0.61 1.22 9.69 Venus 10.2 0.54 1.08 8.60 Mars 4.8 0.12 0 .25 1.91 Table 1.1: Planetary escape velocities and energies. The values given are for an alti-tude corresponding to the exobase level for each planet (see Table 1.2). Neutral atmo-spheric temperatures are taken from Niemann et a l . , 2 4 Banks and Kockarts 2 2 and Fox and Dalgarno 2 5 for Venus, Earth and Mars, respectively. Planet rc (km) Tc {°K) A c Fj/rtc (cm/sec) Earth 500 1000 7.06 7.94xl0 2 Venus 200 275 22.89 1.65xl0- 4 Mars 250 300 4.65 3.39xl0 2 Table 1.2: Planetary exospheric values for hydrogen escape. The exobase or critical level is at the altitude given by rc, the exobase temperature is given by T c , and the escape parameter at the exobase by A c . The Jeans' flux Fj is given with respect to the exobase density n c of hydrogen. Chapter 1. Introduction and Basic Theory 9 Maxwellian distribution function is isotropic, the flux would be zero if the lower hemi-sphere of the velocity space were not artificially removed. This model of atmospheric escape is an oversimplification. The actual distribution function at the exobase is not Maxwellian, because the escape process preferentially removes energetic atoms and non-thermal processes introduce them. The distribution that leads to an outward flow must be anisotropic, that is /(r,c) = r a a ;(r,c) + /™(r , c ,c?) (1.2.11) where 0 is the angle between the velocity c and the outward radial direction. Since the average of the Maxwellian distribution over all velocities is zero, the drift velocity is given by u(r) = | / ™ ( r , c ) c d c (1.2.12) It is reasonable to expect that escaping particles originate not only from the exobase but from a range of radial positions in the vicinity of the exobase. Hence, the differential escape flux found using equation (1.2.12), F(r)dr = n(r)u(r), depends on the radial position. The net escape flux is determined by the integral of F(r)dr over a range of radial distances in the vicinity of the exobase.26 This picture of atmospheric escape, as detailed by Lindenfeld and Shizgal 2 6 and Shizgal and Blackmore, 2 7 which considers both thermal and nonthermal escape as collisionally induced phenomena is discussed in Chapter 4. For rarified regions, such as in the exosphere, the behaviour of gases is best described by the kinetic theory of gases. Of principle interest to most workers in exospheric physics is the description of the bulk physical properties of the exospheric constituents, including the density, temperature, and heat conductivity. An important problem in this calcu-lation of the bulk physical properties is the proper treatment of the transition between Chapter 1. Introduction and Basic Theory 10 the collision dominated lower atmosphere and the nearly collisionless upper exosphere. A quantity which categorizes kinetic theory problems in the exosphere is the Knudsen number, Kn, defined as the ratio of the mean free path of the gas and the local scale height, For small values of the Knudsen number, Kn 0 ( 3 P) + 0 ( 3 P) + 6.98eV(0.22) 0 ( 3 P) + C f D ) + 5.02eV(0.55) 0 ( 3 P ) + CfS) + 2.79eV C f D ) + 0( 1D) + 3.05eV => C f D ) + CfS ) + 0.83eV It is assumed that the available energy is split between the two product atoms equally. From Table 1.1, the 0 escape energy for oxygen at the Venusian and Martian exobases is approximately 8.6 and 1.9 eV, respectively. Thus, only the first two dissociation channels provide enough energy (3.99 eV and 2.51 eV per oxygen atom, respectively) to permit Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 26 direct escape from Mars, and none of the above channels provide sufficient energy for direct escape of oxygen from Venus. The branching ratios for the O f dissociative recombination process are still only poorly known. 9 2 - 9 4 That is, the relative populations of the various product states of atomic oxygen, in addition to translational energy available to each product atom, are poorly constrained. An estimate of the branching ratios of two channels which give rise to product atoms with sufficient energy to escape from Mars is given in brackets next to the appropriate reaction. 9 0 We assume that all dissociation takes place along the (dominant) second channel only, with a dissociative recombination coefficient ccdr = 3 .0x l0 - 7 cm 3 sec - 1 . 6 ' 8 The source of hot oxygen is then given by S(z) = aDR[Ot][e-} = « D H [ 0 + ] 2 (2.1.2) where the square brackets [A] represent the density of species A. Following prior con-vention, we have assumed that the electron density is equal to the O f density.6 Altitude profiles of the O f density and the background gas densities are taken from available atmospheric profiles collected by planetary missions to M a r s 9 5 ' 2 5 ' 6 and Venus, 9 6 ' 2 4 and are illustrated in Figures 2.2 and 2.3 respectively. The densities are extrapolated beyond the tabulated values assuming that they follow a barometric distribution, as given by equation (1.2.4). A representation of the process of hot oxygen atom creation and subsequent diffusion giving rise to exospheric coronae is given by Figure 2.4. Production of hot oxygen atoms occurs primarily in the lower atmosphere where the density of O f is sufficiently high that dissociative recombination is frequent. As can be seen in Figures 2.2 and 2.3, the density of O f peaks peaks near 140 km altitude and decays rapidly above and below the peak altitude. The translationally energetic product oxygen atoms diffuse upward from their Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 27 Figure 2.2: Observed Martian density profiles. The data are from Viking missions of 1977. The solid curve is the 0^ density profile, the dotted curve is the background oxygen density profile, and the dashed curve the background carbon-dioxide density profile. From Fox and Dalgarno. 2 5 Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 28 Figure 2.3: Observed Venusian density profiles. The data are from observations made by the Pioneer Venus orbiter and probe. The solid curve is the 0^ " density profile, the dotted curve is the background oxygen density profile, and the dashed curve the background carbon-dioxide density profile. From Nagy and Cravens. 7 Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 29 Figure 2.4: Stratification of the atmosphere by collisional property. In the collisionally dominated and transition regions, collisions drive energetic atoms toward equilibrium with the thermal background. Above the exobase, in the 'collisionless' region, particles move in trajectories solely under the affect of the planetary gravitational field. Hot atoms in this region have extremely long mean free paths and so may form an extremely extended corona of energetic atoms. A typical hot particle's path is shown on the right hand side of the figure. altitude of creation until they achieve diffusive equilibrium. However, as the background density is still substantial at lower altitudes, the hot oxygen atoms have a high probability to suffer thermalizing collisions. Above the exobase, densities are sufficiently low that collisions are practically negligible, and atoms are free to move in accordance with their translational energy. The formation of an extended hot oxygen corona is thus a problem of the description of production and transport of hot atoms at low altitudes to above the exobase. This chapter examines of the creation of hot oxygen 'coronae' in the exospheres of Venus and Mars using two new approaches to the production and diffusion of hot oxygen in planetary exospheres. The first model for the formation of the hot oxygen coronae is Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 30 given by solution of the diffusion equation, which describes diffusion of the translationally energetic exospheric oxygen atoms through a cold background gas. The general diffusion equation for a minor species of density n diffusing through a background is given b y 2 2 dn d ~di = ~d~z _.. . f dn n dz H(z)\\ where D(z) is the diffusion coefficient, given b y 2 2 + S(z) (2.1.3) = (2.1.4) and S(z), the source of the hot oxygen atoms due to dissociative recombination of Oj", is given previously by equation (2.1.2). It is assumed that the collisional interaction between the hot oxygen and background gases is described by a hard sphere model. The collision frequency and average relative speed are thus given, respectively, by v(z) = nback{z)ohsVave Vave = ( ^ L + ***S (2.1.5) V \™>hot rnbackJ where riback is the thermal oxygen altitude density profile, Th0t and Tback a r e the hot oxy-gen and thermal oxygen temperature profiles, and ahs is the hard sphere collision cross section. The collision frequency and average velocity in equation (2.1.5) have assumed that the colliding species' distribution functions are Maxwellian. The diffusion coefficient is thus determined by the altitude dependent temperature of the hot oxygen atoms and density profile of the background thermal oxygen. For the first model examined in this chapter, we assume an altitude independent (constant) temperature for the hot oxygen atoms. This is equivalent to assuming that there is no thermalization or energy exchange between the diffusing hot oxygen and the background gas. The simple diffusion model is then extended by incorporating an altitude dependent temperature, which is calculated Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 31 self-consistently using coupled density and energy equations based on standard hydro-dynamic theory. 9 7 , 9 8 The incorporation of an energy equation accounts for the effect of thermalization or energy transfer between the hot oxygen and the background. A model of hot atom densities based on solutions of the local Boltzmann equation is also examined. The velocity distribution function of the hot oxygen atoms are calcu-lated at the altitude of the exobase and hot oxygen density profiles are then generated using Liouville's theorem. 1 6 The Boltzmann collision operator rigorously accounts for the thermalizing effect of (elastic) collisions, and so temperature and density profiles are cal-culated in a self-consistent fashion. As for the other models in this chapter, we assume the use of 'hard sphere' or energy independent collision cross sections. The models introduced in this chapter are relatively simple descriptions of the pro-duction and transport of hot atoms. Our goal is thus to investigate whether such models are sufficient to qualitatively describe the altitude distribution of densities of hot exo-spheric oxygen and to compare with similar studies using more complicated models by other w o r k e r s . 1 1 ' 6 9 , 7 9 ' 7 , 6 ' 3 3 ' 3 2 , 8 There are many uncertainties to consider in contrasting and comparing results between the various models. Poor constraints on the values of the dissociative recombination branching ratios, the dissociative recombination coefficient, neutral temperature and density altitude profiles, and the hot oxygen-background col-lision cross section contribute greatly to discrepancies between the various approaches. They also degrade the accuracy of the predictions of hot atom densities. 2.2 Constant Temperature Model As illustrated in Figure 2.4 the generation of an enhanced population of hot atoms at high altitudes in the exosphere requires transporting hot atoms from their altitude of production. It was of interest to determine whether simple diffusion of energetic atoms Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 32 through a stationary thermal background was capable of producing qualitative agreement with densities inferred from satellite and ground based emission studies and calculated by other workers using Monte Carlo and other techniques. If it is assumed that the density of the minor species has reached a steady-state value, the time derivative on the left-hand side of equation (2.1.3) may be set to zero to yield the steady state diffusion equation d_ dz x \ dn n + S{z) = 0 (2.2.1) dz ' H{z)\\ where D(z) is the diffusion coefficient given previously by equation (2.1.4). Equation (2.2.1) is then integrated in altitude from z0 to z to yield = -W)LS(z)dz'+c- (2-2-2) where C\ is the constant of integration. This equation is a linear first-order differential equation, and may be solved by use of an integrating factor," exp(/i(z)), where The integrating factor is applied to equation (2.2.2), and after some rearrangement and collecting of constants, we find n{z) = e~h^ f eh{2,) {P(z') + C x } dz' + C2e~h^ (2.2.3) where for clarity we have defined P W s -W)Cs(z')dz' (2'2-4) The two constants in equation (2.2.3) may be evaluated easily if boundary conditions are applied at the lower (z = z0) boundary. From equations (2.2.2) and (2.2.3), with z = z0, we thus find dn(z0) n(z0) 1 dz + H(z0) C2 = n(z0)eh^ (2.2.5) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 33 For the case of a constant scale height, H(z) = H(z0), and no source, S(z) = 0, equation (2.2.3) reduces to the barometric solution given by equation (1.2.4), •(z - z0) n(z) = n(z 0)exp (2.2.6) L H{z0) j Similarly, for a constant temperature but radially dependent gravitational acceleration, with the scale height given by kT H(z) = mg(z) we find after some algebraic manipulation that the argument of the integrating factor is given by h(z) = R + z0 R4-H{z0) H(z) where R is the planetary radius. Substituting this result in equation (2.2.3) we arrive at '-(z'-z'0y n(z') = n(2^)exp H(z'0) (2.2.7) which is identical in form to the barometric density equation except it is a function of a 'reduced altitude' (or 'geopotential height') z'-z' (* - *o)l 1 + z - zn R + z0, rather than altitude. To determine the hot oxygen density using equation (2.2.3) it is necessary to provide a temperature profile for the hot oxygen in order to calculate the altitude dependent scale height. For this model, some (arbitrary) constant temperature T0 is chosen. There is no a priori justification for a particular choice of this (constant) temperature, other than the physical argument that it should be between the background thermal oxygen temperature and the temperature equivalent to the initial hot oxygen dissociation energy. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 34 The solution of equation (2.2.3) involves a double integral. The integrations are done using a Simpson's Rule with 150 altitude steps. In order to accurately calculate the integrated hot atom source density, -P(z), we do a 250 point Simpson's Rule integration between each pair of altitude points from the Simpson's Rule integration over z' in equa-tion (2.2.3), and simply add the contributions from each to yield -P(zt')- This is necessary as the altitudes of interest span several thousand kilometers and so a straightforward in-tegration of the source over the entire span of z' requires excessive numbers of integration points (and time). The initial altitude is (arbitrarily) chosen below the exobase level, and is 90 km for Mars and 100 km for Venus. As a check of the accuracy of the numerical solution of equation (2.2.3), as detailed above, it was compared to the analytic solutions given by equations (2.2.6) and (2.2.7) for the case where the source is set to zero, and the scale height is constant or dependent only on the gravitational acceleration, respectively. Once the accuracy of the numerical solution of equation (2.2.3) was verified for the above test cases, the effect of an altitude dependent scale height (via equation(1.2.3)) and/or diffusion coefficient on the density profile of the hot oxygen was examined. A value of 8500 K was chosen as the (constant) temperature for the hot atoms.8 For the background thermal oxygen on Mars, a (constant) temperature of 180 K was inferred from the thermal oxygen density profile of Figure 2.2, given the assumption of baromet-ric behaviour. A temperature of 280 K was derived similarly for the thermal oxygen temperature on Venus. A hard sphere (total) collision cross section of 2 . 0 x l 0 - 1 5 cm 2 is used. The solutions for the hot oxygen density profiles for Mars and Venus are shown in Figures 2.5 and 2.6. The density profiles are physically reasonable given the changes in the scale height and the diffusion coefficient. The solid curve in both figures indicates the case where both the scale height and diffusion coefficient are held constant. The curve marked with open triangles illustrates the case where only the scale height is allowed to vary Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 35 Figure 2.5: Martian hot oxygen density profiles derived from solutions of the diffusion equation. The hot oxygen temperature is taken to be a (constant) 8500 K , and the background temperature is a (constant) 180 K . The various curves depict the behaviour of the hot oxygen density to the variation of gravity g and diffusion coefficient D. The solid curve is for both g and D constant with altitude, the open triangle marked curve is for variation of g (and thus scale height) only with altitude, the filled triangle marked curve is for variation of D only with altitude, and the dashed curve is for the density with both g and D varying with altitude. The open circle symbols indicate the background ('cold') oxygen density. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 36 Figure 2.6: Venusian hot oxygen density profiles derived from solutions of the diffusion equation. The hot oxygen temperature is taken to be a (constant) 8500 K , and the background temperature is a (constant) 280 K . The labeling of the curves is as for Figure 2.5. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 37 with altitude. For this case, the hot atom density is maximal at altitude. Physically, this is expected given that the scale height increases with altitude as the gravitational acceleration decreases, all other parameters held constant. The case where only the diffusion coefficient is allowed to vary with altitude is indicated by the curve marked with solid triangles. The integrated hot atom source density, P(z), given in equation (2.2.4) depends inversely on the diffusion coefficient, which depends inversely on the collision frequency. Since the background density is decreasing with altitude, from equation (2.1.5) we have that P(z) decreases with an altitude dependent diffusion coefficient, and so the density at altitude decreases. The rate of decrease is larger for Venus, where the background density (and hence collision frequency) decreases more rapidly with altitude. If both the scale height and diffusion coefficient vary with altitude, as indicated by the dashed curves in the figure, the density profile lies between the values produced when only one or the other is allowed to vary. The density profiles are shown out to an altitude of 6000 km, well beyond the limited altitude range of input density profiles as given in Figure 2.2. There are several reasons for an interest in hot atom densities for such large altitudes. The first is that the enormous enhancement of the hot atom density over the thermal background does not become apparent until above approximately 500 km on Mars. Additionally, the hot atom densities at high altitude are important in determining the nature and magnitude of interactions with the solar wind. 1 4 The dependence of the hot atom density profile on temperature was also examined, and the results are illustrated in Figures 2.7 and 2.8 for Mars and Venus, respectively. Since we do not have any detailed information on the hot atom temperature profile, a constant temperature profile was used, with the value varied between 500 and 8500 K . A n increase of hot atom density at higher altitudes is visible with an increase in the hot temperature, for both Mars and Venus. The effect becomes less pronounced above 5000 Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 38 Figure 2.7: Martian hot oxygen density temperature dependence. The (constant) hot oxygen temperature was varied between 500 and 8500 K , and the scale height (gravity) and diffusion coefficient were allowed to vary with altitude. The solid, short-dashed, long-dashed, dot-dashed, and dot-dot-dashed curves represent hot oxygen temperatures of 500 K , 1000 K , 2500 K , 5000 K , and 8500 K , respectively. The open circle symbols indicate the background ('cold') oxygen density. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 39 Figure 2.8: Venusian hot oxygen density temperature dependence. The (constant) hot oxygen temperature was varied between 500 and 8500 K , and the scale height (gravity) and diffusion coefficient were allowed to vary with altitude. The labeling is as in Figure 2.7. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 40 Figure 2.9: Hot oxygen densities at large altitudes for Venus and Mars. The solid and dashed lines indicate the hot oxygen densities for Venus and Mars, respectively. The hot atom temperature is fixed at 8500 degrees for both planets. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 41 K , especially for Mars, where the low gravitational acceleration and high temperatures combine to yield scale heights of the order of thousands of kilometers. The effect of the stronger gravitational acceleration on Venus is most apparent for temperatures below 8500 K , where hot atom densities are lower than their Martian counter-parts at altitudes above approximately 1000 km. When comparing the hot oxygen density profiles of the two planets, we see that they are remarkably similar in magnitude at the 6000 km altitude of Figures 2.7 and 2.8. The hot oxygen density is greater below this altitude on Venus than Mars. The stronger dissociative recombination source at low altitudes on Venus appears to offset the effect of a lower gravitational acceleration and hence much larger scale height for the hot oxygen on Mars for altitudes up to approximately 6000 km. The hot oxygen density profiles are plotted together for altitudes up to 12 000 km in Figure 2.9. It is clear that for altitudes greater than roughly the planetary radius of Venus that the hot oxygen corona of Mars decays much more slowly due to the lower gravity, and is thus much more extensive, despite the stronger dissociative recombination source on Venus. 2.3 Non-constant Temperature Model In the previous section we assumed that the hot oxygen atoms were at some (arbi-trary) constant temperature. This is equivalent to neglecting energy transfer between the diffusing hot oxygen and the thermal background. In order to extract information about the altitude dependence of the temperature of the hot oxygen atoms in the planetary coronae, we must account for the partial thermalization of the hot atoms as they collide with the thermal background. To extend the model of the previous section to allow for temperature variation with altitude, we employ a standard hydrodynamic approach com-monly used in analysis of particle transport and temperature distribution in the middle Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 42 and upper atmosphere. As previously, we consider a system of hot atomic oxygen diffusing through a thermal background of oxygen. The thermal oxygen is assumed to be at diffusive equilibrium. We are interested in the transport of the translationally energetic, nonequilibrium hot oxygen, and the steady-state altitude dependence of its density and temperature. Following standard hydrodynamic theory used for the solar w i n d 1 0 0 , 1 0 1 we rewrite equations (1.4.1) and (1.4.2) in terms of the the random velocity, multiply by 1, mscs and | m s c 2 , and then integrate over velocity space to obtain the continuity, momentum and energy transport equations for the hot oxygen (denoted by the subscript s): Continuity: dn* dt + V-(n s u s ) = *f* + N (2.3.1) Momentum: nsms dt Energy: + namaua • V(u,) + V • Ps - n.m.G = ^ + M (2.3.2) ^ + « . - v ( | p . ) + |p.(V-u.) + V - q , + P,:Vu. = + £(2.3.3) In equations (2.3..1)-(2.3.3), ns, us, Ts, qa., and Ps are the hot oxygen density, drift (average) velocity, temperature, heat flow, and (anisotropic) pressure tensor, respectively. We have also defined N, M , E as the (altitude dependent) production terms for the hot oxygen number density, momentum, and energy, respectively. The terms and ^ on the right-hand side of equations (2.3.1)-(2.3.3) correspond to the moments of the Boltzmann collision integral, equation (1.4.2), and represent changes of density, momentum, and energy as a result of collisions. They are not shown explicitly as they depend on the type of interaction between colliding particles (e.g. hard sphere) and on Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 43 the form of the distribution functions of the interacting species. In our coupled system of transport equations, we have ignored external force terms except for that due to the gravitational acceleration, G . At lower altitudes, collisions with background species are plentiful and the velocity distribution functions for diffusing species are generally considered to be close to the equilibrium Maxwellian distribution. As collisions decrease with altitude, larger tem-perature and diffusion velocity differences can be maintained by diffusing species, and the distribution function can depart from the equilibrium Maxwellian distribution of the thermal background gas. At the lowest level of approximation, we assume that the dis-tribution function of the diffusing hot species is completely described by the first five moments of the hot atom distribution function: the density, ns, the three components of the drift velocity u s, and the temperature, Ts. We assume the pressure tensor is isotropic (off diagonal elements are zero), and that the stress tensor, heat flow, and higher order moments are zero. Thus, the system of transport equations reduces to d n ° + V • (nsus) = ^ + N(z) (2.3.4) du dt n , m s — + nsmsus • V(u s) + Wps - n3msG = ^ + M(z) (2.3.5) ^ + U * - V ( J ^ ) + fj>.(V-u.) = + E(z) (2.3.6) The collision terms on the right-hand side of equations (2.3.4)-(2.3.6) may be evaluated rigorously for arbitrary interactions between species for the case where both species' distribution functions are Maxwellian. In this case, for arbitrarily large differences in drift velocities and temperatures, we have 5jt = o ^ = nsmsvst{ut - u s)$ s t SE, _ n s m s l / s t [3k{Tt - Ts)Vst + mt(us - ut)2$st] Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 44 where $ s i and tyst are velocity dependent 'correction factors' given for hard spheres b y 9 7 ' 9 8 $st 8 est + est + (-st 4e3t 2est erf(e s i) + 1 1 + 2e2 St e x p ( - 4 ) erf(e 5 t) + | e x p ( - e 2 i ) and where the subscript t denotes a parameter for the background species, which is taken to be thermal atomic oxygen. In the equations above, erf() is the error function, eri(x) » rx 7T JO 2 ^ 7* dt ust is the momentum transfer collision frequency, _ 8 ntnst Vth&HS 3A/7T ms and Tst and pst are the reduced temperature and mass, respectively, defined as msTt + mtTs Tst = Ust ms + mt msmt ms + mt The hard sphere collision cross section is denoted by crjjS- The parameter est is the ratio of the magnitude of the difference in drift velocities to the average 'thermal' speed, that is u, tst Vt th where Vth 2kT, St Pst is defined as the 'average' thermal speed of the gas mixture. We choose to set the background drift velocity to zero (e.g. a stationary background gas). Our problem is thus reduced to finding the density, drift, and temperature profile of the diffusing species. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 45 If we assume that all motion is along the vertical (z) direction, and look for steady-state solutions, our transport equations reduce to yz(nsus) = N(z) (2.3.7) d d nsmsus-—^- + k— (n3Ts) + nsmsg = -nsmsvstus<&st (2.3.8) oz oz k 4 ( n , r , ) + W ^ = £ ( 2 ) - ! ^ x 2 oz 2 oz mt [3k(Ts - Tt)Vat - mtu2sst] (2.3.9) where we have defined the gravitational acceleration as G E — gk. If it is assumed that the drift speed is much less than the average 'thermal' speed of the gas mixture, so that c, j < 1, it may be shown that $st = 1 Vst = 1 In addition, we must specify the production terms N(z), M(z), E(z). We assume, as in the previous section, that dissociative recombination of is the only source of energetic hot oxygen atoms, so that N(z) = aDR[Ot}2 (2.3.10) where osor « (1 — 3) X 1 0 - 7 cm 3 /sec and the units of ./V are c m _ 3 s e c _ 1 . The dissociative recombination reaction is assumed to be 'single channel', with the product 5.02 eV of energy for the dissociation branch shared equally between the two product oxygen atoms. The product atoms will thus have an initial temperature Tproct given b y 6 7 \kTVT0d = 2.5 eV (2.3.11) Li and that the average energy per atom will be 3 Eav — ~^kTpr0(i (2.3.12) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 46 The energy source term E(z) is then simply given by N(z) x Eav, with units of ergs c m - 3 s ec - 1 . We also assume that the distribution of hot atom trajectories following the dissociation reaction is isotropic; in this case, hot atoms are equally likely to have upward motion as downward, and the net momentum gain is zero (M(z) = 0). As an example of such singularities, and in an attempt to gain insight into their phys-ical meaning, we rewrite our original coupled system of transport equations, equations (2.3.7)-(2.3.9), assumed a 1-D plane parallel atmospheric geometry, steady state, and neglected sources and collisions, yielding ^-{nsus) = 0 (2.3.13) dz nsmsus-^- + k— (nsTs) + nsmsg = 0 (2.3.14) dz dz \usk^- {nsTs) + l n s k T s ^ = 0 (2.3.15) 2 dz 2 dz If we now assume a constant temperature and gravitational acceleration, we may write equation (2.3.14) as dns ns ~dz~ = ~H-Cy{n*g) ( 2 - 3 - 1 6 ) where we have used, from the continuity equation, that nsus = C, constant, and that dus dns -q— = —C-Q—/n2s. We have also defined a (constant) scale height, H = kTs/msg. To solve for the density profile, we rewrite equation (2.3.16) as H-Cyin^g), -dz = - ' v s y , d n 2.3.17 nig Integrating both sides, J ' * = f n\g dn ~(z-z0)/H = ln(n,/n0) + (\ - -M 2gH \n2s n2J (z — z ) yexp(A/y2) = e x p ( - v 0 ) + A) (2.3.18) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 47 where A = C2 l{2gHn2) = u2/{2gH) and y = ns/n0. The final result is transcendental for y. For fixed values of A, plots of y exp(A/y2) vs. y (or exp( — (*~ °^) + A) vs. z) may be drawn. Examples are given in Figure 2.10. The gravitational acceleration and scale height are held constant at 3.3 m/s 2 and 30 km, respectively, and are chosen to represent typical values near the Martian exobase for atomic oxygen. The flow speed u0 is chosen to be 300 m/s and 30 m/s, respectively, for curves labeled A and B , yielding values of the parameter A equal to 0.454 and 4.54x 1 0 - 3 , respectively. The curve defined by equation (2.3.18) has a minimum at ymin = \/2A. Using this value in equation (2.3.18), we may solve for z = zmax, Zmaoc = [A- m(v/2Ae)] H + z0 Thus, for given values of the hot oxygen density, ns, and A, there are altitudes above which the solution to equation (2.3.18) does not exist, and forcing the numerical algorithm past these altitudes results in meaningless solutions. Graphically, this is illustrated in Figure 2.10 by the continuation of the solid curves (the RHS of equation (2.3.18) below the minimum of the dashed curves (the LHS of equation (2.3.18). We may rewrite our conservation equations for the above case in a slightly different form in order to gain further insight into the source of the singularity. For the case of no sources, no collisions and constant temperature, we had conservation and momentum equations given by : d , — (nsus) = 0 Oz n s m s u 3-7r— + k—(nsTs) = - nsmsg Oz Oz Substituting the continuity equation into the momentum equation, we have kTsC\ dus nsmsus — = - namsg u2 ) dz where we have made use of Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 48 Figure 2.10: Solutions of the collisionless, sourceless momentum conservation equation for a plane parallel geometry. The gravitational acceleration and scale height are held constant at 3.3 m/s 2 and 30 km, respectively, and are chosen to represent typical values near the Martian exobase for atomic oxygen. The flow speed uQ is chosen to be 300 m/s and 30 m/s, respectively, for curves labeled A and B. The solid lines represent plots of exp(— (z~^°> -\-A) vs. z from equation (2.3.18), with the z axis given by the top horizontal axis. The dashed lines represent plots of y exp(A/y 2 ) vs. y from equation (2.3.18), with the y axis given by the bottom horizontal axis. The vertical axis is common to both plots. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 49 nsus - C, ^ = -C^-jnl dz dz After further manipulation, we may rewrite our momentum equation in the form U - — Y i 1 = -«.9 (2-3.19) \ ms J dz The important point to note from equation (2.3.19) is that while the left-handside may vanish for either du/dz = 0 or u2 = kT/m, the right-handside vanishes only if u = 0. That is, the zeroes on the LHS are not balanced on the RHS, and the two sides vanish independently of each other. This imbalance may be resolved by following the approach used in the examination of the solar w i n d . 1 0 2 ' 1 0 0 , 1 0 1 Rewriting our conservation equations in spherical coordinates, we have 1 d / 2 \ r> ~r~7T~ [nsusr I = 0 H dr v ' d U s . u 9 i t \ n s m s u s - r — + k—(nsls) = — namsg — nsmsvstus dr dr 3 , d , m . 5nskTsdusr2 nsmustvst^, / m m^ 2 dr 2 r 2 dr mt For the case of no sources, no collisions and constant temperature, we may rewrite our momentum equation as for the plane parallel case, yielding / 2 kTs\ dus ( 2kT,\ ms J dr s ^ 3 msr J We may now force the left- and right-handsides to vanish identically by requiring that u 2 = kT/m at the point r = This point is known as the 'critical point' of the solution of the flow velocity. As we add more parameters to our system of conservation equations (e.g. temperature, heat flow) we alter the nature of such singularities, removing some and adding others. 1 0 3 The system of equations given by (2.3.7)-(2.3.9) have excluded such physical processes as heat flow (conduction). As the conductivity (and hence the heat flow) in the upper Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 50 atmosphere is large, it may be that heat flow terms cannot simply be discarded from the energy equation. In fact, such terms may dominate the drift terms, and it was decided that a model which incorporated heat flow should be investigated. We choose to close the transport equations including the heat flow with phenomenologic laws rather than with the introduction of another, higher order transport equation. For example, the phenomenologic heat flow equation (Fourier's Law) may be wr i t t en , 9 8 ' 1 0 4 q s = - A s < V T t - X'tsVTs + Rst(us - ut) (2.3.20) where the various coefficients are detailed in Appendix B , and are strongly dependent on the form of the interaction between species in the gas mixture. It has been shown 1 0 5 that for a range of altitudes spanning the thermosphere and above, the Navier-Stokes (Fourier's Law) expression for the heat flow vector q s = - A , t V r f (2.3.21) may be used in place of equation (2.3.20). Inserting this term into the energy equation, and neglecting nonlinear flow terms in the momentum and energy equations as 'small', our system of equations now becomes d dz d (nsus) = N(z) (2.3.22) k-£- (nsTs) = - nsmsg - nsmsvstus (2.3.23) oz = - = ^ 3 ^ . - 7 1 ) + B W (2.3.24) Previously, we looked for a solution of the continuity equation for the density, given fixed temperature, for the case where the flux was given by the expression for a minor species diffusing in a background and neglecting thermal diffusion, that is equation (2.2.1), ds d z = N{z) s = - D s i dns ns (2.3.25) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 51 where Dst is the mutual diffusion coefficient, and Hs is the scale height for species s. It may be shown that when nonlinear terms in the momentum equation are ignored that our coupled system of equations (2.3.22)-(2.3.24) reduce to equation (2.3.25) for the case where the temperature of the diffusing species is constant with altitude. From equation (2.3.8), we have for this case that us = -Dst ns Hs (2.3.26) where the flux is defined as s = nsus, and the mutual diffusion coefficient Dst and scale height Hs are defined as D St msvst H, = «• msg respectively. Multiplying equation (2.3.26) by the density ns yields the result of the previous section, equation (2.3.25). Similarly, from inspection it is seen that in the case where heat-flow is retained and nonlinear flow terms are discarded, that the momentum equation (2.3.23) also reduces to equation (2.3.26). The continuity, momentum, and energy equations given by equations (2.3.7)-(2.3.9) form a coupled system of simultaneous differential equations. As both the scale height and diffusion coefficient are temperature dependent, they become altitude dependent if the temperature is not held constant. We considered solving the system using a simple 'stepping' finite differencing; that is, we rewrite the equations so that a simple first or-der forward-difference gave us the value of either the density, drift, or temperature at the next altitude step in terms of the values of the other parameters at the previous (current) altitude step. Thus, given a set of input values at the 'bottom' of the atmo-sphere, one could simply 'step' up in altitude to get the new values of each parameter, governed by the transport equations. The difficulty in this method appeared to be a Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 52 strong instability/discontinuity at an altitude where the drift velocity became equal to some equation-dependent 'sonic' speed, analogous to the singular point(s) in the solar wind momentum equation. It appears such singularities are inherent to the truncated transport equations used in transport theory, 1 0 3 and that there is some difficulty (in gen-eral) in the determination of the singular points and integrating/stepping a solution past them. In recasting the equations with the inclusion of the heat flow (and exclusion of the nonlinear drift terms), we hoped to avoid the particular singularity which hampered the direct 'stepping' method. The equations given by (2.3.22)-(2.3.24) can be put in the non-dimensional form ^ = N(0 (2-3.27) ON _ A/" msgz0 _ F msu0z0vst _ • N_ ro q oo \ d£ ~ T kT0 T kT0 TV { l - 6 - ^ ] dT q£ = V (2.3.29) OV *kn0zlw« E{z)zl dX dt ~ X(T)mt N { T T t / T ° > ~ X(T)T0 ~ X(T)T°dTs { 2 - 3 M ) where we have used 7Y = ns/n0 U = us/u0 T = Ts/T0 ' t = z/z0 F = N-U V = dT dt N = N{z)z0 u0n0 KT.) = Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 53 We now have the four unknowns F, N, T, V, where F = N * U. It should be noted that we have chosen to solve for F, and not U as it simplifies the form of the equations by eliminating the appearance of derivatives on the R H S of equation (2.3.27). The set of 'heat-flow' equations was integrated using a standard fourth-order Runga Kut ta integrator, with a set of initial conditions specified at the lower boundary only. The lower boundary altitude was selected such that we could assume that there was no production of hot atoms; hence, N, U, F at the lower boundary were set to zero. The temperature T was set equal to the background temperature. The source terms are the same as for the previous coupled transport equation case. When the temperature was set to a constant it was found that the system yielded results for the density profile that were in rough agreement with those found by direct integration of the diffusion equation (as for the constant temperature model results earlier in this chapter). Transport models based on moment equations derived from the Boltzmann equation are regularly used to address flow problems in planetary atmospheres and in the solar wind. These transport models consist of truncated sets of coupled differential equations involving moments of the distribution function of the flowing species. As the number of moments retained in the expansion increases, large departures from equilibrium for the distribution function can be treated more accurately. However, the choice of the level of truncation also determines the number of singular points in the system of equations. These points are of importance in both steady-state and time-dependent systems. 1 0 3 While some of these singular points are well known (e.g the sonic point in the 'solar wind' equations 1 0 0), in general they must be investigated uniquely for each particular system. These singular points cannot be integrated over using standard methods, and require special attention, which again must be applied on a system-by-system basis . 1 0 3 For simpler systems, it may be possible to easily identify the singular points (e.g. equa-tion (2.3.19)), but for more complex system such a simple treatment does not necessarily Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 5 4 stabilize the numerical system. It is clear from the behaviour observed when numeri-cally integrating equations, although the nature of the singularity/singularities and their treatment is not certain. 2.4 The Linear Boltzmann Equation Model The problem of thermalization of translationally energetic particles in a moderator of atomic or molecular constituents is a well known problem in nuclear physics. 1 0 6 A n example is the problem of the moderation of high energy neutrons resulting from fission reactions. Energetic neutrons are created with a narrow range of initial energies, and diffuse through a moderator, losing energy through collisions with the moderator nuclei. If the moderator is composed of nuclei of moderate to high mass, the average spread of final energies between individual neutrons after some number of collisions is relatively small, and it is possible to characterize the movement of a large number of neutrons through the moderator by means of some average velocity. The loss of energy to the moderator is then viewed as a continuous process, since it is assumed that the average energy (velocity) change due to collisions with moderator nuclei at any given time is minimal. The approximation of thermalization as a continuous process is known as the Continuous Slowing Down Approximation, or CSDA. Application of the CSDA is not limited to neutron moderation; its validity has been in explored in electron degradation problems as w e l l . 1 0 7 ' 1 0 8 We may model the diffusion of hot oxygen atoms resulting from the dissociative re-combination of ions in a manner similar to the CSDA for neutrons in a moderator. The height distribution of the hot oxygen atoms and their thermalization through colli-sions with the ambient background atmosphere is analagous to the 'depth of penetration' Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 55 and rate of thermalization of energetic neutrons moving through a moderator. The den-sity and temperature of the hot oxygen atoms may be self-consistently calculated as moments of the distribution function of the hot oxygen atoms. The hot oxygen distribu-tion function is given by the Boltzmann equation, equation (1.4.1). In the current work, we choose to neglect external force fields (such as magnetic and gravitational fields). In addition, with the CSDA the velocity appearing in the second (or diffusion) term on the left-hand side of equation 1.4.1) is a assumed to be constant, with a value equal to the average speed of the hot atoms. The time and diffusion terms may thus be combined to yield an equation of the form % = -![/] + S (2.4.1) where S is the (steady) source of the translationally energetic atoms, J is the Boltzmann collision operator, given by equation (1.4.2), describing elastic collisions between the hot and (single) background species, and / — f(c,t) is the velocity distribution function for the hot species. Equation (2.4.1) is a local model for the time evolution of the distribution function at some (fixed) altitude. We have assumed that the moderating background gas is in equilibrium and is distributed with a Maxwellian distribution of velocities. In applications of the CSDA to nuclear physics, hot particle energies are of the order of keV-MeV, and the distribution function of the hot atoms remains sharply peaked about the average velocity. For the present case, the hot particle energies are of the order of eV, and so we expect that there will be some broadening of the hot atom distribution function as a result of collisions with the background. The Boltzmann equation given by equation (2.4.1) yields the hot oxygen distribution function as a function of time, at some particular altitude. This altitude is a parameter in the model, and dictates the strength of the source of hot oxygen atoms, the background densities, and other physical parameters for the model. If we assume, as in the CSDA, Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 56 that there is limited thermalization of the hot oxygen atoms and that their distribution function remains relatively narrow, then the average velocity of the hot oxygen atoms is approximately constant. Thus, the time dependent distribution function may be trans-formed to an 'altitude dependent' distribution function through the linear transformation given by z = vhott (2.4.2) where Vhot is an average velocity characterizing the motion of the hot atoms. With this transformation, we may consider the time variation of the distribution function as an altitude variation. However, since the CSDA requires that the hot atom distribution function not thermalize (or broaden) 'too much', we are limited to altitudes 'close' to the true altitude chosen to parameterize the local Boltzmann equation. In order to extend the model to yield the altitude distribution of hot atom densities high into the exosphere, we adopt the method used in many collisionless models of the exosphere,16 detailed shortly. The form of the hot oxygen source distribution is largely unknown, 6 ' 8 and is assumed to be a narrow Gaussian centered about the O2 dissociative recombination energy of 2.5 eV per atom. The altitude dependent source strength is given by equation (2.1.2). The normalization of the source function is chosen equal to the source strength. The initial (t — 0) hot atom distribution function is set equal to a Gaussian of the same form as the source Gaussian, and is normalized to 1. With our steady source of hot atoms, we expect that any final steady state solution should exhibit a departure from Maxwellian. We may rewrite the collision operator J in equation (2.4.1) in its kernel representation,109' yielding ^ - J~dx'{K(x,x')f(x',t)}-Z(x)f(x,t) + S(x) (2.4.3) where x = rrihotc2/2kTb is the hot atom energy normalized to the background energy, and K(x,x') is the non-symmetric collision kernel (that is, K(x',x) ^ K(x,x')). The Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 57 (elastic) collision frequency between hot atoms and the background, Z{x), is given by Z{x) = / K(x,x')dx' (2.4.4) Jo or by Z(x) = JJ fa(g,n)gdfldc in order to ensure that detailed-balance or conservation of particle number is satisfied. Integrating equation (2.4.1) over x, we find that (2.4.5) where /•oo / 'S(x) dx = S Jo Integration of equation (2.4.5) yields n(t) = St (2.4.6) This model gives a density that is linearly increasing with time, with a rate of increase directly proportional to the strength of the source at the altitude of interest. The den-sity may also be calculated directly from the normalization of the velocity distribution function, roo n(t) = / f(x,t)dx (2.4.7) Similarly, it may be shown that the average reduced hot atom energy is given by /•oo E(t) = / xf(x,t)dx (2.4.8) Jo Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 58 2.4.1 Discretization Procedure and Solution We solve the Boltzmann equation using the quadrature discretization method (QDM), the basis of which is described in detail elsewhere111 . The method involves the determination of the velocity distribution function, f(y), at a discrete grid of points, {yi}. The grid points correspond to a set of quadrature points which are the roots of a set of orthogonal polynomials based on the weight function w(y) = y2exp[—y2]. We may thus write N / e-yy2g(y)dy = Yw™9(ym) m=l or J fOO g(y) dy 0 N Y w m g ( y m ) m=l where Wm = wmev™/y2m We rewrite equation (2.4.3) as a function of reduced speed y = yfx, noting the change in the volume element from dx to 2y dy, and then apply our QDM discretization to yield QJ^L = £ Wm2ymK(y2m,y2)fi(y2,t) - Z{y2)U{y2, t) + S(yf) (2.4.9) O Z m=l If the grid of points y\ at which the distribution function is evaluated are chosen to coincide with our grid of speed points ym then equation (2.4.9) is reduced to a set of linear equations which may be represented in vector form as 8t where (9f -w: = M-f + S (2.4.10) M M , , = Wm 2ym Kn\y2m,y2) - Z(yi) 6(ym - yi) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 59 For a hard sphere collision model, Kns(x, x') is the well k n o w n 1 0 9 , 1 1 2 non-symmetric Wigner-Wilkens kernel, Kns(x, x') = ^-AQ^ieviiQVx1 + Ry/x) + ex~x' evi{Ry/x~' + Qy/x) LJ V X ±[evi(Qy/x1 - Ry/x) + ex~x' ed(Ry/x~' - Qv*)]} (2.4.11) where the +(—) sign refers to x > x' (x < x'), and Q = R = x / 7 - V 7 7 = mback/mhot A = crtotnba, l2kTbai ck\ ck with erf(x) the standard error function. It can be shown (see Appendix B.2) that the non-symmetric kernel in equation (2.4.11) may be defined in terms of a symmetric kernel K'(x,x% Kn'(x,x') = Transforming to speed, we then have Kns(y2m,y?) = \J x% Ks(x,x') yi (2.4.12) (2.4.13) To complete the symmetrization of M requires that we make a few further changes to the form of equation (2.4.10). These are detailed in Appendix B.3, with the result that equation (2.4.10) becomes § - = B-f + S OT (2.4.14) where Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 60 and B is a symmetric matrix defined by equations (B.3.11) and (B.3.13). In equation (2.4.14), we have defined a dimensionless 'time' r = At and where / W fi = \l^ffi ( 2 - 4 - 1 5) The solution of equation (2.4.14) is relatively straightforward, and is detailed in Appendix B.4. The final result is fi{r) = Uitl [^(0) + Qir] + pUit3 (FJ(0) e~x>T + Q , - J l l i — J . \ (2.4.16) where U" 1 = U T Q = U TS A, > 0 Fi(r) = U T / 4 (r ) and equation (2.4.15) is used to transform between / and / . As mentioned previously, to determine the altitude density profile, as given by equa-tions (2.4.17)-(2.4.19), we require the Maxwellian distribution function for the hot atoms at the exobase level. Hence, we need to determine the change in the distribution func-tion between the altitude of production of the hot atoms and the exobase. To accomplish this, we make a transformation from time to altitude in equation (2.4.16). A n initial altitude near the production peak and below the exobase is selected to parameterize the local Boltzmann equation. Equation (2.4.16) is then used to find the time dependent velocity distribution function for the hot atoms produced at this altitude. The timescale Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 61 over which we follow the development of the distribution function is dependent on the distance between the production altitude and the exobase (that is, we assume that using our transformation from time to altitude approximates following hot atoms from their production at some initial altitude to the exobase altitude). The distribution function at the exobase level is fitted to a Maxwellian distribution function parameterized by some density and temperature, which are used to calculate the density profile as a function of altitude above the exobase as follows. In the absence of collisions, particles moving through an exosphere follow trajectories determined by the direction and magnitude of their velocities and the strength of the gravitational field. Based on the various types of possible trajectories, classes of particles may be defined. Of interest in our particular case are those classes of particles which can give rise to an extended population of hot atoms. The exospheric density profile is calculated by partial integration of the hot atom Maxwellian distribution function at the exobase. The integration is carried out over the permitted ranges of velocity and angle corresponding to those particle classes which lead to extended hot atom populations. The details of the particle classifications and the calculation of the densities is discussed by Chamberlain 1 6 ' 4 3 and reviewed in detail by Fahr and Shizgal,1 and is detailed in Appendix A. The resulting expressions for the densities of the various exospheric components are nb = -^=nbar 17(3/2, Xcy) - y/l - ?/2exp ne = ^ { [ T ( 3 / 2 ) - 7 (3 /2 , Xcy)] ~-Xcy2 W i 7(3/2, Acy/(l + y)) 1(2.4.17) 1 + 2/ - \Jl - y2exp i + y F(3/2)-7(3/2, Acy/(l + y))H (2.4.18) nbs = 4=^r7(3/2, Acy) (2.4.19) y/TT where y = rc/r, the ratio of exobase position to the exospheric radial position being con-sidered, A c = Eescj{kTc) is the exobase escape parameter given previously by equation Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 62 (1.2.10), ribar — nc exp [—Ac(l — y)] is the barometric density distribution given previ-ously by equation (1.2.4), and T and 7 are the standard gamma and incomplete gamma functions, respectively. The subscripts b, s, e refer to ballistic, satellite, and escaping components, respectively. In order to calculate a time dependent hot oxygen density and temperature, we solve for the hot oxygen distribution function, from equation (2.4.16), and then take the appro-priate moments of the distribution function (equations (2.4.7) and (2.4.8), respectively). The integrations are done numerically using the Q D M points and weights mentioned previously. A hard sphere (total) collision cross section of 2 . 0 x l 0 - 1 5 c m 2 is used. These calculations are carried out for both Venus and Mars. Figure 2.11 shows the time de-pendent hot oxygen density and temperature for Mars, as calculated for two different altitudes, 150 and 220 km. The time dependent hot oxygen density and temperature for Venus, as calculated for the altitudes of 150 and 180 km, are shown in Figure 2.12. For both cases, we note that the density is a linear function of time, as predicted by equation (2.4.6) for a steady source. We also note the sharp exponential decay of the temperature from the initial value of 2.5 eV toward the background oxygen value of 180 K and 280 K for Mars and Venus, respectively. In both cases, the rate of decrease of the temperature of the hot atoms is greater for the lower altitude. Physically, this conforms to our expectation that the background atoms thermalize the hot atoms via elastic collisions much more quickly at lower altitudes, where the background density is higher. Note that while the initial drop in temperature for the hot oxygen atoms is large, the temperature does not fully equilibrate with the background because energy is continually added to the hot oxygen population via the steady dissociative recombination source. The distribution functions are shown at both altitudes for several times in Figure 2.13 for Mars and Figure 2.14 for Venus. The total number of particles will grow (linearly) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 63 0 200 400 600 800 1000 1200 t (sec) (b) Figure 2.11: Time dependence of hot oxygen density and temperature on Mars. The graphs are for altitudes of (a) 150 km and (b) 220 km. The time t is calculated from t = T/A, where A = 7ro?2n;,aefcyj^ m™'* 15 altitude dependent. The value of A and the source strength are (a) 0.0183 & 1199.2 and (b) 0.00157 k 7.08, in units of sec"1 and c m " 3 sec - 1 , respectively. The solid and dashed curves indicate the time dependent temperature and density, respectively. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 64 (a) 20000 8 8. E (b) Figure 2.12: Time dependence of hot oxygen density and temperature on Venus. The graphs are for altitudes of (a) 150 km and (b) 180 km. The time t is calculated from t = T/A, where A = ird2nhacksJ2'^^ is altitude dependent. The value of A and the source strength are (a) 0.25 k 28747.7 and (b) 0.037 k 2187.5, in units of sec"1 and c m " 3 sec - 1 , respectively. The solid and dashed curves indicate the time dependent temperature and density, respectively. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 65 with time due to the steady hot oxygen source. The much higher 0^ density on Venus results in a significantly higher rate of production of hot oxygen particles; however, the background thermal oxygen density is also larger, and so thermalization is also increased. The distribution function near the equilibrium (background) energy grows rapidly on both planets, once production and thermalization of hot oxygen atoms reaches steady-state. The steady-state distribution function will thus be bi-modal, with a peak about the production energy and a peak about the thermal energy. The degradation in energy of atoms from their initial hot production energy to a peak about the thermal (background) energy is clearly visible. The shape of the distribution function as it evolves in time is similar for the two altitudes illustrated. However, the rate of evolution of the hot oxygen distribution function and the number density of hot oxygen atoms (as given by the area under the distribution function curve) varies significantly between the two altitudes for both planets, as is clear from Figures 2.13 and 2.14. Previous workers 1 1 ' 7 ' 6 ' 3 3 ' 8 employed Monte Carlo simulations and their energy distri-bution functions are plotted as histograms for energy 'bins' of discrete widths of 0.05 eV. We compare with their results by taking the value of our distribution (as a function of energy) at the midpoint of a given energy bin and multiplying by the bin width. This yields the hot atom density per bin width at the exobase, as shown in Figures 2.15 and 2.16. It should be noted that a residual 'peak' remains at the dissociation energy of the hot oxygen, 2.5 eV. This is because we consider a constant source of hot atoms. A Maxwellian (or sum of Maxwellians) is fitted to the energy distribution, yielding an exobase temperature and density. The best fit for the altitude of 150 km for both Mars and Venus is shown in Figures 2.17 and 2.18. For Mars, the fit is comprised of two Maxwellians, a cold one with density 3 500 c m " 3 and temperature 800 K , and a hot one with density 800 c m - 3 and temperature 7 500 K . For Venus, the fit is also comprised of two Maxwellians, a cold one with density 10 000 c m - 3 and temperature 800 K , and a Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 66 (b) Figure 2.13: Time evolution of the Martian hot oxygen distribution function from an initial Gaussian distribution at 2.5 eV due to a steady source at that same point. The graphs are for altitudes of (a) 150 km and (b) 220 km. The solid, short dashed, long dashed, and dot-dashed curves correspond to dimensionless times of r = 0.1, 0.2, 0.4, and 1.0 in (a) and r = 0.1, 0.2, 0.4, and 1.9 in (b). The time t is calculated as described in Figure 2.11. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 67 " i — 1 — r (a) (b) Figure 2.14: Time evolution of the Venusian hot oxygen distribution function from an initial Gaussian distribution at 2.5 eV due to a steady source at that same point. The graphs are for altitudes of (a) 150 km and (b) 180 km. The solid, short dashed, long dashed, and dot-dashed curves correspond to dimensionless times of T = 0.1, 0.2, 0.4, and 1.0 in (a) and r = 0.1, 0.2, 0.4, and 1.9 in (b). The time t is calculated as described in Figure 2.12. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 68 0 1 2 E ( e V ) (b) Figure 2.15: Energy distribution of hot oxygen on Mars. The graphs are for altitudes of (a) 150 km and (b) 220 km. The dotted, dashed, and solid curves are for dimensionless times r = 0.4, 1.0, and 1.9, respectively. The results of Ip 6 and Lammer and Bauer 8 are shown for comparison. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 69 Figure 2.16: Energy distribution of hot oxygen on Venus. The graphs are for altitudes of (a) 150 km and (b) 180 km. The dotted, dashed, and solid curves are for dimensionless times T = 0.4, 1.0, and 1.9, respectively. The result of Ip 6 is shown for comparison. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 70 hot one with density 1 200 c m - 3 and temperature 7 500 K . We use equation (2.4.17) to calculate the hot atom density profile for altitudes above the exobase. The results for the ballistic component are shown in Figure 2.19 for Mars, together with the results of Ip 6 and Lammer and Bauer, 8 and in Figure 2.20 for Venus, together with the results of Ip 6 and Nagy and Cravens. 7 The difference in predicted hot atom densities between Venus and Mars originates in the escape parameter, A c , which is the controlling parameter in the collisionless exospheric model, equations (2.4.17)-(2.4.19). From the Maxwellian fits to the energy distributions, we extracted a hot component temperatures of 7500 K for hot oxygen on Venus and Mars, respectively. Referring to Table 1.1 for escape energies for oxygen from Venus and Mars, we thus find that the escape parameter, A c = Eesc/(kTc), is approximately 13.3 for Venus and 3.0 for Mars. This is crucial since the exospheric density profile is approximately given as a barometric distribution with an exponential decay parameterized by A c Thus, even though the hot oxygen density on Venus is higher at the exobase level, it decays much more rapidly with altitude due to the much higher escape parameter, yielding the observed discrepancy between predicted hot oxygen densities for the two planets. This discrepancy would be expected to grow significantly larger at higher altitudes. For both Mars and Venus there is a discrepancy between the currently calculated exospheric densities and the profiles of previously published results. The discrepancy is on the order of 2-10 times smaller than other densities for Mars, and of the order of 2-30 times smaller for Venus. It should be noted, however, that order of magnitude differences exist between models of previous workers. In order to assess possible reasons for the discrepancies between the current model results for the hot oxygen densities and those of previous workers, it is first important to compare and contrast the assumptions and input parameters for each. The input density profiles are identical for each, as is the choice and magnitude of a hard sphere collision cross section between the hot oxygen atoms Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 71 10 1 1 1 i I i I 0 1 2 3 E(eV) Figure 2.17: Maxwellian fit to the hot oxygen energy distribution function at 150 km. The fit, indicated by the plus symbols, is composed of two separate Maxwellians, a cold one with density 3 500 c m - 3 and temperature 800 K , and a hot one with density 800 c m - 3 and temperature 7 500 K . The results of Ip 6 and Lammer and Bauer 8 are shown for comparison. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 72 E(eV) Figure 2.18: Maxwellian fit to the Venusian hot oxygen energy distribution function at 150 km. The fit, indicated by the plus symbols, is composed of two separate Maxwellians, a cold one with density 10 000 c m - 3 and temperature 800 K , and a hot one with density 1 200 cm and temperature 7 500 K . The result of Ip 6 is shown for comparison. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 73 Figure 2.19: The hot oxygen density profile for Mars derived from the ballistic component calculated from the Maxwellian fits to the energy distribution function at 150 km. The densities resulting from the 'hot' and 'cold' Maxwellian fits are given by the short and long dashed curves, and their sum, the total density of the fit, is given by the solid curve. See Figure 2.17 for details on the fit parameters. The results of Ip 6 (IP) and Lammer and Bauer 8 (LB) are shown for comparison. Chapter 2. Diffusion and Boltzmann Equation Models o f Hot Oxygen Coronae 74 1600 1200 6 T 3 3 800 400 h 1000 Density (cm"3) 10000 Figure 2.20: The hot oxygen density profile for Venus derived from the ballistic compo-nent calculated from the Maxwellian fits to the energy distribution function at 150 km. The densities resulting from the 'hot' and 'cold' Maxwellian fits are given by the short and long dashed curves, and their sum, the total density of the fit, is given by the solid curve. See Figure 2.18 for details on the fit parameters. The results of Ip 6 (IP) and Nagy and Cravens 7 ' 3 2 (NC) are shown for comparison. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 75 and the background. The exobase altitudes are identical. The Monte Carlo and present models both approximate the source as fixed at a single altitude, with a dissociation energy of 2.5 eV available to each hot product oxygen atom. The current model assumes a Gaussian spread of dissociation energies about the dissociation energy while the Monte Carlo work utilizes a delta function at the dissociation energy. The current model utilizes the Boltzmann collision operator to account for thermalization of the hot oxygen atoms by the thermal background; in the Monte Carlo work, collisions and the amount of energy exchanged between the hot atom and the background per collision are treated as stochastic processes, and are calculated using random number generators and altitude dependent background density profiles. In neither the work of Ip 6 or Lammer and Bauer 8 is there a detailed description of the statistics used in generation of their final energy density distribution functions and densities, so it is not possible to comment further on the specifics of their simulations. In the current model, as well as in the two stream models of Nagy and Cravens7 and Monte Carlo models of Ip 6 and Lammer and Bauer, 8 the exospheric density profile for the hot oxygen atoms is assumed to be given by the ballistic component of the collisionless exospheric density model, equation (2.4.17). The density and temperature of the hot oxygen distribution function fit to a Maxwellian distribution at the exobase level is required in order to calculate this ballistic density profile. Aside from the present work, only Lammer and Bauer explicitly indicate fitting the exobase distribution function with a Maxwellian distribution in order to extract an altitude density profile. Their work is also the only work to give the exobase altitude density and temperature used as input for the collisionless ballistic density profile, citing a two component fit consisting of a hot component (2 000 c m - 3 and temperature 7 500 K) and a cold component (15 000 c m - 3 and temperature 850 K) . The current model calculates the collisionless exospheric densities using equation Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 76 (2.4.17). The accuracy of this formulation was checked numerically for Earth against previously published results in Fahr and Shizgal 1 and Chamberlain, 1 6 and excellent agree-ment was attained. We were also able to reproduce the results of Lammer and Bauer 8 for the Martian hot oxygen density profile using their exobase density and temperature fit parameters. It was thus concluded that the calculation of the exospheric densities from equation (2.4.17) was not in error in the current work. The most probable remaining source of error is the collision model used to describe the thermalization of the hot oxygen atoms. The probabilistic approach of Monte Carlo is replaced in our kinetic theory method by the time dependent Boltzmann equation. We generate the exobase altitude energy distribution function without need for fitting. The model calculates the time dependent variation of the velocity distribution function at a single, fixed altitude and uses the linear transformation given by equation (2.4.2) to ap-proximate altitude variation of the distribution function. The effect of thermalization is exaggerated since the background density remains fixed at the chosen altitude, ignoring the decrease in background density with altitude. While the Boltzmann equation model does account for collisional thermalization of the hot atoms with the background it does not properly describe diffusive transport of hot atoms over the range of altitudes over which production is non-negligible. The true exobase velocity distribution function of the hot atoms includes the effect of hot atom production from a range of altitudes. In addition, since our transformation between time and altitude is linear, with a proportion-ality constant given by the (constant) average velocity, we require that the hot atoms not thermalize 'too much' (and thus change their average velocity). This limits the model to regions 'close' (in a mean free path sense) to the exobase. As seen from Figure 2.11 for Mars, and especially in Figure 2.12 for Venus, even with these limitations the ther-malization (and thus change in the average velocity) of the hot atoms is not negligible, and so the simple linear transform of equation (2.4.2) does not adequately describe the Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 77 diffusion process of the hot oxygen atoms. A significant approximation in both this model and in the Monte Carlo simulations is that the atmosphere above the exobase level is collisionless, and that equations (2.4.17)-(2.4.19) adequately describe the altitude distribution of hot atom densities in the ex-osphere. The validity of this approximation is examined in the following section by a simple comparison of timescales describing the rate at which diffusion and thermalization proceed. 2.5 Collisional and Diffusional Timescales In order for the hot atoms to form a planetary corona, they must travel from their place of creation (generally low in the atmosphere) into the upper parts of the exosphere. This requires transport from a generally collisional regime near the production peak to a (nearly) collisionless one in the upper exosphere. Diffusion acts to transport hot atoms away from the source region while collisions thermalize the hot atoms. The relative strength of these processes can be compared by examining timescales representative of the speed at which they take place. A collisional timescale is given by where v, the hot atom-background collision frequency, is defined previously in equation (2.1.5). A measure of the 'diffusive timescale' is given by As v increases, TCOU(Z) decreases, but D{z) decreases, and Tdiff(z) increases. A n example of how these processes compare is given in Table 2.3 for Mars and Table 2.4 for Venus. The first column of the tables shows the altitude, the second and third the timescales Tcoll{z) = 1/v (2.5.1) rdlff(z) = H2(z)/D(z) (2.5.2) Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 78 Altitude (km) Tcoii (sec) T d i f f (sec) Tcoll/Tdiff 150.0 3.94e+00 9.64e+04 4.08e-05 200.0 2.49e+01 1.61e+04 1.55e-03 250.0 1.28e+02 3.30e+03 3.90e-02 300.0 7.31e+'02 6.12e+02 1.19e+00 350.0 4.16e+03 1.14e+02 3.66e+01 400.0 2.36e+04 2.11e+01 1.12e+03 450.0 1.34e+05 3.90e+00 3.45e+04 500.0 7.65e+05 7.23e-01 1.06e+06 550.0 4.35e+06 1.34e-01 3.25e+07 600.0 2.47e+07 2.47e-02 1.00e+09 650.0 1.41e+08 4.57e-03 3.08e+10 700.0 8.01e+08 8.44e-04 9.49e+ll 750.0 4.55e+09 1.56e-04 2.92e+13 800.0 2.59e+10 2.87e-05 9.02e+14 850.0 1.47e+ll 5.29e-06 2.78e+16 900.0 8.38e+ll 9.75e-07 8.59e+17 950.0 4.77e+12 1.80e-07 2.65e+19 Table 2.3: Diffusional versus collisional timescales for Mars. The values are for 8500 K oxygen diffusing through a background of 180 K oxygen. A hard sphere (total) collision cross section of 2 . 0 x l 0 - 1 5 cm 2 is used. The Martian exobase for oxygen is located at approximately 250 km. for collisions and diffusion, respectively, and the final column shows the ratio of the timescales. From the tables it is clear that at low altitudes collisions are extremely important in altering the hot atom distribution before there is much chance to diffuse from the source region. For both Mars and Venus, the ratio of the timescales for collisions is less than or equal to the the diffusional timescale for altitudes exceeding the exobase level. We would thus expect that thermalization of the hot atom population would continue beyond the exobase level, and that the final hot atom energy distribution should be somewhat more thermalized than the collisionless exosphere models would predict. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 79 Altitude (km) Tcou (sec) Tdiff (sec) Tcoll/Tdiff 150.0 3.55e-01 1.74e+05 2.03e-06 200.0 7.75e+00 8.25e+03 9.40e-04 250.0 1.34e+02 4.93e+02 2.71e-01 300.0 2.00e+03 3.41e+01 5.86e+01 350.0 3.03e+04 2.32e+00 1.30e+04 400.0 4.59e+05 1.58e-01 2.91e+06 450.0 6.96e+06 1.07e-02 6.47e+08 500.0 1.05e+08 7.31e-04 1.44e+ll 550.0 1.60e+09 4.97e-05 3.22e+13 600.0 2.42e+10 3.38e-06 7.17e+15 650.0 3.68e+ll 2.30e-07 1.60e+18 700.0 5.57e+12 1.56e-08 3.57e+20 750.0 8.45e+13 1.06e-09 7.97e+22 800.0 1.28e+15 7.20e-ll 1.78e+25 850.0 1.94e+16 4.89e-12 3.97e+27 900.0 2.94e+17 3.32e-13 8.86e+29 950.0 4.46e+18 2.25e-14 1.98e+32 Table 2.4: Diffusional versus collisional timescales for Venus. The values are for 8500 K oxygen diffusing through a background of 280 K oxygen. A hard sphere (total) collision cross section of 2 . 0 x l 0 - 1 5 cm 2 is used. The Venusian exobase for oxygen is located at approximately 200 km. Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 80 The strength of the dissociative recombination reaction is directly proportional to the density of ambient Of . Since the O f density peaks relatively low in the both the Venusian and Martian atmospheres (about 140 km for both) and decays rapidly with altitude, the production of hot oxygen is non-negligible only over a very small region of altitudes. Thus, for both Venus and Mars we would expect that the predominant source of exospheric hot oxygen from dissociative recombination would be near 140 km (see Figures 2.2 and 2.3). However, it is important to consider the collisional and diffusional timescales when examining production of hot atoms. If the bulk of the hot atoms are produced over altitudes where collisions are dominant (or TCOU is small), we would expect relatively few hot atoms to escape full or partial thermalization. This would result in a much reduced coronal temperature and extent. 2.6 Summary The formation of hot oxygen coronae in the atmospheres of Venus and Mars via the nonthermal process of dissociative recombination of O f was examined using both hydrodynamic and kinetic theory approaches. Both methods predicted an extended hot oxygen corona. The constant temperature hydrodynamic model predicted hot oxygen densities in the order of 3xl0 3 -4xl0 5 c m - 3 at an altitude of about 1500 km on Mars for hot product temperatures between 1000-8500 K . For a similar altitude on Venus, the predicted hot oxygen densities are of the order of 10-2xl0 6 c m - 3 . Variation of the temperature, especially below 5000 K , produced sizable differences in the predicted hot atom density profile, and showed how sensitive the final density profile was to the choice of hot oxygen temperature. A time dependent Boltzmann equation was used to include the effect of thermaliza-tion of the hot oxygen atoms through collisions with the cold background and include Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 81 the hot oxygen temperature in a consistent manner. The time dependent hot oxygen velocity distribution function was then calculated, with energy degradation described by the Boltzmann collision operator. A simple linear transformation was used to derive the hot atom distribution function as a function of altitude. The distribution at the exobase level was fitted using hot and cold Maxwellian distributions, which yielded exobase tem-perature and densities. These parameters were used in the collisionless theory of Jeans and Chamberlain to calculate hot atom densities at altitudes above the exobase. The present calculations yielded hot oxygen densities of the order of 102 at an altitude of about 1500 km on both Mars and Venus. The hot product temperature for both planets was 7 500 K . The predicted hot atom densities calculated in this way differ from Monte Carlo derived values by a factor of approximately 2-10 on Mars, and 2-30 on Venus, al-though the predicted hot oxygen desnsities are still well above measured thermal oxygen densities for both planets. The current and previous models use common input profiles, collision cross sections, source distributions, exobase altitudes, and a collisionless model to generate exospheric densities, so these factors cannot account for the discrepancy. It should be noted, however, that there are differences of up to an order of magnitude be-tween the results of previous workers due to small changes in model parameters, and so the magnitude of the deviation of the current results are not unreasonable. The C S D A , which requires a roughly constant average velocity in order to legitimize the use of a linear transformation, is less applicable on Venus than on Mars, where there are higher background densities and hence higher rates of thermalization of the hot oxygen atoms. This is particularly true for our model, which fixes the background density at the source altitude despite the fact that we transform from time to space in our local Boltzmann equation. The Monte Carlo simulations, which follow particle motions in altitude directly, account for the change in background density, and so would be expected to yield lower levels of thermalization (that is, a higher population of energetic oxygen Chapter 2. Diffusion and Boltzmann Equation Models of Hot Oxygen Coronae 82 atoms). This is confirmed for both Mars and Venus through comparison with results of previous workers. Chapter 3 Nonthermal Production of Energetic Hydrogen and Deuterium 3.1 Introduction The present state of our understanding of the distribution of hydrogen in planetary exospheres is determined predominantly by measurements of the emissions of exospheric hydrogen. These include Lyman-a and Lyman-/? emissions of atomic hydrogen at 121.6 nm and 102.6 nm, respectively. For the terrestrial exosphere, hydrogen densities have typically been inferred by fit-ting observed emission rates to calculated intensities based on a spherically symmetric exospheric density distribution. 6 6 Information about hydrogen densities near the exobase have been provided by measurements of Lyman-a by low-altitude satellites such as OSO-4, O G O - 4 , 1 1 3 ' 1 1 4 and OGO-6 . 1 1 5 There have also been ground based resonant fluorescent measurements of the geocoronal Balmer-a spectral line at 656.3 n m , 1 1 6 although there has been some difficulty in filtering out contamination in the measured emissions from multiple scattering. 1 1 7 He et a l . 1 1 8 have used a radiative transfer model to attempt to reconcile theoretical models of the thermospheric and exospheric hydrogen density with these data. In all cases, it has been determined that an extended population of hot hydrogen exists about the Earth. The first measurements of Lyman-o; emissions which indicated the presence of a hot hydrogen corona around Venus were made by the ultraviolet photometer on Mariner 5 1 3 83 Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 84 during the flyby of 1967. The data indicated two, distinct, exospheric scale heights. Ini-tially it was thought that the scale heights were indicative of two different species, such as atomic and molecular hydrogen or atomic hydrogen and deuterium. This explanation was reconsidered by several workers 2 , 5 who suggested instead that the Mariner 5 data were more consistent with a two-temperature hydrogen exosphere. Photodissociation of molecular hydrogen was suggested as the source of the hot hydrogen atoms. Analysis of data collected by the ultraviolet spectrometer on Mariner l f j 8 7 , 8 8 during its Venus flyby in 1974, in addition to data obtained by the Lyman-a photometers on the Venera 9 and 10 orbiters 3 in 1975 and Venera 11 and 12 8 9 in 1978, supported the theory of a two temper-ature hydrogen exosphere. High quality data available from the ultraviolet spectrometer aboard the Pioneer Venus Orbiter, 4 which reached Venus in 1978, confirmed the result. It appears that more likely sources of the hot hydrogen are the nonthermal interactions of ionospheric protons with exospheric oxygen and hydrogen or elastic collisions of hydrogen with hot oxygen. 3 2 Although data is less extensive for Mars, Lyman-a emissions similar to those mea-sured on Venus have been recorded in the Martian exosphere by ultraviolet spectrometers aboard Mariner 6 and Mariner 7 1 1 9 - 1 2 0 in 1969 and aboard Mariner 9 1 2 1 in 1971. Calcu-lations of the hot hydrogen distribution about Mars indicate that the hydrogen corona is much less extensive than that of Venus, 3 2 about a factor of 100 times smaller. While such calculations are limited by the lack of available information on ionospheric densities and temperatures involved in the production of the hot hydrogen, the lack of a strong Lyman-a signature for hot hydrogen about Mars is consistent with the relatively large thermal hydrogen population inferred from available d a t a . 1 2 0 ' 9 5 These observations of the exosphere together with insitu mass spectrometric measure-ments provide density and temperature profiles of neutral and charged constituents. For example, the mass spectrometer on the Pioneer Venus Orbiter measured an enrichment Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 85 of the deuterium to hydrogen ( D / H ) ratio in the exosphere of Venus by a factor of 100 relative to the terrestrial value. This enrichment of deuterium relative to hydrogen is believed to arise from the enhanced escape of hydrogen due to nonthermal processes. Nonthermal processes refer to collisions between exospheric species and translationally energetic species (both ions and electrons), generally of ionospheric or plasmaspheric origin. This includes processes such as the collision of hot plasmaspheric protons with exospheric hydrogen, H+* + H -> H* + H+ (3.1.1) which effectively converts the energetic protons to translationally excited hydrogen atoms, H * , with speeds in excess of the escape speed. This nonthermal process is of great impor-tance in reconciling discrepancies between predictions of hydrogen escape fluxes utilizing the thermal Jeans' flux and observations. Models of the chemistry and transport pro-cesses of the terrestrial mesosphere Liu and D o n a h u e , 3 8 - 4 0 Hunten and Strobel 4 1 and Maher and Tins ley 7 7 demonstrated that the flux of hydrogen in all forms is equal to approximately 1.5-1.8xl0 8 c m - 2 s - 1 , independent of the exospheric temperature. How-ever, the Jeans escape flux calculations predict that the hydrogen escape flux should increase with an increase in the hydrogen temperature at the exobase, Tc. The nonther-mal charge exchange process given by equation (3.1.1) was first suggested by C o l e 1 0 in order to reconcile this discrepancy. It has since received a great deal of attention by many o t h e r s . 1 2 2 ' 4 2 ' 4 3 ' 1 2 3 ' 2 6 ' 4 4 ' 5 ' 4 5 ' 7 8 The energetic protons for this process are produced by photodissociation of hydrogenous compounds in the ionosphere. It has been verified by several workers 4 3 ' 4 5 ' 6 6 that the charge exchange induced escape flux decreases with the hydrogen exobase temperature, Tc, and that the sum of the nonthermal and Jeans' fluxes is constant, independent of exospheric temperature and consistent with the meso-spheric models. This has been demonstrated very clearly by Shizgal and Lindenfeld, 4 5 Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 86 who employed a simple hard sphere collisional model to obtain an analytic expression of the charge exchange induced flux. They showed that the charge exchange induced escape, Fce, is given by Fce = nH+yj2kTc/imiHe-x<[{l + r) - V l T ^ ' l / V (3.1.2) = 2kcenH+nH where the important parameter is r = TC/TH+ — 1, and nn+ is an effective density. For sufficiently large TH+/TC temperature ratios, equation (3.1.2) gives a rate coefficient of the form kce = 3.6xl(n 6s-7(l - T C / T H + ) (3.1.3) For temperature values of T#+ = 4000 K and TC — 1000 K , this gives a value of 4 . 8 x l 0 _ 6 s _ 1 close to the estimate obtained by others. 4 2 Despite the approximate na-ture of equation 3.1.2, this result has still been found very useful in the interpretation of ground based Balmer-a observations of geocoronal h y d r o g e n . 1 2 4 , 1 1 8 ' 6 6 It is of considerable interest to extend the work of Shizgal and Lindenfeld to more realistic charge exchange cross sections. In the Venusian exosphere, the charge exchange process is essentially the dominant escape mechanism since the thermal escape is very slow because of the low exospheric temperature (see equation (1.2.9) and Table 1.2). Some important nonthermal charge-exchange processes of importance, in addition to reaction (3.1.1) are reactions involving deuterium ions, H+* + D -> H * + D+ (3.1.4) D+* + H -> D* + H + (3.1.5) It is believed that the enhanced escape of hydrogen due to these and other nonthermal processes may have had significant bearing on the history of possible water loss on Venus. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 87 Current interpretations of the deuterium/hydrogen ratio measured by the Pioneer Venus Orbiter appear indicative of the relative importance of nonthermal loss of hydrogen on Venus in the pas t . 5 8 ' 5 1 ' 5 2 A difficulty with interpretations suggesting the loss of an 'Earth equivalent ocean' of water is that mechanisms yielding enhanced escape of hydrogen do not wholly account for removal of water; 3 7 oxygen must also be removed at a rate con-sistent with the stoichiometry of water. This is particularly difficult on Venus where the escape energy of such massive constituents is high (see Table 1.1). Mechanisms includ-ing incorporation of excess oxygen into crustal rocks or early m a g m a s , 1 2 5 ' 3 7 ' 1 2 6 enhanced nonthermal escape of hot oxygen, 1 2 7 and tectonic resurfacing 5 8 have been suggested, al-though at this time it is still unclear whether any or all of these processes could have been sufficiently vigorous over geologic timescales to absorb the required excess oxygen. Another important process in the energization of hydrogen is the 'resonant' charge exchange process between oxygen and hydrogen, H + * + 0 _> H * + 0 + (3.1.6) and the equivalent process for deuterium, 0 + * + D —> 0 + + D* (3.1.7) Both reaction (3.1.6) and (3.1.7) have been considered as important mechanisms for the enhancement of hydrogen escape. 3 5 Reaction (3.1.6) is able to proceed rapidly in either direction with the overall direction of reaction depending strongly on the relative rate of production of H + and 0 + at the altitude of interest. It is also possible to produce hot hydrogen and deuterium through momentum transfer collisions with hot exospheric oxygen. Such energization of atomic hydrogen via elastic collisions with translationally energetic oxygen is of the form 0* + H - » 0 + H * (3.1.8) Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 88 where the hot oxygen may or may not become thermalized completely in the collision. This process also acts as a sink for hot exospheric oxygen. The analogous process for deuterium would be of the form 0* + D 0 + D* (3.1.9) This energy transfer process was suggested in the early 70's by B r i n k m a n n . 1 2 8 McElroy et a l . 3 5 employed this mechanism to demonstrate that H would readily escape from Venus whereas D would not, thus providing the enhancement of D over H consistent with the Pioneer Venus mass spectrometer measurements. Cooper et a l . 7 5 employed the quantum mechanical differential elastic cross sections for O - H collisions to calculate the energy transfer and the velocity distributions of the product hot hydrogen. They also determined the escape fraction of hydrogen by this process and demonstrated that the angular distribution of the differential cross section is important. They showed that the fraction of product H atoms above the escape speed is 8.5% with the quantum cross sections and 18.8% if isotropic scattering is assumed. The value of this fraction calculated by McElroy et a l . 3 5 for isotropic scattering is 15%. Hodges 1 2 9 recently recalculated the differential cross sections and found some errors in the scattering calculations of Cooper et al.. Hodges did not report any results for escape for either H or D . Gurwell and Y u n g 7 6 reconsidered the enhanced escape of H over D with this same mechanism. Instead of the actual quantum cross sections, they employed Henyey-Greenstein 1 3 0 analytic fits to the cross sections. They calculated product velocity distribution functions for both H and D , and calculate that the fraction of H atoms above the escape energy for isotropic scattering is 15.8%, in agreement with McElroy et al. However, they calculate that the fraction of H atoms above the escape energy with the actual differential cross sections is 4.6% in disagreement with Cooper et al.. They did their calculations with the background gas at both OK and 300K. For the cold gas Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 89 case, the product distributions can be determined analytically whereas with the thermal background, they carry out their calculations with a Monte Carlo simulation. Nonthermal processes provide an important escape mechanism for heavier species, such as oxygen, nitrogen and carbon, for which the thermal escape is negligible. A n important nonthermal process in the Martian exosphere is the dissociative recombination of O f with electrons, mentioned in Chapter 2 by equation (2.1.1); that is, 0+ + e" -* O* + 0* with the product oxygen atoms translationally excited. This process produces a large population of energetic (or 'hot') atoms, and there is evidence of such hot coronae of atomic oxygen in the Martian and terrestrial exospheres. 1 1 ' 3 1 ' 7 ' 3 2 ' 6 ' 3 3 ' 3 4 If the hot oxygen atoms do not escape and are ionized in the upper atmosphere, they may in turn excite other species via charge exchange reactions. An example of such resonant charge exchange reactions of hot ionospheric atomic oxygen with neutral atomic hydrogen and oxygen are given by 0 + * + H -» 0* + H + (3.1.10) 0+* + 0 -> 0* + 0+ (3.1.11) Reaction (3.1.10) may also act as a sink for both thermal and hot hydrogen. Dissociative recombination of Of, as given previously above, is believed to be the main mechanism for the production of a hot oxygen corona on Mars , 6 - 8 although with its high atmospheric concentration of carbon-dioxide, the dissociative recombination reaction C O f + e ^ C O + O* (3.1.12) may also be important in the formation of hot oxygen on Mars, 3 2 although its effectiveness in producing hot oxygen may be limited by the reduction of C O f via reactions such as C O f + 0 Of + CO Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 90 The representation of the dynamics of nonthermal processes is given by the collision cross sections used to approximate the physics of such processes. In many past works examining the role of charge exchange in the production of hot species , 7 7 ' 4 3 ' 5 ' 4 4 ' 4 5 the so called 'Linear Trajectory Approximation' (or L T A ) has been used. This approximation assumes that exospheric charge exchange takes place with no momentum transfer, in effect simply exchanging the incoming velocities of the two collision partners to arrive at their post-collisional values . 1 3 1 ' 7 3 The difference between the actual charge-exchange process and the L T A are illustrated schematically in Figure 3.21. While the L T A of the charge exchange cross section has made detailed (and even analytic) analysis of this process possible, it has been shown that the approximation appears to overestimate the effect of charge exchange on the exospheric escape flux of the hot product . 1 3 1 More accurate and realistic representation of the collision dynamics between ions and atoms in the charge exchange process is of great interest, and has been examined in some detail for certain systems, such as H - H + and D - H + . 7 1 ' 7 3 ' 7 4 ' 7 2 In addition, the velocity distribution function ( V D F ) of the product hot atom, and hence the escape flux, is sensitive to the form of the collision cross section, it is vital that it accurately reflect the physics of the process, and that these cross sections be included in models describing the escape process. Given the importance of nonthermal collisional processes in the exosphere, it is clear that the standard collisionless models have to be reconsidered. For this reason, our goal is to calculate the non-equilibrium and nonthermal V D F of the product hot atoms resulting from nonthermal collisional processes. A consistent description of these product velocity distribution functions ( P V D F s ) is required to, ultimately, accurately estimate escape fluxes resulting from nonthermal mechanisms. This is done in Chapter 4. The current chapter is concerned with the calculation of the P V D F s for several atomic Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 91 u | (b) Figure 3.21: Schematic of the charge-exchange and LTA. The charge-exchange collisional process is illustrated in (a), where an electron is transferred from the neutral to the ion, and momentum is transferred from the ion to the neutral. The LTA is shown in (b), where again an electron is transferred from the neutral to the ion, but no momentum is transferred between the ion and neutral. The result is that the collision appears simply to have switched the initial velocities of the ion and neutral, u = v' and v = u', with no change in direction. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 92 systems and nonthermal processes. In section 3.4 we investigate a kinetic theory descrip-tion of the production of nonequilibrium, hot products resulting from nonthermal pro-cesses. To calculate these P V D F s , we require an accurate representation of the energy and angle dependence of various nonthermal collision processes. As mentioned previously, this entails the calculation of quantum mechanical differential and total collision cross sections. We detail the calculation and form of realistic quantum mechanical collision cross sections in section 3.2. The systems examined include direct-plus-exchange ( D P E ) elastic collisions, direct elastic (DIR) collisions, charge exchange (CE) and L T A collisions for H - H + and D - H + , and direct elastic collisions for 0 - H and 0 -D. The resulting P V D F s are used to calculate energy exchange rate coefficients, which are compared to calcula-tions using a simple hydrodynamic approach. The energy exchange rate coefficients are used to examine the time evolution of the average energy for hot products resulting from nonthermal processes as well as providing a check on our kinetic theory approach of the same calculation. The P V D F s are also used to calculate the escaping fractions of the hot products, and compare these results with those of previous w o r k e r s . 1 3 1 ' 7 1 ' 7 3 ' 7 4 ' 7 6 3.2 Collision Cross Sections As mentioned in the introduction to this chapter, early work in the field approximated the charge exchange mechanism by considering the transfer of momentum during the charge exchange collision to be minimal (i.e. the L T A ) . If the pre- and post- collisional velocities of the ion and atom are indicated by the pairs (c-,c^ ) and (c,-,ca), this is equivalent to the simple interchange of particle labels, that is c'- = ca and = c,-. Alternately, this is also equivalent to a charge exchange cross section which is a delta function in the backward scattering direction, 0 = ir, and that there is no contribution from elastic scattering during the collision. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 93 In order to more accurately represent the details of the kinematics of the collision process, it is necessary to calculate the collision cross section using quantum mechanical methods. The importance of the angle and energy dependence of real, quantum cross sections was noted by Cooper et a l . 7 5 in their study of O H . This has been demonstrated by S h i z g a l 1 3 1 ' 7 1 for hydrogen-proton and hydrogen-deuterium systems and was recently re-confirmed by Hodges and B r e i g 7 3 ' 7 4 and Hodges. 1 2 9 3.2.1 Quantum Mechanical Scattering The procedure for the calculation of quantum mechanical (QM) differential cross sections is well known. A detailed description of these calculations is given in Appendix D . l , standard textbooks on quantum mechanics, 1 3 2 " 1 3 4 and in the l i t e ra ture ; 1 3 5 ' 1 3 6 ' 7 3 ' 7 2 here we simply summarize the important results. The collision of a pair of particles of masses m i , m 2 may be reduced to the motion of a single particle of reduced mass p = raira2/(rai + m 2 ) moving in a potential due to the scattering center, given by V(r). In quantum mechanics, this scattering problem involves the solution of the Schrodinger wave equation. The amplitude of the scattered wave, denoted by J-', may be related to the differential collision cross section by (see Appendix D.2) - ^ ( 2 / + l)P,(cosr;)e^ sin<5, (3.2.1) and so, from equation ( D . l . 14), the differential elastic cross section is _1_ oo 2 o(E,0) = y^(2/ + l)P,(cos0)e i 5' sin5, (3.2.2) From equation (CO.5), the total elastic cross section is given by oUE) = — £ ( 2 / + l ) s i n 2 £ , (3.2.3) Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 94 Equations (3.2.2) and (3.2.3) above assume only single channel scattering (that is, a single interaction potential describes the particle interaction). For most real systems of atoms and molecules, there are several elastic scattering channels for any given energy. In such cases, the elastic differential scattering cross sections are given by combining the contri-butions of each channel according to their quantum mechanical statistical we ight . 1 3 7 ' 1 2 9 The momentum transfer cross section is defined as a weighted integral of the differ-ential cross section over all angles, or amt(E) = J a(E, 6) (1-cos 6) dSl (3.2.4) It appears in the final form of the energy exchange rate coefficient, discussed later in this chapter. In addition, the momentum transfer collision cross section is often used in aeronomical calculations. 1 3 7 The linear trajectory approximation (or L T A ) cross section is often used to approxi-mate the charge-exchange process. As was mentioned briefly in the introduction to this section, this model cross section is the equivalent of a charge exchange cross section which is a delta function in the backward scattering direction, 9 = TT, and assumes that there is no contribution from elastic scattering during the collision. The definition for the L T A cross section i s 7 2 , l t M = & « L ± ^ > ( , 2 . 5 ) 2ir sm(7r — e) where e is a positive number between 0 and 7r, and o^s j s s o m e (arbitrary, energy independent) total hard sphere cross section; the L T A approximation is given in the limit e —> 0+. Another definition of the L T A , which incorporates energy dependence, is given b y 7 3 1.0 0.5 0.0 -0.5 J I L J I I I L J I L 0 1 2 3 4 5 6 7 8 9 10 Internuclear separation, r (aQ) Figure 3.22: Interaction potentials for H+-H. The solid curve is for the gerade or lSa state, and the dashed curve is for the ungerade or 2Pau state. The unit of energy for V(r) is the ionization energy of atomic hydrogen, while the unit of internuclear separation is the Bohr radius, a0 = 0.5292 cm. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 98 e (radians) e ( radians) (a) (b) Figure 3.23: Differential cross sections for H + - H . The differential scattering cross sections, in square Angstroms, are shown for fixed energy as a function of the scattering angle, 9. In (a), the energy is 0.01 eV; the top part of the figure shows the differential cross sections for direct (solid line) and charge exchange (dashed line) cross sections, and the bottom part of the figure shows the direct plus exchange ( D P E ) differential elastic cross section. In (b), the energy is 1.0 eV; the top part of the figure shows the differential cross sections for direct elastic scatter, the middle portion shows the differential elastic cross section for the charge exchange cross section, and the bottom figure shows the direct plus exchange ( D P E ) differential elastic cross section. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 99 6 ( r a d i a n s ) 6 ( r a d i a n s ) (a) (b) Figure 3.24: Differential cross sections for D + - H . The differential scattering cross sections, in square Angstroms, are shown for fixed energy as a function of the scattering angle, 8. In (a), the energy is 0.01 eV; the top part of the figure shows the differential cross sections for direct (solid line) and charge exchange (dashed line) cross sections, and the bottom part of the figure shows the direct plus exchange ( D P E ) differential elastic cross section. In (b), the energy is 1.0 eV; the top part of the figure shows the differential cross sections for direct elastic scatter, the middle portion shows the differential elastic cross section for the charge exchange cross section, and the bottom figure shows the direct plus exchange ( D P E ) differential elastic cross section. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 100 10000 3 1000 F 10000 1.0E-4 1.0E-3 1.0E-2 1.0E-1 E n e r g y ( e V ) 1.0E+0 1.0E+1 1.0E-4 1.0E-3 1.0E-2 1.0E-1 E n e r g y ( e V ) 1.0E+0 1.0E+1 (a) (b) Figure 3.25: Total elastic collision cross sections. Figure (a) is for H + - H , Figure (b) is for H + - D . In both (a) and (b), A denotes total cross sections for direct plus exchange ( D P E ) and B denotes total cross sections for charge exchange ( C E ) . The solid lines are for the standard total elastic cross section and the dashed lines are for the momentum transfer cross sections for each respective type of cross section ( D P E or C E ) . Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 101 = ^E(2;+1){SIN2^.+ + sin 2 <$,_} (3.2.10) k 1=0 and crcteot(E) = 2TT ['ace(E,6) sm 6 d6 Jo C O = I a E ( 2 f + l ) s m a ( * I l + - h-) (3.2.11) K 1=0 respectively, where again + or — refers to a value calculated using either the gerade or ungerade potential. Some examples of differential cross sections for hydrogen-proton scattering calculated as described here are illustrated in Figure 3.23, with total cross sections shown in Figure 3.25. We assume that the potentials of Figure 3.22 and equations (3.2.8)-(3.2.11) are ap-plicable to deuterium-proton scattering when the reduced mass for the D + - H system is used in the calculation of the phase shifts. This is equivalent to ignoring the 'coupling coefficient' between the lowest two electronic states of D + - H , and appears justified given that the coupling potential is four orders of magnitude smaller than the other interaction potentials. 7 4 Differential cross sections for deuterium-proton scattering are illustrated in Figure 3.24 and total cross sections are shown in Figure 3.25. 3.2.4 Cross Sections for 0 - H and O-D As discussed previously, an important nonthermal mechanism on Mars and Venus is the production of translationally energetic oxygen atoms from dissociative recombination and the subsequent energy transfer to H and D , equations (3.1.8) and (3.1.9). The present work considers the use of quantum mechanical collision cross sections for 0 - H and O-D and the determination of the enhanced product velocity distributions of hot H and D . The collision cross sections were evaluated as detailed previously in this chapter, with the primary difference being that in this case there are four interactions Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 102 Figure 3.26: Interaction potentials for 0 - H . The labels A , B , C , and D denote the curves for the 4n, 2 E , 4 E , and 2 n states, respectively. The unit of energy for V(r ) is the ionization energy of atomic hydrogen, while the unit of internuclear separation is the Bohr radius, a0 = 0.5292 cm. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 103 potentials, or four scattering channels, corresponding to the electronic states 4II, 2 S ~ , 4 E ~ and 2II. These potential energy curves were taken from Figure 1 of Hodges , 1 2 9 and are shown in Figure 3.26. We assume that deuterium-oxygen collisions are governed by the hydrogen-oxygen interaction potentials, with the only change being the reduced mass of the system. The potentials for the 2 n and 2 E - are also listed in Table 1 of van Dishoeck et a l . 1 4 0 , 1 4 1 The bound state potential 2 n was fitted to a Morse potential of the form, V(r) = De[l -exp(-Bx)}2 where x = r — re and 8 = 30(l + Xix + A 2 x ) . 1 4 2 ' 1 4 3 We chose the values De = 4.582, r e = 1.838, B0 = 1.239, Ai = 0.001 and A2 = 0.037. The other three potentials were digitized from the work of Hodges, 1 2 9 and a spline fit was used to interpolate in the table. A n important part of the potentials is their asymptotic forms for large r. These are of the form 9.13/r 6 (in au) for the S states and 9.22/r 6 for the n state, 7 5 and represent the long range interaction of the two polarizable atoms as dipoles. The way in which the tabulated potential functions are made to go over to these long range forms is somewhat arbitrary. It is not known at what radial positions these asymptotic forms become valid and it appears that the tabulated potentials do not extend into this region. Consequently a rather arbitrary merging of the tabulated potentials with the asymptotic portion is made. Hodges alludes to this difficulty in attempting to reconcile some differences between his cross sections and those reported by Cooper et al. The actual forms of the potentials in the intermediate range (6-10 au) will have to await further quantum mechanical calculations. At short range, the potentials are fitted with polynomials of the form a/r + b + cr, as in Hodges. 1 2 9 These potential functions were used to calculate four sets of phase shifts and four sets of cross sections. The cross section for O - H and O-D collisions are the sum of these cross Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 104 Figure 3.27: Total elastic collision cross sections. Figure (a) is for O - H , Figure (b) is for O - D . The labels A , B , C , and D denote the total cross sections corresponding to the 4 n , 2 S , 4 S , and 2 n states, respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 105 Figure 3.28: Total elastic collision cross sections. The labels A and B for the solid curves denote the total elastic cross sections for O-D and 0 - H , respectively, as given by the weighted average of the total cross sections for each of the 4 states of O - D / O - H . The dashed curve indicates the analogous elastic momentum transfer cross section; at this scale, it is identical for 0 - H and O-D. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 106 100000 i 1 1 r I I I i I I I . I , I , I 0I i I i I , I , I , I , I 0 30 60 90 120 150 180 0 30 60 90 120 150 180 8 (degrees) 8 (degrees) (a) (b) Figure 3.29: Differential cross sections for (a) 0 -H and (b) O-D. The differential scatter-ing cross sections, in square Angstroms, are shown for fixed energy as a function of the scattering angle, 6. The differential cross sections are shown for an energy of 1.0 eV, are are a weighted sum of the four states for each system. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 107 sections weighted by the statistical weights, W{, determined from the degeneracies of the electronic states, that is, o(E) = J2WME) i where i denotes one of the four channels. The energy variation of the total cross sections for O - H and O-D collisions are shown in Figure 3.27. While similar, the results are slightly different for the two systems. The oscillatory curve is for the 2 n state, while the smooth curves correspond to the 4 n , 2 E ~ , and 4 E ~ states. Hodges found that the cross sections for the 2 E ~ and 4 E ~ states (the two lowest curves in the figure) are indistinguishable, whereas here we find a slight difference. Cooper et a l . 7 5 reported a significant difference between these two cross sections. It would appear that the manner in which the potentials at small r are merged to the asymptotic form is responsible for slight differences in the cross section but probably not the large difference noted in the results of Cooper et al., as discussed by Hodges. We agree with his conclusion that the large difference that they report is unreasonable. The statistically weighted cross sections for H-0 and D-0 collisions are shown in Figure 3.28, along with the momentum transfer cross section for elastic scatter. There is a small but real difference between the weighted total elastic cross sections, in contrast with Hodges' supposition that these cross sections are the same. This is important since the object of these calculations is to study the enhancement of escape of H over D. Differential cross sections for both types of collisions at 1.0 eV are illustrated in Figure 3.29. The large forward scattering peak is evident in both cross sections. 3.3 Energy Exchange Rate Coefficients A n aspect of interest in the production of hot atoms involved in the nonthermal escape process is the rate at which energy is transferred between various atmospheric Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 108 constituents. Calculation of the energy exchange rate coefficients enable us to compare the relative efficiency with which various nonthermal processes operate. We may also use them to study the time evolution of the average energy of between species of differing (initial) temperatures. Energy exchange rate coefficients are reported by Hodges and Bre ig , 7 3 who examine charge exchange and momentum transfer effects in hydrogen-proton collisions. They carry out quantum mechanical calculations of hydrogen-proton cross sections for both processes. However, they choose to use Monte Carlo integration techniques to treat the Boltzmann collision operator, citing complexities in the relationships of pre- and post-collisional encounter velocities produced by wide-ranging variations in the collision cross section. We show in this section that such Monte Carlo techniques are not required for the calculation of the energy exchange rate coefficients. We follow the method and notation of Shizgal and F i t z p a t r i c k 1 4 4 - 1 4 6 and Shizgal . 1 4 7 The derivation of the energy exchange rate coefficient involves the calculation of the energy moment of the Boltzmann equation. This yields the rate of change of the average energy of particles of species 1. In general, the energy exchange rate coefficient depends in a complicated way on the details of the differential collision cross section and the form of the velocity distribution functions of both species involved in the collision process. This involves solving an equation containing the Boltzmann collision operator, as detailed in equation (1.4.2). It was demonstrated by Shizgal and F i t z p a t r i c k 1 4 4 ' 1 4 5 that it is possible with some simplifying assumptions to reduce the problem of the calculation of the energy exchange rate coefficients to a single integral over the total momentum transfer collision cross section. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 109 3.3.1 Theory The rate at which species 1 transfers energy to the background, species 2, for a uniform system in the absence of external forces is given by the energy moment of the collision operator in the Boltzmann equation, that is E =~dT = J 2m^-mdCl = ni(r)n2(r) JJJ {/i(r,ci)/2(r,c'2) - /i(r,ci)/2(r,c2)} x -micl ga(g1fl) dfl dc2 dci (3.3.1) where the Boltzmann collision term was given by equation (1.4.2). We assume that both velocity distribution functions are Maxwellian distributions at different temperatures, f M r M ( m i V2 f i ( C i ) = \2MJ 6 X P m,c;-/2 2W If primed/unprimed velocities in the collision integral in equation (3.3.1) are interchanged, we have that / / /1/2 ^ m i c i V dto dc2 dci = j j / i / 2 ^rriic'fgo- dtt dc2 dci where 2 9 g'(r(g', ft') dfl' dc[ dc'2 = gcr(g, ft) dfl dc\ dc2 and where we have also used microscopic reversibility; that is, we assume that every forward collision has a corresponding inverse collision. Since summing over all possible collisions is equivalent to summing over all possible inverse collisions, the value of the integration is unchanged by the interchange of primed and unprimed variables. Equation (3.3.1) can thus be written RE = nxn2 JJJ / f / f ^ m a (cf - c 2) gadVldc2dCl (3.3.2) Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 110 W i t h the transformation to dimensionless velocity coordinates, * = (0*C (3.3.3) equation (3.3.2) becomes ^ = j]J 9°Md*^ (3.3.4) Hodges and B r e i g 7 3 evaluate equation (3.3.4), an eight dimensional integral, using Monte Carlo integration. They cite, as justification for resorting to Monte Carlo techniques, that direct integration of equation (3.3.4) is 'impractical'. However, the techniques for the reduction of equation (3.3.4) to a one dimensional integral over a total cross section are well k n o w n , 1 4 4 - 1 4 6 and are outlined below. The center of mass velocity and relative velo'city are given, respectively, b y 2 9 g = c2 - Ci G = M1c1 + M2c2 (3.3.5) where Mi = m , / m 0 m0 = m x + m 2 We have standard relations 2 9 between particle velocities and the center of mass and relative velocities, namely ci = G - M 2 g c2 = G + M l g c't =' G ' - M 2 g ' c/2 = G' + M l g ' (3.3.6) Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 111 We may define a reduced center of mass velocity, S, and relative velocity £ , £ = MU2-MU, (3.3.7) where MXT2 M i T 2 + M2Tx M2TX MXT2 + M2TX We may relate S and £ to G and g by 144 (mxl2kTx + m2/2kT2)^ [G + M 1 M 2 ( J 1 1 - T 2 ) / ( M 1 T 2 + M2Ti)*\ 2kT( S (3.3.8) where p = m i r a 2 / m 0 is the normal reduced mass, and T e / / = (M\T2 + M2Tx)jm0 is the effective temperature. Inverting equation (3.3.7), we have i2 = MIE + MU (3.3.9) We note that equations (3.3.7) and (3.3.8) are valid for the pre-collisional (primed) ve-locities by simply replacing variables with their primed counterparts. For elastic collisions, we also have 2 9 g g - 2 ( g - k ) k (3.3.10) where the unit vector k is the 'apse line' or line joining the two particles at the point of closest approach during a collision; it corresponds to the external bisector of the angle between the relative velocities before and after the collision, g' and g (see Figure E.59). Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 112 Equation (3.3.10) is a result of the principles of conservation of momentum and energy. From equations (3.3.6) and (3.3.10) we have c'x = G - M 2 g + 2 M 2 ( g - k ) k (3.3.11) where by conservation of momentum G ' = G . Using the primed analogue to equation (3.3.9), and utilizing equation (3.3.11), we then have & = MfE-M%{l + a)t' + M%at (3.3.12) where a = Mi(Ti/T2 — 1). After some algebra (see Appendix E), we find g - & = 4 ( M 1 ^ ) l ( l + « ) ( ^ k ) ( S - k ) + 4 M 2 « ( l + a ) ( $ . k ) 2 (3.3.13) From conservation of energy (see Appendix E), we have that d^di2 = dEd£ (3.3.14) and so equation (3.3.4) becomes RE = ^kT^JJJ [ 4 ( A * i A * 2 ) * ( l + « ) ( £ • k ) ( S - k ) + 4A4 2 a ( l + « ) (£ • k) 2] x e-~-? (^J11^2 (adCldSdZ (3.3.15) The first term in equation (3.3.15) is odd in the cosine of the angle between k and S, and so vanishes when integrated over dS; that is, if we choose our coordinate system such that the polar direction is along 2, then dS = sinO d9 d(f> E2 dE, and so k • S = kE cos 9, which means that / kE cos 9d(cos 9) = 0. Thus, equation (3.3.15) becomes RE = ^kT^M2a(l + a) ^f^j J e~* dE> J J (£-kfe^H^ d£l df-= ^ k T l AM2a(l + a) ' 4 . jT e~* E2 dE x JJ ( £ - k ) V * 2 iadVld^ Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 113 0 0 —2 3 47T / e~~ E2 dZ = 71-2 The integral over dZ may be done analytically. i and we have RE = ^ ^ 4 ^ 2 0 ( 1 + OL) {^y^j 4TT yy {lcozx)2e-?eodSldi where we have defined x a s the angle between k and £ . Wi th the transformation (see Figure E.59) cos 2 x = ^ ( l - c o s 0 ) and the definition of the momentum transfer cross section, equation (3.2.4), we thus have /•oo , RE = N{TUT2) e-U6. T 2 , AE* < 0, and energy is transferred from species 1 to species 2, the background, and vice versa. As well, we note that when T i = T 2 , AE* = 0, and there is no further energy exchange (equilibrium). Thus, depending on the initial temperatures (energies) of the interacting species, we can have either relaxation (cooling) or excitation (heating) of species 1. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 116 3.3.3 Results and Discussion The calculation of the energy exchange rate coefficient, using equation (3.3.16), involves a single integral over the reduced relative velocity £. From a computational standpoint, this one-dimensional integral is much more efficient and accurate than than the original triple integral of equation (3.3.1). The momentum transfer cross section(s) required for the £ integration are calculated at 91 energies, from 1 x 1 0 - 4 to 10 eV, and at 361 angles (half degree increments), and are illustrated in Figure 3.25. The integral over d^ is done using a transformation to reduced energy, z = £ 2 , or f°° 2 1 f°° / e-^ £5 omt d£ = -I e~ V omt dz Jo 2 Jo and applying an interpolation scheme with z points located at the 91 energies used when calculating the differential cross section (the interpolation scheme is detailed in Appendix E) . We compute values of in units of c m 3 s _ 1 , for several momentum transfer cross sections, corresponding to differential cross sections for elastic direct (DIR), charge ex-change ( C E ) , and the linear trajectory approximation (LTA) for the charge-exchange cross section. The temperatures of the ion-neutral pairs are chosen to approximate the observed or estimated range of temperatures for exospheric conditions on the terrestrial planets (Venus, Earth and M a r s ) . 5 , 1 4 8 For hydrogen-proton collisions, Table 3.5 summarizes the results for the DIR, C E , and L T A cross sections for several combinations of hydrogen and proton temperatures. We note that the L T A cross section overestimates the energy exchange rate coefficient for charge exchange by approximately 30% for all combinations of proton and hydrogen temperatures. We observe a uniformly monotonic increase in the energy exchange rate coefficient for all cross sections as the hydrogen temperature is increased for a fixed proton temperature except for the direct (DIR) elastic interaction, which remains roughly Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 117 Cross T[H] T[H+] Section 500 1000 1500 6000 DIR 2.30 2.29 2.28 C E 8.87 9.19 9.51 L T A 11.12 11.43 11.74 4000 DIR 2.27 2.29 2.30 C E 7.47 7.84 8.19 L T A 9.67 10.06 10.43 2000 DIR 2.07 2.14 2.19 C E 5.82 6.27 6.69 L T A 7.80 8.33 8.81 Table 3.5: Energy exchange rate coefficients for several H+-H cross sections. The tem-peratures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of the energy exchange rate coefficients are in units of 1 0 - 9 c m 3 s- . Cross T[H] T[H+] Section 500 1000 1500 6000 DIR L T A 20.78 22.07 20.48 21.27 19.90 20.99 4000 DIR L T A 17.67 18.58 17.77 18.16 17.50 18.08 2000 DIR L T A 13.23 13.71 13.70 14.02 16.02 16.72 Table 3.6: Energy exchange rate coefficients for several H + - H cross sections as reported by Hodges and Bre ig . 7 3 The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of the energy exchange rate coefficients have been converted to units of 1 0 - 9 c m 3 s - 1 for comparison with Tables 3.5 and 3.12. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 118 constant for all temperatures considered. The results for hydrogen-proton collisions may also be compared with the work of Hodges and Bre ig , 7 3 which is summarized in Table 3.6. Comparing tabulated values for the L T A we note that there is a discrepancy between their results and the present work, on the order of a factor of roughly 1.6 to 1.9. It should also be noted that while Hodges and Breig designate the elastic cross section as the regular 'direct' elastic collision cross section, it appears that they have mislabelled it, and in fact those results are for the D P E collision cross section. If we sum the DIR and C E rate coefficients in Table 3.5 to yield the D P E result, we get a similar discrepancy as with the L T A . The calculation of the energy exchange rate coefficients is repeated in an independent calculation later in this chapter, and the same discrepancy is observed; it is unclear as to its origin. It may be a systematic statistical problem with their Monte Carlo integration scheme, but without repeating their calculation exactly this is difficult to determine absolutely. The cross sections from both works have been compared and show excellent agreement, so that cannot be the source of the discrepancy. In addition, their L T A results exhibit strange variation as the hydrogen temperature is changed for a fixed proton temperature. For a proton temperature of 2000 K , the L T A rate coefficient appears to increase with increasing hydrogen temperature; however, for a proton temperature of 4000 K , the rate coefficient appears to increase then decrease, and for a proton temperature of 6000 K the rate coefficient decreases. The results of the current work show a uniform increase in the energy exchange rate coefficient as the hydrogen temperature is increased, regardless of the proton temperature or the cross section used, with exception of the DIR interaction. In order to determine which behaviour was correct, we fitted the momentum transfer cross section for the L T A with a power law of the form omt(E) = aEp Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 119 where E is the energy in eV, and crmt(E) is in A 2 . Several sets of parameters were assessed to fit amt(E), and are illustrated in Figure 3.30. The solid curve is the true L T A momentum transfer cross section. The short dashed, long dashed, and dot-dashed curves represent values of a and p corresponding to 83/-0.28, 80/-0.35, and 40/-0.45 respectively. These fits attempt to account for the sharp increase in the L T A cross section at low energy. A fourth fit, given by the dot-dot-dashed curve and corresponding to fit parameters of 91/-0.12, respectively, fits only the high energy portion of the cross section. W i t h this form of the momentum transfer cross section, the integral in equation (3.3.19) may be carried out analytically to yield an energy exchange rate coefficient of the form k E = *!L J ^ k l A { * W r ( p + 3) 7rm0 V n where T is the standard gamma function. The energy exchange rate coefficients as a function of the effective temperature, T e / / , are plotted in Figure 3.31. Regardless of the fit, it appears that for the L T A cross section we find that the energy exchange rate coefficient increases monotonically with an increase in the effective temperature, T e / / , verifying the pattern seen in Table 3.5. The dependence of the rate coefficients with the bath and test particle temperatures, T2 and T i , is shown in Figure 3.32. The solid curve is for the D P E cross section and the dashed curve is for the C E cross section. A (constant) value of 7\ = 1500 K is used for both curves. The values of ^ ( T i , T 2 ) are normalized by the value of & E ( T 2 = 2~i), and so the relative magnitudes of the D P E and C E curves cannot be compared directly. From the figure we see that the energy exchange rate coefficients increase almost uniformly. The disparate masses of the two collision partners makes the identification of the hot atom and background important for deuterium-proton collisions. The energy exchange rate coefficient, as given by equation (3.3.19), is explicitly dependent on the effective Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 120 Figure 3.30: Fits to the L T A momentum transfer cross section. The solid curve is the true L T A momentum transfer cross section. The short dashed, long dashed, dot-dashed and dot-dot-dashed curves represent fit parameters a and p corresponding to 83/-0.28, 80/-0.35, 40/-0.45 and 91/-0.12 respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 121 Figure 3.31: Temperature dependence of kE for the L T A . The solid curve is the true L T A momentum transfer cross section. The short dashed, long dashed and dot-dashed curves represent fit parameters a and p corresponding to 83/-0.28, 80/-0.35 and 40/-0.45 respectively. The plus symbols correspond to a fit with parameters 91/-0.12. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 122 H II H M H Figure 3.32: Dependence of energy exchange rate coefficients on bath temperature for H+-H. The value of kE{Tx,T2) is plotted relative to the value of kE(T2 = Tx) as a function of T2 for fixed values of T i . The solid curve is for the direct-plus-exchange (DPE) cross section, dashed curve for the charge exchange (CE) cross section. A constant temperature T i = 1500 K is used for both curves. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 123 temperature, T e / / through the factor \ J 2 n k ^ e f f . In addition, k g is implicitly dependent on the effective temperature through the definition of £, as mentioned in connection with the momentum transfer cross section appearing in equation (3.3.19). For hydrogen-proton collisions, we had merely to consider a single process for the formation of hot hydrogen, as given by equation (3.1.1), H+* + H D P E ' C E , ' L T A H* + H + However, to produce hot deuterium, we have from processes given in equations (3.1.4) and (3.1.5) that H+* + D D P ^ I R H+ + D* C E ^ L T A H » + D + ^ D+* + H D P J ^ m D+ + H* C E ^ L X A D * + H + ( H ) where the labels D P E , DIR, C E , L T A indicate the types of processes which produce the products on the right-hand side of the reaction. In each case, we consider the neutral reactant as species 1 (test particle) and the ion reactant as species 2 (bath). Thus, process I above can have two different sets of products, depending on the type of cross section considered; the same is true for process II. We tabulate the energy exchange rate coefficients for the above processes in Table 3.7. Columns labeled T refer to the reaction given by the processes labeled I above, while columns labeled 'IF refer to the reaction given by the processes labeled II above. Thus, a DIR reaction which is labeled 'I' refers to the first branch of process I, while a C E reaction labeled TI ' refers to the second branch of process II. It should be noted that entries for DIR 'IF and C E / L T A T are thus representative of the transfer of energy to Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 124 T [neutral] Cross 500 1000 1500 T[ion] Section I II I II I II 6000 DIR 1.741 1.76 1.731 1.78 1.72 1.78 C E 8.026 6.098 8.812 6.482 8.336 6.850 L T A 9.714 7.790 9.861 8.193 10.00 8.568 4000 DIR 1.78 1.69 1.781 1.73 1.78 1.76 C E 6.667 5.259 6.849 5.691 7.027 6.097 L T A 8.382 6.870 8.567 7.351 8.745 7.789 2000 DIR 1.663 1.55 1.691 1.63 1.72 1.69 C E 5.027 4.241 5.258 4.784 5.478 5.258 L T A 6.606 5.693 6.868 6.327 7.116 6.869 Table 3.7: Energy exchange rate coefficients for several H + - D cross sections. The tem-peratures of both deuterium and the protons, in degrees Kelvin, are indicated in the table. Values under columns marked T are for reactions where the test particle is the deuterium atom, with a bath of hot protons; values under columns marked 'IF are for reactions where the test particle is hydrogen, with a bath of hot deuterium ions. The tabulated values of the energy exchange rate coefficients are in units of 1 0 - 9 c m 3 s - 1 . produce hot hydrogen atoms, rather than hot deuterium atoms, via reactions with hot deuterium ions. Comparing values from Tables 3.5 and 3.7, we observe that energy exchange rate coefficients corresponding to production of hot deuterium via DIR collisions with protons are approximately 80-90% of the corresponding value for hydrogen-proton production of hot hydrogen. For C E / L T A collisions, the deuterium energy exchange rate coefficients are approximately 70-80% of the hydrogen-proton values. We believe this is due to the less efficient transfer of energy between particles of disparate masses. The dependence of the rate coefficients with the bath and test particle temperatures, T 2 and T i , for deuterium and hydrogen is shown in Figure 3.33. The solid curve is for the Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 125 DIR cross section and the dashed curves is for the C E cross section. A constant value of Xi = 500 K is used in the calculation of both curves. The values of kg(Ti,T2) are normalized by the value of kE(T2 = Xi) . From a comparison with Figure 3.32, we can see that the behaviour of the DIR curve is very different from that for the hydrogen-proton system. At first it was thought that the curvature of the DIR curve was due to numerical error, such as poor convergence in the integration of equation 3.3.19. To eliminate this possibility, the details of the integrands and convergence of the integration yielding the rate coefficients were examined. The results are shown in Figures 3.34 and 3.35. Since the integrands are well behaved and the convergence of the integral is smooth, we deduce that the result is in fact physical and not numeric. It may simply be a result of the change in efficiency of energy transfer for colliding partners of disparate mass. The energy exchange rate coefficients for O-H and O-D direct elastic collisions are tabulated in Table 3.8. We can see immediately that the energy exchange is much less efficient for these systems than for hydrogen-proton and deuterium-proton systems, owing to the very large difference in mass between the O-H and O-D colliding pairs. The energy exchange rate coefficients for O-D are larger than those for O-H, reflecting the larger mass of deuterium and the larger reduced mass of the O-D system. We may, as mentioned in the preceding section, calculate the time evolution of the average energy of species 1 using equation (3.3.21). We chose three temperature pairs from those used for the calculation of the energy exchange rate coefficient: TJJ+ = 6000 K/TH = 500 K , TH+ = 4000 K/TH = 1000, and TH+ = 2000 K/TH = 1500. These were chosen as they represented the full range of possible temperature ratios previously considered. In addition, we have examined both the case where, initially, T\ > T2 and where T2 > T\. In plotting the results, we have used dimensionless time units Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 126 0.8 1 1 1 1 1 1 1 1 1 i I i i i i I 0 5 10 15 T 2 / T l Figure 3.33: Dependence of energy exchange rate coefficients on bath temperature for D+-H. The value of kE(T1,T2) is plotted relative to the value of kE(T2 = Tx) as a function of T2 for a fixed value of Tx = 500 K . The solid curve is for the direct elastic (DIR) cross section, the dashed curve for the charge exchange (CE) cross section. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 127 (a) (b) Figure 3.34: Details of the integrand in calculation of the energy exchange rate coefficient for D H + . The neutral temperature, T i , is in both cases fixed at 500 K . In (a), the deuterium ion temperature is fixed at values of 500, 1875, 3250 and 4625 K , corresponding to the labels A - D , respectively. In (b), the proton temperature is fixed at 500, 2333, 4166 and 6000 K , corresponding to the labels A-D, respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 128 160 120 C/3 Figure 3.35: Convergence of the energy exchange rate coefficient for DH+. The solid curve is for the DIR cross section, the dashed curve for C E . The sum of the integration to calculate the energy exchange rate coefficient, as a function of energy, is shown for the highest T 2 temperatures of Figure 3.34. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 129 Cross T[H] or T[D] T[0] Section 500 1000 1500 6000 OH OD 0.179 0.258 0.207 0.284 0.226 0.303 4000 OH OD 0.170 0.241 0.201 0.272 0.222 0.294 2000 OH OD 0.159 0.219 0.194 0.258 0.217 0.284 Table 3.8: Energy exchange rate coefficients for O-H and O-D direct elastic cross sections. The temperatures of O and H / D , in degrees Kelvin, are indicated in the table. The tabulated values of the energy exchange rate coefficients are in units of 10~9 cm 3 s _ 1 . corresponding to those used by Clarke and Shizgal, 7 2 Tca = At where A = ri2cr0yj32TrkT2/m and where cr0 — lcm2. If we choose n 2 = 1 c m - 3 , then for a background of hydrogen at 500, 1000, and 1500 K we find values of A corresponding to 2.04,2.88 and 3.53 x 106 s - 1 ; hence, a unit of r is of the order of sub-microseconds of real time for these choices of parameters. Comparing the definition with that given for T in equation (3.3.21), we immediately see that Tcs = '^-^32TrkT2/m KQ For the case of Ti > T 2 , we choose the protons as species 1, and hydrogen as the bath. The results are plotted in Figure 3.36. For the case of T\ > T 2 , we choose hydrogen Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 130 as species 1, and the proton as the bath; the results are plotted in Figure 3.37. It is clear from the figures that the rate of energy transfer depends on the initial temperature (energy) ratio; as is indicated by the energy exchange rate coefficients of Table 3.5, as the temperature disparity becomes more pronounced, the rate at which energy is exchanged between the components increases, and the rate at which cooling or heating occurs is greater than for smaller temperature differences. As in Clarke and Shizgal, 7 2 we define two timescales to better quantify the time evolution process. For a cooling (relaxation) system, we define T\je as the time required for the energy ratio, E/Ethermah to decay to 1/e of its original value, and T I . I as the time it takes for the energy ratio to decay to 1.1. For a heating (excitational) system, we define T\/e as the time for the energy ratio to increase to e times its original value, and T i . i as the time it takes for the energy ratio to reach 1/1.1 = 0.91. These timescales are tabulated for both cases in Tables 3.9 and 3.10. As a check of our time evolution calculations, we compared energy relaxation times for protons with initial energies of 0.646eV (4993 K) and 1.27eV (9816 K) in a hydrogen bath at 298K, using both the cross sections given previously in this chapter and those used by Clarke and Shizgal. 7 2 The time evolution results are illustrated in Figure 3.38, while a comparison of relaxation times is given in Table 3.11. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 131 0 ' 1 I i I i J i I 0.00 0.05 0.10 0.15 0.20 x Figure 3.36: Time evolution of the average energy in a relaxing H + - H system. Solid curves are for the direct-plus-exchange (DPE) cross section, dashed curves for the charge exchange (CE) cross section. Time is in reduced dimensionless units r, correspond-ing to those used by Clarke and Shizgal. 7 2 The labels A , B, and C correspond to pro-ton/hydrogen temperatures of 6000 K/500 K , 4000 K/1000 K , 2000 K/1500 K , respec-tively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 132 o.o 1—1—1—1—I—1—1 1 1 I i i i i I i i i i I i t i i 1 0.0 0.1 0.2 0.3 0.4 0.5 x Figure 3.37: Time evolution of the average energy in a heating H + - H system. Solid curves are for the direct-plus-exchange (DPE) cross section, dashed curves for the charge exchange (CE) cross section. Time is in reduced dimensionless units r , correspond-ing to those used by Clarke and Shizgal. 7 2 The labels A , B , and C correspond to pro-ton/hydrogen temperatures of 6000 K/500 K , 4000 K/1000 K , 2000 K/1500 K , respec-tively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 133 T[H] Cross 500 1000 1500 T[H+] Section T1/e T J . I r 1 / e T J . I r1/e T J . I 6000 D P E C E LTA 0.0158 0.0197 0.0158 0.188 0.226 0.189 4000 D P E C E LTA 0.0680 0.0864 0.0684 0.162 0.203 0.162 2000 D P E C E LTA n/a n/a n/a 0.0683 0.0904 0.0688 Table 3.9: Timescales for the time evolution of the average energy of H for several H + - H cross sections. The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values are timescales, in dimensionless units, for H being heated by H + (i.e. timescales are for the time evolution of the average energy of the hydrogen). A n entry of 'n/a' means that the particular timescale is not applicable for those particular temperatures. The timescales are defined in the text. The dimensionless time units are those of Clarke and Shizgal, 7 2 T = At, where A is defined in the text and is of the order of 2.04 - 3.53 x 106 s _ 1 . Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 134 T[H] Cross 500 1000 1500 T[H+] Section r 1 / e T1a r1/e T1A r 1 / e T1A 6000 D P E C E LTA 0.0381 0.0500 0.0383 0.206 0.287 0.209 4000 D P E C E LTA 0.0921 0.1228 0.0929 0.180 0.243 0.182 2000 D P E C E LTA n/a n/a n/a 0.0733 0.0977 0.0739 Table 3.10: Timescales for the time evolution of the average energy of H + for several H + - H cross sections. The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values are timescales, in dimensionless units, H + being cooled by H (i.e. timescales are for the time evolution of the average energy of the protons). An entry of 'n/a' means that the particular timescale is not applicable for those particular temperatures. The timescales are defined in the text. The dimensionless time units are those of Clarke and Shizgal, 7 2 T = At, where A is defined in the text and is of the order of 2.04 — 3.53 x 106 s _ 1 . Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 135 Figure 3.38: Comparison of the time evolution of the average proton energy. The average energy per proton divided by the average energy per proton at the bath temperature is plotted against dimensionless time units r corresponding to those used by Clarke and Shizgal. 7 2 In all cases, the initial average energy per proton is taken to be 0.646 eV (approx. 4993 K ) , and the (constant) background temperature is 298 K . The solid, short dashed, and long dashed curves represent the time evolution of the energy ratio for the D P E , C E , and direct (D) cross sections, respectively, as calculated using cross sections as defined in Clarke and Shizgal. The dash-dot and dash-dot-dot-dot curves represent the time evolution for D P E and C E cross sections, respectively, as calculated using cross sections as defined previously in this chapter. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 136 E0 = 0.646 eV E0 = 1.27 eV T l / e T l . l T l / e n . i D P E C E LTA 0.030 0.040 0.031 0.202 0.291 0.205 0.023 0.029 0.023 0.216 0.309 0.219 DPE+ CEt LTA2t 0.016 0.022 0.006 0.103 0.137 0.050 0.012 0.017 0.004 0.111 0.148 0.053 DPE* CE* LTA2* 0.016 0.022 0.006 0.080 0.105 0.029 0.012 0.030 0.004 0.077 0.101 0.024 Table 3.11: Comparison of timescales for the time evolution of the average energy per proton for the H + - H system. The (constant) temperature of the background H is taken as 298 K , and the initial proton energy is 0.646 eV (approx. 4993 K ) . The values reported in the table are the timescales T i / e and T\_\, defined in the text, for the D P E and C E cross sections defined previously. The D P E , C E , LTA entries are for calculations using cross sections from this work, D P E , CE , LTA2 entries marked T are calculations using cross sections from Clarke and Shizgal, and D P E , C E , LTA2 entries marked * are the reported results from Clarke and Shizgal. 7 2 A l l times are reported in the dimensionless time units of Clarke and Shizgal, T — At, defined in the text. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 137 3.4 Product Velocity Distribution Functions 3.4.1 Theory The rate of production, Q, of fast moving particles depends on the rate of binary en-counters between the two species involved in a collision and on the velocity distribution functions of the two species. We therefore write that 2 6 Q ( r , C l ) = m(r)n2(r) Jj /i(r, ci)/ 2(r, c'2)go(g, 9,') dQf dc2 (3.4.1) where g, f i ' and o(g,tt') are as defined previously in Chapter 1. The units of Q(r, Ci), above, are c m - 6 s2. The densities of the collision pair are nx(c) and n2(c) and their distribution functions are /i(r, c[) and / 2 (r , c'2), respectively. Altitude dependence enters through the density profiles n\ and n 2 and the temperature profiles, T(r) , contained in the distribution functions. The differential cross section, cr(g, ft'), is required in equation (3.4.1), and plays an important role in the velocity dependence of Q(r, c\). We call Q(r, ci) the product velocity distribution function (PVDF) for species 1. It should also be noted that the distributions in equation (3.4.1) are for the incoming (or pre-collisional) particles, and the integration is over the velocity of the product ion (for a reaction of the type specified by processes (3.1.1) or (3.1.4), for example). With the kinematics of the collision, the velocities of the products (ci, c2) can be related to the velocities of the reactants (c^ , c2) so that Q(r, ci) is a function of the velocity of the translationally energetic product. The hot atom production is thus due to 'inverse' collisions, with the primes denoting the atom velocities before collision. A schematic of the collision process is shown in Figure 3.39. In fact, it is clear from equations (1.4.2) and (3.4.1) that the rate of production, Q(r, Ci) is-identical to the 'gain' (inverse) term of the Boltzmann collision operator. A difficulty in the calculation of Q(r, ci) is that in general the distribution functions Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 138 Figure 3.39: Schematic of an elastic binary collision. Particles 1 and 2 enter into an elastic collision with velocities u' and v', and leave with velocities u and v, respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 139 of nonthermal species are far from equilibrium and are not Maxwellian. However, in order to render the problem tractable it is assumed that the two species are adequately described by Maxwellian distributions at different temperatures, / t ( c : ) = ( 2 ^ ) 2 e X P 2kTi Substituting this form for the distribution functions in equation (3.4.1), we have, for some fixed altitude r, Q(ci) = n m 2 (^^J J/ e-^c'?-x^g J{g,Cl)dg (3.4.7) J 0 where J(g,Cl) = JJe-^Wc^^o*{g^')d^d9 = f e~aC^J^dn (3.4.8) and Ja, = j e - ^ r g ' + 7 g . g ' ] ( T ( 5 ) n y ^ ( 3 4 J ) In equation (3.4.9), it should be noted that the solid angle dQ' is about the direction of g', whereas dQ, is about the direction of g. The geometry is illustrated in Figure 3.40. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 141 Figure 3.40: Coordinate system for the calculation of collisional production. The direc-tion of gis taken as the polar direction. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 142 Using the above geometry and the addition theory for non-coplanar angles, 1 4 9 cos ip' — cos 8 cos ij) + sin 8 sin ij) cos (f> we may write equation (3.4.9) as Jn, = r J, and made the definition /•27T J, = / e - ^ l 5 sin0 sin*/, cos Jo = 2TT I0(8c\g sin 8 sin^>) where I0(z), the modified Bessel function of order zero, 1 5 0 is given by l r W) = -7T Jo e±zcosed8 If we substitute all of the above intermediate results into equation (3.4.8), and use the definition dQ = sin#^ d8^ d(f>psi, we have J(g,Cl) = {2TT)2 T e-™2™eo(g,8)K^d(cos8) (3.4.11) J 0 where Krp = r e-DlCOS^I0(D2sm4>)d(cosi(;) J 0 Di = agci + figcx cos 8 D2 = flgci sin 8 Dl + Dj = g2c\(a2 + 82 + 2a 8 cos 9) The integral may be done analytically, 1 5 1 yielding K = 2smh(xjD2 + D2) \JD\ + D\ Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 143 If we now combine equations (3.4.12), (3.4.11), (3.4.7) and (3.4.3), we have \ 7T J Ci JO J { 9 ^ ] = 2 ( 2 ^ C 1 = r . ^ - ' . M t ^ ^ T O ^ , ) (3.4.12) Jo <7v a; + pz + zap cos 0 Before implementing the result of equation (3.4.12), we wish to non-dimensionalize our various parameters. We define the velocity v0 — 2kT*/m*, where T*,m* are arbitrary temperature and mass values, and then define dimensionless variables x = ci/v0, z = g/v0, A' = Av20, B' = Bv20, a' = av2, B' = /3v20, 7 ' = -yv20, and a' = a/cr0, where A, 5 , a, / 3 , 7 are defined previously in equation (3.4.6), and C T 0 is some factor chosen to non-dimensionalize the differential cross section. With these definitions we may now write equation (3.4.12) as Q(x) = Vgnin2 i^v7^1^2^ 87r2<70 J z3 e~B'z2 J(z,x)dz J{ztx) = r e-A'x2e-<'z2cose(T'(z,0) x Jo s i n h ( W ^ / 2 + ^ 2 + 2a^cos(J) ^ , „ „ =—-a(cos0) (3.4.13) Since J is dimensionless, all the dimensions for Q(x) are contained in the multiplicative factors in front of the integral; it is easy to verify that the units of Q(x) are s 2 c m - 6 . It is possible to perform the integrations in Q(x) and J analytically for a suitable choice of differential cross section. For the hard sphere cross section (see Appendix C ), (7(2,0) = d2/4:, we have ' v W V 7rd2 Q(x) = v y ^ n i l ^ - ^ j ^={N(r,s,u+)-N(r,S,u.)} (3.4.14) where • Ax = a'2 + P'2 Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 144 A2 = 2a'{3' N(r,s,t) = -1 erf(s) VrTi2 t exp st t = u± u± = { j (a> + P)y/£ for + The details of the derivation are given in Appendix E. 3.4.2 Total Collision Rates The preceding section described the rate of production of fast moving particles of speed c\ by binary encounters between two species; that is, Q(c\) represents the rate at which collisions which create particles of species 1 with speed c\ occur. In order to calculate the total number of collisions between species 1 and 2 we must integrate over all possible speeds c\. It should be noted that, in the absence of any particle sources or sinks, integration of the Boltzmann collision operator over dc\ must be equal to zero because of particle conservation;29 that is, particles of species 1 and 2 may be moved between different elements of velocity space by elastic collisions, but the total number of particles remains constant. We have, from equation (1.4.2), Either term above defines the total collision rate between species 1 and 2. Using the definition of the production rate, equation (3.4.1), and the definition of the collision 0 and so Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 145 frequency, Z{ci) = JJ f2a(g,tt)gdQ,dc2 fi o-tot(g) g dc2 (3.4.15) we may thus calculate the total collision rate using either of Qtot = { J j <5(ci)rfci (3.4.16) ( j f1Z(c1)dc1 These total collision rates may also be compared with the work of Fitzpatrick and Shizgal 1 4 5 for the hard sphere case in order to verify the calculation of the product velocity distribution functions/collision frequencies. 3.4.3 Rate of Production of Escaping Atoms Once we have calculated the product velocity distribution function, Q(ci), we are in-terested in making some statement about the number of hot atoms which are able to escape. A detailed discussion of the calculation of a rigorous, altitude dependent escape flux which includes the moderating effect of the ambient background atmosphere is dis-cussed in Chapter 4. For the moment, we may use the product velocity distribution functions to compare the relative rates of production of escaping atoms. This incor-porates the nonthermal production of hot atoms and the use of realistic collision cross sections. It is assumed that any hot atom produced by nonthermal excitation (via any of the elastic collision mechanisms we have previously detailed in the calculation of the energy exchange rate coefficients: DIR, C E , or the LTA) with velocity greater than the escape velocity escapes from the planetary atmosphere. Thus, in order to calculate the rate of production of escaping atoms, we simply integrate the product velocity distribution func-tion (or the gain term of the Boltzmann collision operator) normalized by the densities Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 146 of the colliding species, over all velocities C\ greater than the escape velocity, Air f°° a e s c = / Qic^cldcr (3.4.17) n\Tl2 Jccsc where cesc = \j2GM/r. It may be easily verified using equation (3.4.1) that aesc has units of a rate, cm 3 s _ 1 . Equation (3.4.17) ignores reverse collisions, which reduce the overall rate at which escaping atoms are produced. We incorporate reverse collisions by adding the loss term of the Boltzmann collision operator to Q{c\) before integrating; that is, aesc = / [Q{ci)- hZ{cx)]cldc1 (3.4.18) nin2 JCesc where Z(ci) is the collision frequency between species 1 and 2, defined previously. We may also utilize the PVDFs in order to estimate the fraction of hot atoms created by a nonthermal process which can escape from the atmosphere. This is given by Air r°° airs:c = a e s c / / Q{c1)c\dc1 (3.4.19) n\n2 Jo 3.4.4 Results and Discussion The calculation of the PVDF ' s in the present work follows closely the presentation and calculations done by Shizgal. 7 1 It is a direct method in comparison to the Monte Carlo integration used by Hodges and Breig 7 3 to calculate the so-called 'differential scatter-ing coefficients'. The calculation of the product velocity distributions (PVDFs), from equation (3.4.13), requires a double integration, over the scattering angle 6 and reduced energy z. This integration requires the differential elastic cross section(s), o(g,9). The cross section(s) are calculated at 91 energies, from 1 x 1G 1 - 4 to 10 eV, and at 361 angles (half degree increments), and are illustrated in Figures 3.23 and 3.24. In the calculation of the PVDFs , we use a Simpson's Rule with one-half degree spacing of grid points for the Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 147 6 integration. An interpolation scheme using z points located at the 91 energies used in calculating the differential cross section is used to perform the z integration (the scheme is detailed in Appendix E). The linear trajectory approximation adopted in the tables and graphs that follow is that given by equation (3.2.6), unless otherwise indicated. Some examples of the product velocity distribution function normalized by the species' densities, Q(ci)/n-iri2, for H + - H are illustrated in Figures 3.41-3.43. We have chosen the same pairs of temperatures as in the calculation of the energy exchange rate coefficients because they approximate ion and neutral temperatures approximating exospheric con-ditions on the terrestrial planets. For comparison with the work of Hodges and Breig, 7 3 we have also calculated and plotted the 'differential scattering rate coefficient' defined as Aire P ( C I ) = —L}\f\ (3.4.20) dec where P(c\) is equivalent to — calculated by Hodges and Breig. In equation (3.4.20), J[/] is the Boltzmann collision operator, as defined by equation (1.4.2). The units of P(c\) and Q ( c i ) / n i n 2 in the figures are 10 1 6 cm 2 and 10 1 4 s2, respectively. In the calculation of the collision operator, we use equation (1.4.2), with the gain term given by the product velocity distribution function, Q(ci), as defined in equation (3.4.1), and the loss term given by the species 1 distribution function times the collision frequency, Z(ci) , which is given by equation (3.4.15). For the collision frequency integration we use 400 Simpson's Rule points to ensure convergence for all values of c\. It should be noted that the vertical scale varies between figures in order to accentuate the details in the form of the PVDFs . Figures 3.41-3.43 have similar features. The plots of P(ci) may be interpreted as indicative of the transfer of energy from the hot protons to the cold hydrogen. The negative portions of the curve are due to the removal of hydrogen for those velocities (the 'loss' term of the Boltzmann collision operator), while the positive portions are due to the addition of hot hydrogen for those velocities (the 'gain' term of the Boltzmann Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 148 > | " " | M I I | . II I | I 1 I I | I I U | I I I I | I 0 5 10 15 2 0 2 5 3 0 3 5 4 0 c. (km/s) (a) (b) 10 15 20 25 C| (km/s) 75 ;l 1 1 1 y~ r<<1111111 1 1 | 1 M 1 | M 1 1. 50 j- I ; "j 25 j- "j 0 i- j , / A - 2 5 J ' ' 1 - i - 5 0 J ' ^ '/ -= -75 100 ' r -z 125 i -j 150 175 200 225 - i 1 . . . .= 15 2 0 25 c, (km/s) (c) (d) Figure 3.41: Product velocity distribution functions for the direct-plus-exchange (DPE)cross section for H + - H . In (a), the hydrogen temperature is fixed at 500 K , and the proton temperature is varied, with values of 2000 K (long dashed), 4000 K (short dashed), and 6000 K (solid). In (b), the proton temperature is fixed at 6000 K , and the hydrogen temperature is varied, with values of 1500 K (long dashed), 1000 K (short dashed), and 500 K (solid). In (c) and (d) we have plotted -P(ci), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 149 (a) (b) 10 15 20 25 30 35 40 Cj (km/s) 50 .' 1 1 1 | 1 ' ' M i l 1 1 j 1 1 1 1 J I 1 1 1 [ 1 1 1 1 j 1 1 I 1 [ 1 1 I T | 25 : / ' V : 0 / '/ — : Y\ \ i •25 Vi 1 I\ i1 -o -50 '\\\ '' -Cu -75 -100 I ^ -125 -j -150 "I IV i i . . . .: 10 1 5 20 25 30 35 40 c, (km/s) (c) (d) Figure 3.42: Product velocity distribution functions for the charge exchange (CE) cross section for H + - H . In (a), the hydrogen temperature is fixed at 500 K , and the proton temperature is varied, with values of 2000 K (long dashed), 4000 K (short dashed), and 6000 K (solid). In (b), the proton temperature is fixed at 6000 K , and the hydrogen temperature is varied, with values of 1500 K (long dashed), 1000 K (short dashed), and 500 K (solid). In (c) and (d) we have plotted -P(ci), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 150 c, (km/s) c, (km/s) (c) (d) Figure 3.43: Product velocity distribution functions for the linear trajectory approxima-tion (LTA) cross section for H + - H . In (a), the hydrogen temperature is fixed at 500 K , and the proton temperature is varied, with values of 2000 K (long dashed), 4000 K (short dashed), and 6000 K (solid). In (b), the proton temperature is fixed at 6000 K , and the hydrogen temperature is varied, with values of 1500 K (long dashed), 1000 K (short dashed), and 500 K (solid). In (c) and (d) we have plotted P(c\), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 151 collision operator). It is clear from subfigure (c) in Figures 3.41-3.43 that for a fixed hydrogen temperature, increasing the proton temperature increases the average energy of the cold hydrogen. Graphically, this is illustrated by the increase in the position of the peak and width of the positive portion of the P(ci) curve, respectively. From subfigure (d) in Figures 3.41-3.43, we see that for a fixed proton temperature, increasing the hydrogen temperature results in a smaller number of hydrogen atoms being energized, but with a greater average increase in energy. Graphically, this is illustrated by the decrease in the depth of the negative portion of the P(c\) curves and the change in position of the peak and width of the positive portion of the P(c\) curve, respectively. The pattern of change with hydrogen and proton temperatures is the same for each of the D P E , C E and LTA cross sections; however, the relative magnitude of the changes varies between them. A comparison of the PVDFs for the three cross sections for fixed hydrogen and proton temperatures is shown in Figure 3.44. For subfigures (a) and (c), the proton temperature is 6000 K and the hydrogen temperature is 500 K , representing the largest difference in colliding partner temperatures. For subfigures (b) and (d), the proton temperature is 2000 K and the hydrogen temperature is 1500 K , representing the smallest difference in colliding partner temperatures. The difference between the large energy and small energy cases is evident by the difference in vertical scale between subfigures (a) and (b) or (c) and (d) and the change in shape of the PVDF ' s and differential scattering coefficients. The D P E includes both direct elastic momentum transfer and charge exchange, and so the result that this cross section is most efficient at producing hot hydrogen is not unexpected. For only charge exchange, the effects are reduced by approximately 30-40% when considering the plots of -P(ci); we may thus deduce that charge exchange is relatively efficient, as compared to direct elastic momentum transfer collisions, since the C E values are 60-70% of the magnitude of the D P E values. This roughly holds even Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 152 c, (km/s) c, (km/s) (c) (d) Figure 3.44: Product velocity distribution functions for the D P E , C E , and LTA cross sections for H + - H . In (a), the hydrogen temperature is 500 K and the proton temperature is 6000 K . In (b), the hydrogen temperature is 1500 K and the proton temperature is 2000 K . For both cases, we show the resulting PVDFs for the LTA (long dashed), C E (short dashed), and D P E (solid) cross sections. In (c) and (d) we have plotted P(c\), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 153 when considering the plots of the PVDFs , Q(ci), when considering the higher energy portion (greater than approximately 10 km/s). The LTA cross section behaves similarly to the C E one, except that it appears to overestimate both the rate at which energy is transferred from the protons to the hydrogen and the change in the average energy of the cold hydrogen. These overestimations are visible as a subtle shift toward higher velocities in the Q(c\) plots, and a change in the width and depth of the positive and negative portions of the P{c\) plots. We also plot P V D F s and differential scattering coefficients for the deuterium-proton system. These are illustrated in Figures 3.45-3.46 for the DIR and C E interactions for a range of ion and neutral temperatures. The features of these curves are very similar to those for the hydrogen-proton curves, and most of the discussion for that system also applies to this one. The disparate masses of the two colliding partners makes energy transfer slightly less efficient for the deuterium- proton or deuterium ion-hydrogen systems than for the hydrogen-proton system. This is reflected in the vertical scale of the P V D F ' s compared to those for the hydrogen-proton system. The PVDFs and differential scattering coefficients for the oxygen-hydrogen and oxygen-deuterium systems are shown in Figures 3.47-3.48. The figures represent results for the DIR interaction for a range of collision pair temperatures. The general shape and features of these curves are very similar to those of the previous systems. The magnitude of the P V D F s decrease as we move from D-H, O-D, and 0-H, and we also note that the peak of the PVDFs shift toward higher velocities. For the case where the colliding partners are close in temperature the differential scattering coefficients, -P(ci), show some structure for velocities below approximately 5 km/s. There are several checks we may perform in order to verify that we are correctly calculating the PVDFs . One good initial check is that the condition given by equation (3.4.16) is satisfied; we have verified that this is correct within 2-7% depending on the Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 154 C| (km/s) c, (km/s) (c) (d) Figure 3.45: Product velocity distribution functions for the direct elastic (DIR)cross section for H + - D . In (a), the deuterium temperature is fixed at 500 K , and the proton temperature is varied, with values of 1000 K (dot-dashed), 2000 K (short dashed), 4000 K (long dashed), and 6000 K (solid). In (b), the deuterium temperature is fixed at 1500 K , and the proton temperature is varied, with values of 1500 K (dot-dashed), 2000 K (short dashed), 4000 K (short dashed), and 6000 K (solid). In (c) and (d) we have plotted P(ci) , as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 155 c, (km/s) c, (km/s) (c) (d) Figure 3.46: Product velocity distribution functions for the charge-exchange (CE)cross section for D + - H . In (a), the hydrogen temperature is fixed at 500 K , and the deuterium ion temperature is varied, with values of 2000 K (long dashed), 4000 K (short dashed), and 6000 K (solid). In (b), the hydrogen temperature is fixed at 1500 K , and the deuterium ion temperature is varied, with values as in (a). In (c) and (d) we have plotted -P(ci), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 156 c, (km/s) C l (km/s) (c) (d) Figure 3.47: Product velocity distribution functions for the direct elastic (DIR) cross section for 0 -H. In (a), the hydrogen temperature is fixed at 500 K , and the oxygen temperature is varied, with values of 2000 K (long dashed), 4000 K (short dashed), and 6000 K (solid). In (b), the hydrogen temperature is fixed at 1500 K , and the oxygen temperature is varied, with values as in (a). In (c) and (d) we have plotted P(c\), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 157 (c) (d) Figure 3.48: Product velocity distribution functions for the direct elastic (DIR) cross section for O-D. In (a), the deuterium temperature is fixed at 500 K , and the oxygen temperature is varied, with values of 2000 K (long dashed), 4000 K (short dashed), and 6000 K (solid). In (b), the deuterium temperature is fixed at 1500 K , and the oxygen temperature is varied, with values as in (a). In (c) and (d) we have plotted P(ci), as defined in the text, for the same temperatures and with the same labeling as for (a) and (b), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 158 temperature ratio chosen. Graphically, from Figures 3.41-3.43, this means that the area under the negative and positive portions of the P(ci) curves must be equal, since particle conservation requires that the integral of the collision operator over all C i must be zero. As a second check of the accuracy of our PVDFs , we use them to calculate the energy exchange rate coefficients calculated in a previous section of this chapter. That is, we calculate We use equations 3.3.17 and (3.3.19) to obtain the energy exchange rate coefficient, kE, from RE as given above. The results of this calculation are given in Tables 3.12, 3.13 and 3.14. By comparison with the earlier results (using equation (3.3.16)), given in Tables 3.5, 3.7 and 3.8, we see that the two very different methods of calculation agree to better than 10% over most of the temperature ranges, and never differ by more than 20%. To further quantify the relative accuracy of the two methods, we compare the results for the energy exchange rate coefficient for the hard sphere cross section (the hard sphere cross section is used because for this choice it is possible to calculate the analytic solution of equation (3.3.16), given by equation (3.3.17)). The results are summarized in Table 3.15 in the form of the ratio of kE from equation (3.3.16) divided by the value from equa-tion (3.4.21). For all pairs of temperatures, equation (3.3.16) reproduced the expected analytic solution given by equation (3.3.17) to better than 1%, and so the deviation from unity in Table 3.15 may be attributed entirely to the accuracy with which we are to calculate the triple integral of equation (3.4.21). We perform another check on the PVDFs by calculating the charge exchange reaction rate. That is, we calculate the rate of both the forward and backward reaction for deuterium and hydrogen, (3.4.21) D+ + H D + H+ Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 159 Cross T[H] T[H+] Section 500 1000 1500 6000 DIR 2.47 2.58 2.72 C E 8.94 9.26 9.57 LTA 11.66 12.02 12.37 4000 DIR 2.42 2.61 2.91 C E 7.78 8.25 8.73 LTA 10.08 10.59 11.10 2000 DIR 2.18 2.50 2.62 C E 5.90 6.48 7.32 LTA 7.90 8.56 9.48 Table 3.12: Energy exchange rate coefficients for several H+-H cross sections, calculated by integration of the product velocity distribution functions. The temperatures of both H and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of the energy exchange rate coefficients are in units of 10~9 cm 3 s _ 1 . Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 160 T [neutral] Cross 500 OT 1500 T[ion] Section I II I II I II 6000 DIR 1.897 1.874 1.987 1.994 2.101 2.145 C E 8.388 6.301 8.577 6.748 8.772 7.170 LTA 10.124 8.015 10.298 8.487 10.476 8.915 4000 DIR 1.953 1.770 2.109 1.950 2.345 2.238 C E 7.026 5.328 7.275 5.851 7.548 6.372 LTA 8.807 6.928 9.059 7.510 9.331 8.066 2000 DIR 1.797 1.670 2.094 1.928 3.065 2.908 C E 5.139 4.233 5.485 4.756 6.052 5.376 LTA 6.735 5.654 7.119 6.236 7.716 6.884 Table 3.13: Energy exchange rate coefficients for several H + - D cross sections, calculated by integration of the product velocity distribution functions. The temperatures of both deuterium and the protons, in degrees Kelvin, are indicated in the table. Values under columns marked 'I ' are for reactions where the test particle is the deuterium atom, with a bath of hot protons; values under columns marked 'II' are for reactions where the test particle is hydrogen, with a bath of hot deuterium ions. The tabulated values of the energy exchange rate coefficients are in units of 10~9 cm 3 s - 1 . Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 161 Cross T[H] or T[D] T[0] Section 500 1000 1500 6000 OH OD 0.184 0.261 0.211 0.287 0.230 0.306 4000 OH OD 0.179 0.247 0.209 0.278 0.229 0.300 2000 OH OD 0.185 0.236 0.223 0.279 0.256 0.317 Table 3.14: Energy exchange rate coefficients for O-H and O-D direct elastic cross sections. The temperatures of O and H / D , in degrees Kelvin, are indicated in the table. The tabulated values of the energy exchange rate coefficients are in units of 1 0 - 9 cm 3 s"1. T[H] T[H+] 500 1000 1500 6000 0.95 0.95 0.94 4000 0.96 0.95 0.92 2000 0.98 0.96 0.88 Table 3.15: Comparison of energy exchange rate coefficients for the hard sphere cross section using equations (3.3.16) and (3.4.21). The temperatures of both H and the protons, in degrees Kelvin, are indicated in the table. The values reported in the table are the ratio of the result from equation (3.3.16) divided by that of equation (3.4.21). Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 162 The calculation is done both by direct integration of the the loss term and direct integra-tion of the P V D F , as given by equation (3.4.16), and is illustrated in Figure 3.49. The results of Hodges 7 4 are also shown in the figure. It is noted that he calculates different rates for the forward/reverse charge exchange reaction, since he does not neglect the cou-pling coefficient in the calculation of the cross sections for D H + , although it is difficult to account for the magnitude of the difference between them given that the coupling coefficients are of the order of 1 0 - 4 the size of the coefficients of the primary channels. As mentioned previously, we have neglected the cross-coupling interaction, and so our rates for the forward and reverse reactions are identical. It is clear, however, that there is good agreement with Hodges' calculated rate coefficients. It should be noted that, as in most of the work done by Hodges, Monte Carlo integrations are used to perform all of the velocity integrals. Our calculation, in the form of a triple integral, does not have to resort to the Monte Carlo methodology in order to extract the rate coefficients. The results of the calculation of the rate of production of escaping hot hydrogen atoms, aesc, are given in Tables 3.16 -3.18. The same set of proton and hydrogen temperatures are used as in the calculations of the energy exchange rate coefficients. Values in columns marked A are calculations using equation (3.4.18) and include reverse collisions; values in columns marked B use equation (3.4.17), ignoring reverse collisions. We assume a hydrogen exobase altitude of 200 km for Venus, 3 5 5 00 km for Earth, 4 4 and 250 km for Mars. 6 The escape velocities for these altitudes are 10.2 km/s, 10.78 km/s, and 4.85 km/s, respectively. We see several interesting trends in the values of the rate of production of escaping hot hydrogen for the terrestrial planets, as given in Tables 3.16 -3.18. The LTA cross section results overestimate the rate of production of hot atoms by 20-50% for Venus and Earth, and 5-30% for Mars, in good agreement with the discrepancy of 30-50% reported by Shizgal. 1 3 1 This is observed for production of hot atoms including and excluding reverse Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 163 Figure 3.49: Dependence of charge exchange rate coefficients on the effective temperature for D + - H . The value of the charge exchange rate coefficients are plotted as a function of the effective temperature, T e / / . The solid curve is from integration of the P V D F , the dashed curve from integration of the Boltzmann collision operator loss term. The circle and plus symbols are the results of fits to Monte Carlo calculations by Hodges7 4 for the rates of D+-H and D - H + , respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 164 T[Hj Cross 500 1000 1500 T[H+] Section A B A B A B DIR 6000 C E LTA DIR 4000 C E LTA DIR 2000 C E LTA 0.850 0.851 3.282 3.282 4.130 4.130 0.475 0.475 1.894 1.894 2.542 2.542 0.092 0.092 0.404 0.404 0.603 0.603 0.934 1.040 3.385 3.426 4.166 4.207 0.544 0.644 1.986 2.024 2.575 2.613 0.120 0.214 0.428 0.463 0.591 0.626 0.966 1.706 3.289 3.572 3.988 4.271 0.574 1.273 1.901 2.166 2.407 2.672 0.116 0.773 0.307 0.552 0.402 0.647 Table 3.16: Rate of production of escaping H on Venus for several H + - H cross sections. The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s - 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 165 T[H] Cross 500 1000 1500 T[H+] Section A B A B A B DIR 0.741 0.741 0.116 0.773 0.854 1.336 6000 C E 3.033 3.033 3.140 3.161 3.104 3.289 LTA 3.859 3.859 3.906 3.927 3.796 3.981 DIR 0.392 0.392 0.449 0.500 0.486 0.942 4000 C E 1.656 1.660 1.753 1.772 1.717 1.891 LTA 2.256 2.256 2.299 2.318 2.193 2.366 DIR 0.064 0.064 0.085 0.133 0.089 0.520 2000 C E 0.297 0.297 0.323 0.341 0.246 0.407 LTA 0.450 0.450 0.450 0.468 0.322 0.484 Table 3.17: Rate of production of escaping H on Earth for several H + - H cross sections. The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s _ 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 166 TJfl] Cross 500 1000 1500 T[H+] Section A B A B A B DIR 2.196 4.311 1.582 8.630 1.211 11.466 6000 C E 4.597 5.385 3.044 5.675 2.072 5.914 LTA 5.107 5.894 3.411 6.042 2.345 6.187 DIR 1.703 3.626 1.179 7.624 0.881 10.310 4000 C E 3.622 4.322 2.308 4.666 1.486 4.953 LTA 4.183 4.883 2.696 5.053 1.756 5.223 DIR 0.971 2.653 0.541 6.248 0.308 8.740 2000 C E 2.171 2.764 1.127 3.154 0.465 3.485 LTA 2.711 3.304 1.427 3.455 0.592 3.612 Table 3.18: Rate of production of escaping H on Mars for several H+-H cross sections. The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s - 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 167 collisions. As may be expected given their similar size and escape velocities, the results for a e s c for Venus and Earth are very similar. From Tables 3.16 and 3.17 the data for both planets indicate that for cases where a large separation occurs between the temperatures of the protons and hydrogen that the LTA cross section yields larger aesc values than either the C E or the DIR cross sections; as the temperature of the hydrogen increases, the temperature discrepancy diminishes, and the LTA predicted aesc values fall below those using the DIR cross section, although they continue to exceed aesc values using the C E cross section. We observe a monotonic increase of the values of aesc with an increase in hydrogen temperatures, for fixed proton temperature, for calculations of aesc excluding reverse collisions. Physically, this is reasonable; imparting some fixed amount of energy to a distribution of particles of higher average energy results in more particles reaching the escape energy. The situation is roughly reversed when we include reverse collisions, where aesc decreases as we increase the hydrogen temperature, for a fixed proton temperature. For low hydrogen temperatures the rate of reverse collisions is insufficient to offset the rate of transfer of nonthermal energy, and so the values for aesc including and excluding reverse collisions are roughly identical. However, as the hydrogen temperature increases, the rate of reverse collisions increases, and so the discrepancy between the values of aesc including and excluding reverse collisions increases. It is important to note that even for the highest value of the hydrogen temperature considered here, 1500 K , the most probable hydrogen speed is still only approximately 5 km/s, less than the escape speed. We would expect that even with the aid of nonthermal energy exchange, the number of hot hydrogen atoms produced which avoid thermalization to escape is relatively modest, given the high escape velocity. The situation is much different for Mars, primarily because of its much lower escape speed. From Table 3.18 the predicted values of aesc are 1.5-3 times that for Venus or Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 168 Earth. When reverse collisions are neglected, we observe a monotonic increase of the values of aesc with an increase in hydrogen temperatures, for fixed proton temperature. However, when we include reverse collisions, the observed values of aesc increase and then decrease with increasing hydrogen temperature for fixed proton temperature. This behaviour is most likely related to a change in the balance between the rate of production, the rate of thermalization, and the value of the mean thermal hydrogen speed relative to the escape speed. As for Venus and Earth, for low hydrogen temperatures, the rate of reverse collisions is insufficient to offset the rate of production, and the values of aesc including and excluding reverse collisions are approximately the same. As the hydrogen temperature increases, the mean thermal hydrogen speed begins to exceed the escape speed, and so we sample more of the hydrogen distribution when we calculate the rate of production of escaping hydrogen. That is, the peak of the hydrogen distribution is very close to the escape speed, and the bulk of the hydrogen atoms are involved in the integration of the reverse collision term. As the hydrogen speed increases further, even more of the hydrogen distribution is sampled, and the mean thermal hydrogen speed is in excess of the escape speed. The rate of reverse collisions is now calculated integrating over the bulk of the hydrogen distribution (as opposed to the tail of the distribution when the mean thermal speed of the hydrogen is less than the escape speed, as for Venus and Earth), and so the rate of reverse collisions increases faster than the rate of production of hot atoms, slowing the overall rate of production of hot hydrogen atoms. For a species with a much heavier mass, such as oxygen, we would expect that the mean thermal speed would be much lower for a given temperature, and that the values of aesc would be similar to that of Venus and Earth for hydrogen. We may again compare our results with those of Hodges and Breig, 7 3 whose results for the production of hot hydrogen for Earth are summarized in Table 3.19. Once again we see a (variable) discrepancy between our results and their work, similar but not identical Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 169 Cross T[H] T[H+] Section 500 1000 1500 6000 DIR LTA 6.79 7.34 6.60 6.88 6.26 6.43 4000 DIR LTA 3.56 4.18 3.54 3.86 3.25 3.49 2000 DIR LTA 0.58 0.78 0.55 0.67 0.43 0.50 Table 3.19: Rate of production of escaping H on Earth for several H + - H cross sections taken from the results of Hodges and Breig. 7 3 The temperatures of both hydrogen and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s - 1 . to the discrepancy in the values of the energy exchange rate coefficients. As before, it is difficult to comment on the source of this discrepancy. The same calculation of rates of production of escaping particles is done for deuterium-proton systems. These results are reported in Tables 3.21-3.22. The interpretation is similar to those done for hydrogen-proton collisions. We can see from the tables that it is more difficult to produce escaping deuterium than escaping hydrogen, regardless of the type of interaction considered. The values in the tables of rate of escape production of hot deuterium by collisions with protons are roughly a factor of 2-10 times smaller than those for hydrogen- proton collisions. In addition, direct elastic collisions are much less efficient at production of hot neutrals than charge exchange reactions. This is because although some momentum is transferred in each charge exchange collision, the pre-collisional ion keeps most of its energy when it becomes a post-collisional neutral. The same calculation of rates of production of escaping particles is done for oxygen-hydrogen systems, and oxygen-deuterium systems. These results are reported in Tables Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 170 T [(neutral) Cross 5 0 0 - iOOfj 1500 T(ion) Section A B A B A B 6000 DIR C E 0.117 1.340 0.117 1.340 0.145 1.370 0.145 1.370 0.176 1.394 0.191 1.400 4000 DIR C E 0.042 0.450 0.042 0.450 0.055 0.465 0.055 0.465 0.071 0.484 0.085 0.489 2000 DIR C E 0.184 0.018 0.184 0.018 0.310 0.020 0.333 0.021 0.520 0.021 1.800 0.026 Table 3.20: Rate of production of escaping D on Venus for several H + - D cross sections. The temperatures of both deuterium and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of a e s c are in units of 1 0 - 9 cm 3 s - 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. T(neutral) T(ion) Cross Section A 500 B A 1000 B A 1500 B 6000 DIR 0.072 0.072 0.096 0.096 0.122 0.127 C E 1.069 1.069 1.102 1.103 1.132 1.135 4000 DIR 0.023 0.023 0.033 0.033 0.044 0.050 C E 0.314 0.314 0.330 0.330 0.350 0.351 2000 DIR 7.3x10- 4 7.3xl0" 4 1.3x10" 3 1.4xl0" 3 2.5x10" 3 7.7xl0- 3 C E 0.009 0.009 0.010 0.010 0.011 0.013 Table 3.21: Rate of production of escaping D on Earth for several H + - D cross sections. The temperatures of both deuterium and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s - 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 171 T(neutral) Cross 500 1000 1500 T(ion) Section A B A B A B 6000 DIR C E 1.45 4.51 1.61 4.57 1.39 3.95 3.49 4.73 1.17 3.11 5.90 4.89 4000 DIR C E 1.04 3.23 1.18 3.29 0.99 2.77 2.92 3.47 0.82 2.07 5.17 3.65 2000 DIR C E 0.43 1.38 0.56 1.43 0.38 1.04 2.09 1.63 0.25 0.50 4.11 1.85 Table 3.22: Rate of production of escaping D on Mars for several H + - D cross sections. The temperatures of both deuterium and the protons, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s - 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. 3.24-3.25. The values here are smaller by several orders of magnitude than rate of pro-duction for hydrogen-proton or deuterium-proton systems, mainly due to the disparate masses of the collisional pair. It should be noted that for all of the terrestrial planets, the production of deuterium is only approximately l/10th that of hydrogen due to direct elastic collisions with hot oxygen. While energy transfer overall is more efficient for col-lision partners with similar masses, the lighter hydrogen has more particles created with velocities in the high energy 'tail ' of the P V D F , which is the portion that contributes to the rate of production of escaping hot product atoms. We may utilize equation (3.4.19) in order to approximate the fraction of escaping hot atoms produced by a given nonthermal process. This estimate does not account for thermalization effects or reverse collisions. We assume that the hot oxygen atoms have an initial velocity of 5.6 km/s, corresponding to approximately 2.5 eV, and that the hydrogen Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 172 T[H] fc T[D] Cross 500 1000 1500 T[0] Section A B A B A B 6000 OH OD 0.0203 0.0047 0.0204 0.0047 0.0598 0.0101 0.0916 0.0102 0.0906 0.0196 0.3159 0.0240 4000 OH OD 0.0060 0.0006 0.0061 0.0006 0.0290 0.0024 0.0606 0.0025 0.0470 0.0067 0.2704 0.0111 2000 OH OD 0.0006 6.0xl0- 5 0.0007 6.0xl0- 5 0.0069 0.0001 0.0382 0.0002 0.0082 0.0007 0.2257 0.0050 Table 3.23: Rate of production of escaping H on Venus for direct elastic collisions between O-H and O-D. The temperatures of the hydrogen, deuterium, and oxygen atoms, in degrees Kelvin, are indicated in the table. The tabulated values of a e s c are in units of 1 0 - 9 cm 3 s _ 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. T[H] fc T[D] Cross 500 1000 1500 T[0] Section A B A B A B 6000 OH OD 0.0128 0.0027 0.0129 0.0027 0.0425 0.0060 0.0592 0.0060 0.0728 0.0122 0.2228 0.0140 4000 OH OD 0.0033 0.0003 0.0033 0.0003 0.0193 0.0012 0.0359 0.0012 0.0367 0.0037 0.1856 0.0055 2000 OH 0.0002 0.0003 0.0043 0.0207 0.0063 0.1518 OD 1.4xl0- 5 1.5xl0~5 1.5xl0- 5 1.5xl0~5 0.0003 0.0021 Table 3.24: Rate of production of escaping H on Earth for direct elastic collisions between 0 - H and O-D. The temperatures of the hydrogen, deuterium, and oxygen atoms, in degrees Kelvin, are indicated in the table. The tabulated values of a e s c are in units of 1 0 - 9 cm 3 s _ 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 173 T[H] fe T[D] Cross 500 1000 1500 T[0] Section A B A B A B 6000 OH OD 0.330 0.250 0.781 0.283 0.189 0.256 1.789 0.714 0.081 0.194 2.497 1.262 4000 OH OD 0.226 0.152 0.668 0.185 0.117 0.164 1.677 0.612 0.036 0.117 2.327 1.162 2000 OH OD 0.102 0.052 0.533 0.084 0.034 0.057 1.489 0.493 0.001 0.024 1.956 1.031 Table 3.25: Rate of production of escaping H on Mars for direct elastic collisions between O-H and O-D. The temperatures of the hydrogen, deuterium, and oxygen atoms, in degrees Kelvin, are indicated in the table. The tabulated values of aesc are in units of 1 0 - 9 cm 3 s - 1 . Columns marked ' A ' are for production of hot atoms including the effect of reverse collisions; columns marked ' B ' neglect reverse collisions. is at a temperature of 0, 100, 200 and 300 K . We summarize the results of the production of hot hydrogen from elastic O-H collisions on Venus in Table 3.26, along with results calculated by Cooper et a l . , 7 5 McElroy et a l . 3 5 and Gurwell and Yung. 7 6 It should be noted that the results of McElroy et al. assume isotropic scattering in their calculations. Cooper et al. perform calculations for both 'real' and isotropic scattering systems, and assume a Maxwellian distribution for the hot hydrogen product. Gurwell and Yung use a Henyey-Greenstein (HG) functional fit to the integrated angular distribution functions derived from the cross section calculations of Cooper et al. for a relative impact energy of 0.15 eV. They assume that the same fit is applicable to the OD system. They modify the anisotropy parameter in the H G function to approximate forward and isotropic scattering. They use several different approximations to the branching ratios for the 0^ dissociative recombination process, noting that while there is fairly large variation in the values of Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 174 the branching ratios the average energy of the final hot oxygen distribution is generally fairly uniform. Their final calculation of product hot atom velocity distributions is done by Monte Carlo. Our calculation assumes that the hot oxygen distribution function is Maxwellian, with a peak about the energy of the main dissociation branch (+5 eV). We also assume that the product hot atoms are Maxwellian. We utilize the (anisotropic) quantum mechanical differential cross sections directly for our calculation. As can be seen from the table, our calculated escape fractions lie between those of Cooper et al. and Gurwell and Yung, using anisotropic scattering, and those of McElroy et al. and Gurwell and Yung using isotropic scattering. Gurwell and Yung note the discrepancy between their results and those of Cooper an co-workers, especially for the case of anisotropic scattering, especially since they used the cross section data of Cooper et al. in the H G fit. However, as has been mentioned previously, published values for O-H cross sections in Cooper et al. are most likely in error, and it would therefore seem likely that calculations based on them would also be in error. The product velocity distribution functions Q(c\) for both O-H and O-D for various initial hydrogen and deuterium temperatures are shown in Figure 3.50). The results for both systems show a distribution of production of hot product speeds dominated by a peak at lower speeds and with a long tail which slowly decays with increasing speed. Their form is remarkably similar to the results of Gurwell and Yung, as given in their Figure 5. The same calculations may be performed for O-D. Escaping production fractions for both O-H and O-D direct elastic collisions for all three terrestrial planets are given in Table 3.27. McElroy et a l . 3 5 predicted that deuterium escape resulting from collisions with hot oxygen was negligible compared to hydrogen escape resulting from the same process. This was because their scattering model resulted in a maximum scattered prod-uct velocity which was less than the escape velocity on Venus for all branches of the O2 Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 175 1.6 0.4 0.0 I I I — ^ \ 1I "'/(• -1 -1 10 15 Cj (km/s) (a) 10 15 C, (km/s) (b) Figure 3.50: Product velocity distribution functions for 0 -H and O-D due to direct elastic scattering. The 0 -H results are illustrated in (a), O-D results in (b). Initial oxygen most probable velocity is 5.6 km/s. In both figures, the hydrogen/deuterium temperatures are 100 K (long-dashed), 200 K (short-dashed), and 300 K (solid), respectively. Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 176 T[H] (K) Present Work C l C2 Mc GY1 GY2 0 9.73 6.0 100 9.91 5.1 11.6 200 10.16 6.9 15.5 300 10.41 8.5 18.8 15.0 4.6 15.8 Table 3.26: Escape fraction for hot H from O-H elastic collisions on Venus. The results of Cooper et a l . 7 5 for forward and isotropic scattering are given in C l and C2; results of McElroy et a l . 3 5 in Mc; and results of Gurwell and Young 7 6 for forward and isotropic scattering, in GY1 and GY2, respectively. The pre-collisional hot oxygen atom is assumed to have a velocity of 5.6 km/s. Venus Earth Mars . T[H] or T[D] (K) H D H D H D 0 9.73 6.91 8.38 5.95 27.62 22.71 100 9.91 6.97 8.62 6.01 28.36 23.02 200 10.16 7.04 8.84 6.06 29.54 23.37 300 10.41 7.12 15.0 9.07 32.00 23.82 Table 3.27: Escape fraction for hot H and hot D from elastic collisions with 0 . Columns marked H are results for hot product hydrogen, D are results for hot product deuterium. The pre-collisional hot oxygen atom is assumed to have a most probable velocity of 5.6 km/s. • Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 177 dissociative recombination reaction. 7 6 From the table it is clear that while the produc-tion of deuterium with velocities at or above the escape velocity for any of the terrestrial planets is not as efficient as for hydrogen, it is hardly negligible, especially for Mars where the escape energy for deuterium is only 0.25 eV. We agree with the assessment of Gurwell and Yung that the differential production of hot H and D on Venus is insufficient to account for the observed fractionation observed by the Pioneer Venus Orbiter. The magnitude of the D / H fractionation on Venus would require a mechanism which is much more selective in removal of hydrogen. 3.5 Summary Product velocity distribution functions describing the rate of production of hot atoms for a variety of nonthermal processes, including direct-elastic and charge-exchange colli-sions were calculated using realistic, quantum- mechanical collision cross sections. The atomic systems H - H + , D - H + , O-H, and O-D were examined. The product velocity dis-tribution functions for the various systems and nonthermal processes were compared and contrasted. The product velocity distribution functions were used to yield a direct estimate of the escaping fractions of H and D as a result of nonthermal direct elastic energization by hot oxygen atoms. These kinetic theory calculations were compared to work done by other workers. Some discrepancies were found in the magnitude and form of the collision cross sections used in this work and that of some other workers. It was found that the current calculations, incorporating quantum mechanical cross sections, were in reasonably good agreement with predictions of other workers based on anisotropic cross sections. They agreed with more the more recent works in that the fraction of hot deuterium produced via direct energization by hot oxygen is not negligible, though less than the escape fraction Chapter 3. Nonthermal Production of Energetic Hydrogen and Deuterium 178 of hot hydrogen produced in the same way. In the current calculations, no fitting was required in the extraction of the escape fractions, nor was it necessary to fix the relative impact energy or approximate the differential collision cross section using functional fits. Energy exchange rate coefficients were calculated using both a direct hydrodynamic approach and integration of the product velocity distribution functions calculated previ-ously. The results of these calculations for H - H + were compared with results obtained using a more complicated Monte Carlo approach. It was found that while roughly simi-lar in magnitude, the current approach correctly predicted the overestimation of energy exchange by the hot ion in the linear-trajectory approximation to the charge-exchange process. A similar calculation of charge exchange rate coefficients was carried out for D - H + , and it was found that the current results were in good agreement with the Monte Carlo results of other workers, but required no fitting or smoothing of statistical fluctu-ations. Chapter 4 K i n e t i c Theory Calculat ions of Nonthermal Escape Fluxes The exospheric escape problem has been examined by many previous workers for a variety of different species and for a large number of different planetary bodies. The standard collisionless approach 9 , 1 6 to the escape of atmospheric species assumes that an equilibrium Maxwellian velocity distribution function exists at the exobase. This distri-bution is integrated over all upward directed velocities of magnitude equal to or greater than the escape speed, as illustrated in equation (1.2.8), in order to yield the escape flux. The result is the well known Jeans' flux, given by equation (1.2.9). Other similar approaches to the calculation of the escape flux include examinations of the steady-state, collisionless solution (6f /6t = 0, J = 0) of the Boltzmann equation, reviewed by Chamberlain 1 6 and others.9'1 The problem of non-Maxwellian effects associated with thermal escape was examined by Lindenfeld and Shizgal . 2 6 ' 1 5 2 They studied the escape of hydrogen from Earth, and found that the perturbation of the escaping minor constituents distribution function leads to an escape flux somewhat less than that predicted by Jeans theory. The steady-state case of the Boltzmann equation, where the gravitational term is explicitly ignored (but incorporated as a boundary condition) was examined by Shizgal and Blackmore. 2 7 They used a plane-parallel model for the atmosphere, and performed an altitude dependent kinetic theory calculation for escape of hydrogen and helium from Earth, and hydrogen escape from Mars, and found good agreement with Monte Carlo calculations. Lindenfeld and Shizgal 4 5 employed a simple collisional model to obtain an expression of the charge 179 Chapter 4. Kinetic Theory Calculations of Nonthermal Escape Fluxes 180 exchange induced flux that can be compared with the thermal escape flux. They cal-culated a charge exchange rate coefficient for Earth with a value of 4 . 8 x l 0 _ 6 s - 1 , close to the estimate of 5 x l 0 - 6 s - 1 obtained by Bertaux, 4 2 and showed that the total of the charge exchange and Jeans' hydrogen escape flux was roughly constant. This was con-sistent with the results of Liu and Donahue, 3 8 - 4 0 who demonstrated that the escape flux should be approximately 1.8xl08 c m _ 2 s - 1 . A more realistic model which includes the actual density and temperature profiles was carried out by Maher and Tinsley. 7 7 Their estimates of hydrogen escape fluxes for low and middle latitudes on Earth are 1.5xl08 c m - 2 s - 1 . A comprehensive assessment of terrestrial hydrogen and deuterium fluxes was carried out by Yung et a l . 7 8 They used a one-dimensional photochemical model extending from the middle atmosphere to the exobase to model sources and sinks of both hydrogen and deuterium. They calculated net and charge-exchange escape fluxes for hydrogen of 3.02xl0 8 c m _ 2 s _ 1 and 1.41xl08 c m - 2 s - 1 . The same calculations for deuterium yielded values of 3.5xl0 4 c m - 2 s - 1 and 5.4xl0 3 c m - 2 s - 1 , respectively. It has been estimated 3 7 that the average escape flux for proton charge-exchange pro-duced hot hydrogen on Venus is 1.7xl07 c m - 2 sec - 1 , although it may have been as much as 1000 times larger in an early atmosphere much richer in hydrogen. Deuterium escape is thought to be approximately one-tenth as efficient. McElroy et a l . 7 9 ' 3 5 examined the nonthermal energization of hydrogen and deuterium via the dissociative recombination processes involving using a radiative transfer formulation. They found fluxes of ap-proximately 8xl0 6 c m - 2 sec - 1 for hydrogen, while concluding that the analogous process for deuterium was negligible. Kumar et a l . 1 5 3 estimated that the hydrogen escape flux due to charge exchange with hot protons was 1.2xl06 c m - 2 sec - 1 . A Monte Carlo simulation of hydrogen escape on Venus by Hodges and Tinsley 1 5 4 estimated the planetary average escape rate due to charge exchange was 2.8xl0 7 c m - 2 s - 1 . Shizgal 7 1 performed a kinetic Chapter 4. Kinetic Theory Calculations of Nonthermal Escape Fluxes 181 theory calculation of escape rates for nonthermal hydrogen on Venus using quantum colli-sion cross sections for hydrogen-proton charge exchange collisions. He found escape rates of 1.1 and 6.2 x 10 7 cm~ 2 s _ 1 for the dayside and nightside exosphere. These values were smaller than those reported by Kumar et a l . , 1 2 3 who performed a simpler rate calculation to find escape fluxes of 2.7 and 13.3 x 107 cm~ 2 s _ 1 , respectively. Calculations of hydrogen escape for Mars are few, as temperature and density data are poor for hydrogen. Using Mariner observations, Barth et a l . 1 2 1 and Anderson 1 2 0 es-timate a value between (1-2) x 108 cm~ 2 s - 1 . Lammer and Bauer 8 employed a Monte Carlo technique to simulate the thermalization and transport of hot oxygen atoms pro-duced by dissociative recombination of on Mars. They estimated an escape flux of 6x 10 6 c m _ 2 s _ 1 , corresponding to a mass loss of oxygen atoms at a rate of 0.14 kg/s. Other estimates of the oxygen escape rate are given by Fox 8 0 at 3x 10 6 cm- 2 s - \ and by McElroy 8 1 at (6-7) x 106 c m ^ s " 1 . In Chapter 3, we considered the result of the nonthermal charge exchange process, employing realistic differential cross sections. Others 7 1 , 7 6 have also considered the pro-duction of hot atoms from charge exchange processes. However, it is notable that dis-tributions of hot atoms calculated in this way for Venus do not yield temperatures in tremendous agreement with the observed two temperature exosphere.71 A n important consideration may involve the partial thermalization of the energetic atoms produced owing to collisions with other ambient species, something not incorporated into the mod-els of Chapter 3. Due to thermalization, it is expected that the escape rate will be reduced, bringing the theory into better agreement with observations. To account for the competition between thermalization and reactive processes in a rigorous fashion re-quires an altitude dependent kinetic theory treatment of the Boltzmann equation, and is currently beyond the scope of this work. However, it is possible to make approximations of the thermalization process in order to yield escape fluxes which reflect the effect of Chapter 4. Kinetic Theory Calculations of Nonthermal Escape Fluxes 182 the ambient neutral atmosphere on hot atom escape rates. This is addressed in further detail in the theory section of this chapter. It is also clear that the theoretical structure of the exosphere must be reconsidered in order to include the effect of escape due to collisional nonthermal processes. The collision-less picture is not a true reflection of the actual exosphere, although one could consider the exosphere as almost collisionless with infrequent nonthermal collisional processes as a first order deviation from the collisionless state. 1 2 2 ' 5 However, in order to properly de-scribe the transition between the collision dominated and collisionless regimes properly, the notion of a sharp discontinuity at the exobase that divides the atmosphere between hydrodynamic kinetic regimes has to abandoned. Instead, the escape process should be considered as occurring from a range of altitudes above and below the exobase.2 6'1 This concept originates in the older works of Biutner, 1 5 5 Jensen 1 5 6 and Jockers, 1 5 7 and has been considered by Shizgal and Lindenfeld 2 6 ' 4 5 and Shizgal 7 1 and more recently by Johnson. 1 5 8 ' 1 5 9 It forms the basis for the treatment of the escape process used in this chapter. 4.1 Theory The basic concept of the model proposed by Lindenfeld and Shizgal 2 6 is that transla-tionally energetic species are produced by elastic or reactive collisions. These energetic particles are moving in all directions but only those moving radially outward can escape, provided they do not suffer further collisions that change the direction and magnitude of their velocity so as to prevent it. Since the density decreases approximately baro-metrically (see equation (1.2.4)), the rate of collisions decreases with increasing altitude, and hence the population of translationally energized particles decreases. However, with the decrease in the density the mean free path increases and the probability of escape Chapter 4. Kinetic Theory Calculations of Nonthermal Escape Fluxes 183 increases. Thus, the rate of escape as a function of altitude attains a maximum value occurring in the vicinity of the exobase. The net escape flux may be written as an integral over radial position, The units of Fnet are those of a flux, c m - 2 s _ 1 . In equation (4.1.1), ro is some reference level deep enough in the atmosphere that the escape probability is zero while ri is the altitude above which the atmosphere can be considered 'collisionless' and there is no collisional production of hot, nonthermal atoms. F(r), the rate of production of energetic atoms per unit volume escaping from radial position r, is given by where Q ( r ? c i ) is the velocity dependent production rate of energetic particles at radial position r and p(r, ci), is the probability that a particle with speed c\ at position r will escape. The explicit form of p(r, Ci) exhibits a very weak dependence on the speed.2 6 The probability for escape, p(r,c\), increases to approximately 1/2 with altitude while Q(r, Ci), the production of energetic particles through collisions, decreases with altitude. The units of F(r) are c m - 3 s _ 1 . It is assumed that Fnet represents the net escape flux from the entire range of altitudes, but only along the direction dr. If it is further assumed that the flux from the atmosphere or system is isotropic, then the overall planetary escape rate is found by multiplying Fnet by the area of the spherical shell defining the uppermost boundary of the system. We consider a system with two species. Particles of type 1 are the escaping, colli-sionally energized species, while particles of type 2 comprise the (thermal) background. Collisions occur only between the two different types; non-linear self collisional types of processes are not considered. If we choose a 'plane parallel' atmosphere, with geometry as shown in Figure 4.51, the probability of escape of a particle with speed Ci, moving in a (4.1.1) (4.1.2) Chapter 4. Kinetic Theory Calculations of Nonthermal Escape Fluxes 184 direction making an angle 9 with the zenith direction and which suffered its last collision at an altitude r is given b y 1 6 0 , 2 6 where Z(ci) is the dimensionless collision frequency for particles of species 1 with the where A is the non-dimensionalizing factor for Z(c\), g is the relative velocity, fl is the solid angle about the scattering direction, and a is the differential scattering cross section. If the total elastic collision cross section for collisions between the two particle types is assumed to be the hard sphere cross section, nd2, where d is the mean diameter of the colliding pair, d = (d\ + d2)/2, it can be shown (see Appendix F) that the collision frequency is given by (4.1.3) background, n^r') is the density of the background gas, and the range of 6 for escaping particles is taken as 0 < 0 < TT/2. The velocity dependence of equation (4.1.3) arises through the collision frequency, Z(c\), which is defined as (4.1.4) where A = d2 (^p -J 2 , and where we have defined a dimensionless energy x = the mass ratio 7 = m 2 / m i , and used the standard definition of the error function erf(x). The escape probability p(r, Ci) from equation (4.1.2) is given by the average of P(r, cx,Q) over all possible escape angles 6, that is (4.1.6) Chapter 4. Kinetic Theory Calculations of Nonthermal Escape Fluxes 185 where we have defined dt /CO e-stC and where y is the vertical optical depth, y = Trd2J\2(r')dr' (4.1.7) Equation (3.4.1) gives the rate of production of particles of type 1, with speed ci , due to collisions. Only those particles with speeds exceeding the escape speed will per-manently escape the planetary atmosphere. Thus, in calculating the net escape flux, equation (4.1.1), we want to include only those particles of type 1 which have speeds exceeding the escape speed. Thus, the rate of production of hot atoms with speeds exceeding the escape speed becomes26 Q(r,d) = H{Cl - C e , c )n i ( r )n 2 ( r ) Jj / i ( r , c ' ^ r , c'2)ga dft dc2 (4.1.8) H(C\ — Cesc) where 0 C\ < Cesc 1 Ci > C e sc is simply the standard Heavyside function. If it is assumed that the velocity distribution functions of both the background and the escaping species have reached equilibrium, we get the equilibrium rate of production of hot atoms, Qeq(r,Cl) = H(c1-Ce8C)n1(r)n2{r)JJf»f( 0; that is, they have translational energy exceeding the gravitational potential 'well' of the planet. Those that are non-escaping have E < 0, with the value E = 0 marking the division between the two classes. If we define Xt = Xsm9 Xr = X cos 9 we see from equation (A.0.1) that the boundary E — 0 is equivalent to X2 = X2 + X2 = y or a circle of radius ^/y in (Xr,Xt) phase space. Particles with (Xr,Xt) inside this circle execute closed, elliptical orbits while particles with (Xr,Xt) outside the circle execute open hyperbolic ones. The other important concept in the division of particle classes is whether or not a particle trajectory intersects the exobase level, r c . A particle whose trajectory touches the exobase at only one point (a "glancing" trajectory) will have 9C = 7r/2; that is, from equation (A.0.2), Xc is perpendicular to r. Any trajectory crossing r c must have 9C ^ ir/2. If we substitue 9C = ir/2 into equation (A.0.2), square it, and then substitute for X2 from equation (A.0.1), we have after some simplification that X ' ~ ] 4 ^ X r = xl (A.0.3) Appendix A. Density Profiles in a Collisionless Exosphere 201 where 1 + 2/ Equation (A.0.3) defines the particle trajectories which do not cross the exobase; it is of the form a2 b2 which defines a hyperbola with foci at c = a2 + b2 = ±Xb/y and eccentricity e = c/a = 1/y. If we define the polar angle of the hyperbola as 0P, we have for any point on the hyperbola that Xt = X sin 8P XT = X cos 0P or, substituting into equation (A.0.3) and simplifying, smOp = ^-{X2 + l-y)12 (A.0.4) Ji. The various classes of particles are listed in Table A.28, and may be graphically illustrated for fixed values of y, as in Figure A.56 (which corresponds to a radial position giving y = 1.5). The densities of particular classes of particle at altitudes above the exobase are cal-culated by assuming that the distribution function at the exobase is a Maxwellian. From equation (1.2.7), we have = / fM(cc)dc J partial / m. \ " — 2-jrn TIL ,c partial mc2 2kTc c2 dc sin 0 dO (A.0.5) Appendix A. Density Profiles in a Collisionless Exosphere 202 i b M j p i n j : ri\- / \ \ by / \X$C* Y Satellite Xt 1 Satellite A, Ballistic / \ Fly-by Capture Xr Figure A.56: A graphical illustration of exospheric particle classes in the phase plane of (Xr,Xt). The altitude used for this particular figure corresponds to a value of y = 1.5. The quantities y, Xb and 0P are defined in the text of this Appendix. The escaping and ballistic particle regions are shaded to clarify the regions of phase space reprsented by those particle classes. Based on figures by Fahr and Shizgal.1 Appendix A. Density Profiles in a Collisionless Exosphere 203 Particle Class Orbit X limits limits Region of (Xr, Xt) (B)allistic k (S)atellite elliptic (b) 0^y 0 < 0 < 7T outside large circle (B)allistic elliptic (b) and cross rc 0 < X < Xb Xb^y o < e < ep outside large circle k inside hyperbola Table A.28: Exospheric particle classes. The incoming particles are sometimes divided into a fly by and capture component. The abbreviations 'b' and 'nb' under the 'Orbit' column refer to bound or non-bound orbits of the type specified. Based on Fahr and Shizgal. 1 If we use the definition of the escape parameter from equation (1.2.10), E \ -'-'esc c = TFe GMm kTcrc change to dimensionless variables y,X as previously, and use equation (A.0.1), equation (A.0.5) becomes 3 nk = 2Trnbar(—)[[ e~XcXX2 dX sin0 + Qy/x~) Zi V cc ±[exi{Q\/x~> - Ry/x~) + ex~x' erf (Ry/x1 - Qy/x~)}} where the +(—) sign refers to x > x' (x < x'), and (B.2.1) Q = 1 2 R = 1 2 7 = mback/mhot A = hs SkTback V Kmhot This is multiplied by the symmetrizing factor Jx (x-x')/2 y^Le-(x-x<) so that Ks(x,x') = Fs-Kns(x,x') Appendix B. Details from Chapter 2 209 or Ka{x,x') = -AQ2J^L={e(x'-x)l2 [erfYQVz7 + Ry/x) ± e r ^ Q v x 7 - Ry/x)] 2 y \Jx' x + e (*-* ')/2 [ e r f(i2v^ + Qy/x~) ± ed(Ry/x1 - Qy/x~)]} (B.2.2) is symmetric in its arguments. The +(—) sign refers to x > x' (x < x'), or in terms of the magnitude and order of the arguments, +(—0) where the first argument is larger (smaller) than the second. If Ks(x,x') is symmetric in its arguments then Ks{x,x') = Ks(x',x) Exchanging the arguments of equation (B.2.2) yields K'(x',x) = l-AQ2.[-^={eix-x')/2 [erf(Qv^ + RVx') ± ed{Qy/x - Ry/x')] 2 y \Jx x' + e (* '-*)/2 [erf'R^ + Q^i) ± erf(Ry/x~ - QVx1)]} (B.2.3) If we recognize that the error function is anti-symmetric in its argument, i.e. erf(—z) = —erf ( 2 ) , and then rearrange equation (B.2.3), we have K'(x', x) = I ^ Q 2 . G={e(*'-*>/2 [ed^Qy/x1 + Ry/x) ± -er{(Qy/x~' - Ry/x)] 2 y \Jx x' + e ( * - * ' ) / 2 [erf(i2vV + Qy/x~) ± -ed{Ry/x~' - Qy/x~)]} (B.2.4) where +(—) sign refers to x < x' (x > x'), or in terms of the magnitude and order of the arguments, +(—) where the second argument is smaller (larger) than the first This is exactly opposite to what we had in equation (B.2.2). To synchronize the sign convention between the two equations we choose the convention that the +(—) sign refers to x > x' (x < x') irrespective of the actual order of the arguments. With this convention equation (B.2.4) becomes K'(x',x) = \AQ2J-^={e^x'-^12 [eri(Qy/x~' + Ry/x) ± erf(Qvx7 - Ry/x)] 2 y yjx x' + c (*-* ')/2 [erf(i?V^+ Qs/x") ± eri{Ry/x~' - Qy/x')]} (B.2.5) Appendix B. Details from Chapter 2 210 and comparison with equation (B.2.2) verifies that Ks as defined in equation (B.2.2) is indeed symmetric in its arguments. The non-symmetric kernel Kns may thus be replaced by Kns(x,x') = 1 \fx or, re-written in terms of speed quadrature points, X e(*-*') Ks(x,x') KnS(ylyl) = Jy-e^)K°{ylyl) v y> B .3 Symmetrization of the Discretized Collision Operator From Chapter 2, we had the equation Yt = M-f + S (B.3.6) We first define a dimensionless 'time' T = At where A is the quantity defined in equation (2.4.11), yielding ^ = M-f + S (B.3.7) OT where M and S are simply M and S divided by A. The eigenvalue problem corresponding to the equation (B.3.7) is MV»(n) = A„V(n) (B.3.8) or /"CO / dx M(x,x')^n\x) = Xn^n\x') (B.3.9) J 0 Appendix B. Details from Chapter 2 211 where A /"CO M = dx M(x,x') Jo M(x,x') EE [Kns(x,x')-Z{x')6(x-x')]/A Our goal is to transform this eigenvalue problem into one of the form B^n> = \n{n) (B.3.10) where B is a real, symmetric matrix. We know the symmetrizing factor for the kernel from equation (2.4.12), and so we define \ £ 2 f - e * - * ' B(x,x') M(x,x') and then substitute this into equation (B.3.9) to yield rOO Jo d X \ v ' cx—x' {^e*-*' B(x,x')il>W(x) = \ni>{n)(x') 2 We then transform to a quadrature in speed, with x' = y2 and x = y2, 3=1 V y i or N 3=1 \ Wie~y'. Rearranging, we finally have N I 3 = 1 WJ . / » '-51 *!•> e y> where (B.3.11) (B.3.12) BtJ = 2y/ywWiW3 Bjti (B.3.13) Appendix B. Details from Chapter 2 212 is, by inspection, symmetric in (since B is a symmetric matrix). Comparing equa-tions (B.3.10) and (B.3.12), we see that *!B) = V^^B) (B.3.14) If we redefine our original eigenvectors according to equation (B.3.14), / Wi U = ^f, (B.3.15) we find after some manipulation that equation (B.3.7) becomes df * . ~ = B-f + S (B.3.16) where B is a symmetric matrix defined by equations (B.3.11) and (B.3.13), and V e~y2 B.4 Solution of the Eigenvalue Problem From equation (2.4.14), we had = B-f + S (B.4.17) df dr This is analagous to the eigenvalue problem of the f o r m 1 6 2 , 1 6 3 BU = A -U (B.4.18) where B is a symmetric matrix and U is the matrix of eigenvectors of B. If we make the definition f = U F in equation (B.4.17), we find that Appendix B. Details from Chapter 2 213 or after multiplication by U 1 from the right-hand side OF ^ ^ „ — = D - F + Q OT (B.4.19) where Q = U _ 1 S D = U _ 1 B U The matrix D is a diagonal matrix of the eigenvalues of B. It can be shown that all the eigenvalues are inherently negative, except for A i which we expect to be zero since it represents the equilibrium eigenfunction. We choose to redefine the eigenvalues with the negative sign explicitly shown, i.e. A = — A. We then have Ai = 0 < A 2 < A 3 < . . ' . < A n With this convention we may write the jth component of equation (B.4.19) as dFj dr If we use an integrating factor, this may be rewritten as .d_ dr ex>TFj = Qj or integrating both sides dr Carrying out the integration on the right-hand side and rearranging yields "1 _ e-A>i F3(T) = FMe-^ + Qj Appendix B. Details from Chapter 2 214 N With our previous definition of f = U • F , or fj = ^VijFj, we finally arrive at 3 = 1 N Mr) = Uitl [*i(0) + Q a r] + [ Fj(0) + Q, 3=2 1 - e~A> .A,-where we have used L'Hopital's Rule to account for the indeterminacy caused by the zero eigenvalue in the j = 1 term. Appendix C The Differential Collision Cross Section The scattering cross section is of great importance in any detailed analysis involving the Boltzmann equation. It appears explicitly in the definition of the collision operator in equation (1.4.1), and describes the physics of the binary collisions for the system. The physical meaning of the cross section is made clear by considering a typical scattering configuration, such as is illustrated in Figure C.57, where a uniform beam of particles is incident on a (repulsive) scattering cemter. For the particular case of a two-body collision, both the incident beam and the target become single particles, but the description of the process remains the same. The incident beam is characterized by its intensity 7, also called the incident flux density, which gives the number of particles per unit time which cross a unit area per-pendicular to the incident direction. There is a finite region of interaction over which the incident beam is deflected (or scattered) from its incident direction after which the effect of the scatterer diminishes and the scattered beam continues linearly in the scattered direction. Referring again to Figure C.57, we may quantify this process as fol lows. 1 3 2 ' 1 3 3 For a central force type of potential (i.e., dependent only on radial distance), there is symmetry about the axis of the incident beam, and the element of solid angle may be written as dii = 2TT sin(0) d9 (C.0.1) The angle 0 is the angle between the incident and scattered directions, and is not supris-ingly known as the scattering angle. For any given incident particle, the amount of 215 Appendix C. The Differential Collision Cross Section 216 Figure C.57: A typical scattering configuration. An incident beam of particles interacts with a scattering center, and part of the beam is scattered through an angle 9 into a solid angle dQ. scattering is determined by the impact parameter b (the perpendicular distance between the direction of incidence and the scattering body) and the energy E of the particle. If the incident intensity or flux is given as 7, the number of particles per second scattered into the solid angle dQ (that is, between 6,9 + d9) will be proportional to both 7 and dQ, AN = aldQ (C.0.2) where a is a constant of proportionality, and is called the differential cross section. We assume that different values of b cannot lead to the same scattering angle (that is, that 9 is single valued function of b). Thus, the particles scattered into 9,9 + d9 are those incident particles with an impact parameter between b, b + db, AN = 2%Ibdb (C.0.3) Equating equations (C.0.2) and (C.0.3), and making use of the definition of the solid Appendix C. The Differential Collision Cross Section 217 angle from equation (C.0.1), we find that the differential cross section is given by \db o = « » ( » ) ' • » ' ( C ' a 4 ) The absolute value signs are used in equation (CO.4) because b and 6 may vary in opposite directions despite the fact that particle number A i V must always be positive. Physically, the differential cross section is the collisional obstacle or area presented by the scatterer to the incident beam. The term 'cross section' is appropriate since the dimensions of o are length-squared or area. The total scattering cross section is simply the integral of equation (CO.4) over all directions, &tot — Jodtl = 2TT r o sin(0) dd (C.0.5) Jo A special form of the differential cross section that is regularly used in kinetic the-ory is the hard-sphere cross section, often also referred to as the "billiard ball" model. It approximates each scattering particle as a rigid sphere of a fixed radius r0 with an interaction potential defined as constant > 2r0 V(r) = I oo < 2rG The hard-sphere differential cross section arising from this potential is adiff (2r 0) 2 4 A n important feature of the hard-sphere differential cross section is its independence of both energy and angle. The total cross section is found by integrating over all solid angles, °tot ~ J adiffdil r sin(0)i Jo = 2-K J m(e)de = ird2 Appendix C. The Differential Collision Cross Section 218 where d = 2r 0, the hard-sphere diameter. Though not a particularly realistic model of particle interaction, the hard-sphere cross section is useful. It is often used to reduce the complexity of the collision operator in the Boltzmann equation, and is also useful in the examination (analytically) of cross section dependent properties that would be extremely difficult or impossible to calculate with a more realistic model. Finally, it can be useful as a test of computer codes, since there are many results for a large class of problems utilizing the hard-sphere cross section. Appendix D Scattering Theory D. l Quantum Mechanical Scattering The description of the scattering (or collision) process is well described in most stan-dard texts on quantum mechanics. 1 3 2 - 1 3 4 In this section we give an outline of the theory describing Q M scattering, and the procedure for the calculation of the phase shifts and hence the differential and total elastic collision cross sections. We are interested in finding the wave function ij>(r) for a particle with a given en-ergy E moving in a potential field V(r). This problem involves solving the Schrodinger equation, 1 3 2 ' 1 3 3 HV>(r) = E^[v) (D.l.l) where H , g + is the Hamiltonian operator and represents the total energy. If we substitute for the linear momentum operator, p = — ihW, and use the definitions k ~ A A s * P p = kh, 219 Appendix D. Scattering Theory 220 1 2' E = ~-pg2 = p2 h2k 2 1,2 2p 2p the Schrodinger equation can be rewritten as 2p •V2 + V(r) %2k2 (D.1.2) In equation (D.1.2), p is the reduced mass and k is the wave number. If we assume that the potential V(r) is spherical, and define U(r) = pV(r) (D.1.3) then equation (D.1.2) becomes W2 + k2-U(r) iP(r) = 0 (D.1.4) At large radial distances we expect the potential V(r) to fall off to zero, and the effect of a scattering body on another particle is assumed to be negligible. If the potential term of equation (D.1.4) is set to zero we have the 'free-motion' Schrodinger equation V2 + k2 V(r) = 0 (D.1.5) A solution of equation (D.1.5) is ^in(r) = Ae >kr (D.1.6) which is a plane wave (that is, k • r, the phase of the wave, is a constant). We expect that in the absence of a scatterer that the wave function would be of this form, and that the inclusion of a scattering potential simply adds an additional component to the incident wave function. We approximate the effect of a collision by examining the change in the incident particle wave function after interaction with the potential field of the central scattering body. This is similar to the treatment of diffraction effects in optics or Appendix D. Scattering Theory 221 acoustics. 1 3 2 As the incident particle, in the form of a plane wave as given by equation (D.l.6), approaches the scatterer it begins to feel the potential field of the scatterer. We assume that scattering from a point target produces a spherical wave of the form eikr ^c(r) = F{B)-— (D.1.7) r where the J-(9) is the scattering angle dependent scattering amplitude, and allows for angular anisotropy in the scattered wave function. The total scattered wave function is a superposition of the wave functions in equations (D.l.6) and (D.1.7), ^ ( r ) = A e ' k - r + j r ( 0 ) - — (D.1.8) r Since the scattered spherical wave decays as 1/r as we move away from the collision interaction region, the total wave function asymptotically approaches a plane wave again. The effect of a the scattering body on the incident wave function is illustrated by Figure D.58. We may now relate the quantum mechanical wave function terminology with the classical derivation of the collision cross section given in Appendix C. The intensity of the incident beam is defined as = | A | 2 ^ (D.1.9) and that of the scattered beam as he = \ipsc\ 2g = ( D . U 0 ) where g = p/p is the relative velocity of collision. The number of particles per second scattered through an angle 6 into an area dS (perpendicular to the scattered beam) is AN = LJS Appendix D. Scattering Theory 222 Figure D.58: Wave function scattering by a central potential. A n incoming plane wave is partially transmitted/reflected as another plane wave and partially scattered as a spherical wave by the potential field due to the scattering body. = mt*JidS (D. i . i i ) H fj. If the incident beam is normalized such that \A\2 = 1, then Iin = — (D.1.12) If we put equation (D.1.12) into equation ( D . l . l l ) , and use the definition that dS/r2 = dfl, we find that AN = \F{6)\2 IindVl (D.1.13) Comparing this with the form of equation (C.0.2), we immediately find that the quantum mechanical differential cross section is given by a = \f(0)\2 (D.1.14) Appendix D. Scattering Theory 223 D.2 Calculation of the Phase Shifts To calculate the scattering amplitude J7, we must match equation (D.l.8) to the asymptotic form of the solution of the Schrodinger equation including the interaction potential, equation (D.l.4), which was written as [V2 + k2 - U(r)] ii>(r) = 0 The solution for the wave function may be written as an expansion of the form oo tf(r) = £ C W , M ) z=o where we assume that 1>i(r,8) = Ri(r)Pi(coaO) The /th term in the expansion of ip is known as the /th partial wave. In general, the solu-tion of equation (D.l.4) is a three dimensional problem, where ij) = ij)(r,0, as one looks along the incident beam axis), and so the problem reduces to a two dimensional one. We may thus write }_d_ 28Ri Pl + ^ L 2 P , r 2 r 2 dr \ dr where L is the quantum mechanical angular momentum operator 1 3 3 with the = ^ H w ) = -i{i+i)p< If this is substituted into equation (D.l.4), the original problem reduces to determining the i?/'s from the radial Schrodinger equation, 1 d ( 23Ri r r2 dr \ dr ) + ( f c2 - c / e / / ( r ) ) / i , ( r ) = 0 (D.2.15) Appendix D. Scattering Theory 224 where we have defined an 'effective potential' Ueff = U(r) + If we consider the'free motion' (U(r) = 0) version of equation (D.2.15), we may rewrite it in the form of Bessel's equation of half-integral order, r ' | £ + 2 r ^ + ( * V + .-(/ + l ) ) « ( r ) = 0 The general solution is of the form of a linear combination of two independent solutions, Ri(r) = Cnji(kr) + C2mi(kr) where j\ is the spherical bessel function, n; is the spherical Neumann function, and Cu and C21 are coefficients to be determined by boundary conditions. Since the spherical Neumann functions diverge at the origin, and we require that Ri be finite everywhere, we set C21 = 0. Our 'free motion' solution is thus = Cuji(kr)Pi(cose) (D.2.16) It can be shown 1 3 4 that this is exactly equivalent to an expansion in Legendre polynomials of our plane wave 'free motion' solution of equation (D.1.6). Since we are interested in the form of tj) as r —• oo (that is, far from the interaction region), we examine the asymptotic behaviour of the radial portion of equation (D.2.16), ,-w \ ^ sm(kr — lir/2) / T ^ ^ , „ x #,(»•-> 00) = Cu—-—: — (D.2.17) kr We now consider the case where the potential in equation (D.2.15) is not zero. If we restrict ourselves to the region far from the interaction region, r >^ 0, and assume Appendix D. Scattering Theory 225 that the potential U(r) is negligible in this region, then we find that equation (D.2.15) again reduces to Bessel's equation of half-integral order. The general solution is a linear combination of spherical Bessel and Neumann functions, as previously. However, we cannot simply exclude n/ as before since the origin is no longer included in the region of interest. Setting Cu — 7; cos Si Cii — —7; sin Si we have the solution Ri(r) = 7; [ji(kr) cos Si — nj(fcr) sin S{\ (D.2.18) For r —• oo, this takes the form sin(fcr - /7r/2 + SA Rt(r - c o ) = 7 i and thus the asymptotic wave function is given by sm(kr — fa/2 + Si) „, „s ,^ n ^s jjji = 7 ; — i ' T u PticosO) (D.2.19) kr Comparing equations (D.2.17) and (D.2.19), we see that the inclusion of the potential causes a small shift Si in the phase of the radial part of the asymptotic wave function. This phase difference between asymptotic solutions of the radial wave equation with and without the potential is known as the phase shift. In general, the phase shifts are calculated by finding the nodes r0 of Ri from equation (D.2.18), yielding tan 6, = ^rH (D-2.20) ni(kr0) In order to finally find the scattering amplitude we now match the expected form of the asymptotic wave function, equation (D.l.8), to equation (D.2.19) . 1 3 3 ' 1 3 4 If we express Appendix D. Scattering Theory 226 the sine function as a sum of exponentials, sin z = -(eiz-e~iz) 2^v ; this yields Cu - - • • " i k r S 27^ [eikre~ilT/2 ~ e-ikreil*'2] P,(cos 9) + F{9) • &— = OO 2ib [JkTz~iiV'2ziSl ~ e-tkreil*/2e-iS'} P,(cos 9) 1=0 oo 1=0 Since e~lkr and e , f c r are linearly independent, the coefficients of these terms must sepa-rately be equal. For the terms involving e~lkr, we have oo oo Y,CuPi(coS9) = £ 7 , P , ( c o s 0 ) e - i 5 < f=0 1=0 If we multiply by P;»(cosfl), integrate over d(cos9), and use the orthogonality condition for the Legendre polynomials, [+i P,Pi,d(coS6) = ^-rS„, J-i 21 + 1 we find 7/ = CueiS> (D.2.21) If we make use of the fact that the asymptotic expansion of the plane wave in equation (D.1.6) is given by equation (D.2.16), we have ^ = &n(r) = e i k r = J2C^Ji(kr)Pi(cos9) i If we multiply by P//(cosfl) and integrate over d(cos0), again using the orthogonality condition for Legendre polynomials, and then make use of the definition for the spherical Appendix D. Scattering Theory 227 Bessel functions, we have Clljl^2J+T6u' = / - i eikrcos9M^se)d(cos9) = 2ilji(kr) or Cu = (2l + l)il (D.2.22) We are now (finally!) ready to solve for the scattering amplitude itself. Equating the coefficients of the e%kr terms yields oo oo Cu P/(cos 0)e-il7r'2 + 2ikF{6) = £ 7, P,(cos 6) e-'(W2-««) If we substitute for Cu and 7/ from equations (D.2.22) and (D.2.21), and simplify as for the e~tkr case, we find 1 00 HO) = ^ £ ( 2 / + l)P,(cos0) 1 After taking a factor elSl out of the brackets and using the definition of the sine function, we finally arrive at the form F(9) = - £ ( 2 / +l)P,(cos0)e i 5 z2((a'2 + f3'2y ra'+P' Ja'-P 1 exp -iz2 2a! B' ' u2 \ 2a>8>) x smh(xzu) du (E.5.8) This is of the form e~Au2 s i n h ( 5 « ) du (u = J e~Au2e-Budu Kt - K: (E.5.9) If we complete the square in the argument of the exponential of equation (E.5.9) we have for the first term Kt = \eB2'iA j e-(^+Bu+By4A) d u Appendix E. Details of the Production of Hot Atoms 235 where the final form follows directly from integral tables. 1 5 0 Similarly, B \ K: \fAe^eri(u^A- 2VX From equations (E.5.8) and (E.5.9) we may identify A B B2/AA B/2y/A uVA jz2 2a'8' xz a''F . = x a'8' ~2~y UZ\ 2a'(3' Putting it all together, we thus may write equation (E.5.8) as h{z) = N8 K{yjA1 + A2) - Kiy/A! - A2) where we have collected coefficients and defined K(u) 1 1 M r _eC2e-y'z2A1/A2 u = 2z*]l A2l> exi(uz'+ C) - eri{uz - C) u c 4 y Substituting equation (E.5.10) into (E.5.6), h = \{^'eC2 L z e -(~i'AilA2-B>)zi K{sJAx + A2) - K(\/AI - A2) dz (E.5.10) (E.5.11) (E.5.12) Appendix E. Details of the Production of Hot Atoms 236 With equation (E.5.11) giving the form of K(u), equation (E.5.12) thus reduces to the form (ignoring factors in front of the integrals) where I(jAx + A2) - I(y/A1 - A2) f°° -? 2 I(ri±) = / yz [eri{n±y + 8)- eri(r)±y - 8)] dy Jo (E.5.13) and V = < y/Ax + A2 = a' + 0' for + y/Ax - A2 = ct'-B' for -Since the two terms in I(r]±) differ only by the sign of the argument of the erf() function, we really only have to do the integral = / ye erf(r/y + 8)dy Jo We do this integral by parts; let u = evi(ny + <5), dv = ye ^ so du = ^ 7 ? e ^ w + ^ 2 dy and v = - e " ^ 2 / 2 ^ , and thus h(r}) = -—-ev{(ny + 8) +-= exp - ((£ + vV + Wy + dg c\ = c2 • c2 = g2 + cl + 2gCl cos 8 (F.0.4) Appendix F. Derivation of the Hard Sphere Collision Frequency 242 Making the change of variable of integration from c2 to g, and using the results in equation (F.0.4), we now have v{cx) = N2nd2JJJexp -\(g2 + c\ + 2gcx cos 9) vth g3 sin 9 d9 d dg Since there is azimuthal symmetry, the integration over (j> simply yields a factor of 2 i r . If we make the change of variable p = cos we have u(ci) = 2N2n2d2 f Jo exp vth g / exp 7 Jth {2gc1[i) dp dg The integral over p is simply an exponential and may be done easily, yielding i/(d) 2N2TT2d2v2h r°° 2 7 C l Jo exp exp 2*/gci Jth J — exp -2jgcx Jth 92dg If we make the change of variable from velocity to reduced speed, V = g/vth we then have v{x) = N27T2d2vfh / o o 7 £ 'th f00 T 1 — / exp [~7(x2 + y2) (exp [2~jxy] - exp [-2 7zy]) y2 dy J 0 We may rewrite this as N2ir2d2vK r~ v(x) = ^ V t h (exp [ - 7 (s - y)2] - exp [ - 7 (x + y)2}) y2 dy (F.0.5) The integral over y in equation (F.0.5) may be performed using ???, yielding v{x) N2*2d2v}h 7 £ xe 1 X ( 1 x2 1 where erf(x) is the standard error function, erf(x) -- -^= / e s 2 ds y/TtJo Appendix F. Derivation of the Hard Sphere Collision Frequency 243 Substituting back in the definition of N2, and simplifying the result in the brackets, we finally have v(x) \/^d2vth x/7 --yx + 2 \^fx + 2x^/7 \eri(xfyx) (F.0.6) If we wish instead to express the collision frequency as a function of reduced energy, x = - x mlc\l2kT2 2 / 2 Cl/Vth 2 we have (F.0.7) where we have defined A = y/ird2vth It should be noted that since the background distribution is normalized to 1, the units of the collision frequency v are not sec - 1 , but rather cm 3/sec. A dimensionless collision frequency may be readily defined as v{x") Bibliography [1] H.J . Fahr and B. Shizgal. Modern exospheric theories and their observational rele-vance. Rev. Geop. & Sp. Phys., 21:75-124, 1983. [2] D.E . Anderson. The Mariner 5 U V photometer experiment: A n analysis of H Lyman-a data. J. Geophys. Res., 81:1213-1216, 1976. [3] J .L. Bertaux and et al. Lyman-o: observations of Venera 9 and 10. i . the nonthermal hydrogen population in the exosphere of Venus. Plan. Sp. Sci., 26:817-831, 1978. [4] A.I . Steward, D.E. Anderson, L .W. Esposito, and C A . Barth. Ultraviolet spec-troscopy of Venus: Initial results from the Pioneer Venus orbiter. Science, 203:777-779, 1979. [5] R.R. Hodges and B . A . Tinsley. Charge exchange in the Venus ionosphere as the source of the hot exospheric hydrogen. J. Geophys. Res., 86:7649-7656, 1981. [6] W . H . Ip. On a hot oxygen corona of Mars. Icarus, 76:135-145, 1988. [7] A . F . Nagy and T.E . Cravens. Hot oxygen atoms in the upper atmosphere of Venus ' and Mars. Geop. Res. Lett., 15:433-435, 1988. [8] H . Lammer and S.J. Bauer. Nonthermal atmospheric escape from Mars and titan. J. Geophys. Res., 96:1819-1825, 1991. [9] J .H. Jeans. The Dynamical Theory of Gases. Cambridge University Press, fourth edition, 1925. [10] K . D . Cole. Theory of some quiet magnetospheric phenomena related to the geo-magnetic tail. Nature, 211:1385-1387, 1966. [11] J .H. Yee. A Theoretical and Experimental Study on the Atomic Oxygen Corona in the Earth's Atmosphere. PhD thesis, University of Michigan, 1980. [12] A . E . Hedin. Hot oxygen geocorona as inferred from neutral exospheric models and mass spectrometer measurements. J. Geophys. Res., 94:5523-5529, 1989. [13] C A . Barth and et al. Ultraviolet emissions observed near Venus from Mariner 5. Science, 158:1675-1678, 1967. [14] G. Kotova et al. Study of the solar wind deceleration upstream of the Martian terminator bow-shock. J. Geophys. Res., 102:2165-2173, 1997. 244 Bibliography 245 [15] B . D. Shizgal and G. Arkos. Nonthermal escape of the atmospheres of Venus, Earth, and Mars. Rev. Geophys., 34(4):483-505, 1996. [16] J .W. Chamberlain. Planetary coronae and atmospheric evaporation. Plan. Sp. Sci., 11:901-960, 1963. [17] D . M . Hunten and T . M . Donahue. Hydrogen loss from the terrestrial planets. Ann. Rev. Earth & Plan. Sci., 4:265-292, 1976. [18] B . A . Tinsley. Effects of charge exchange involving H and H + in the upper atmo-sphere. Plan. Sp. Sci., 26:847-853, 1978. [19] D . M . Hunten. Thermal and nonthermal escape mechanisms for terrestrial bodies. Plan. Sp. Sci., 8:773-783, 1982. [20] D . M . Hunten. Escape of atmospheres, ancient and modern. Icarus, 85:1-20, 1990. [21] K . K . Mahajan and J. Kar. A comparative study of Venus and Mars: Upper atmo-spheres, ionospheres, and solar wind interactions. Indian Journal of Radio & Space Phys., 19:444-465, 1990. [22] P . M . Banks and G. Kockarts. Aeronomy, Parts A & B. Academic Press, 1973. [23] T . M . Donahue. Fractionation of noble gases by thermal escape from accreting planetesimals. Icarus, 66:195-210, 1986. [24] H.B. Niemann, W.T. Kaspizak, A . E . Hedin, D . M . Hunten, and N .W. Spencer. Mass spectrometric measurements of the neutral gas composition of the thermosphere and exosphere of Venus. J. Geophys. Res., 85:7817-7827, 1980. [25] J .L. Fox and A. Dalgarno. Ionization, luminosity, and heating of the upper atmo-sphere of Mars. J. Geophys. Res., 84:7315-7333, 1979. [26] M . J . Lindenfeld and B . Shizgal. Non-maxwellian effects associated with the thermal escape of a planetary atmosphere. Plan. Sp. Sci., 27:739-751, 1979. [27] B . Shizgal and R. Blackmore. A collisional kinetic theory of a plane parallel evap-orating planetary atmosphere. Plan. Sp. Sci., 34:279-291, 1986. [28] P .K. Kundu. Fluid Mechanics. Academic Press, 1990. [29] S. Chapman and T .G . Cowling. The Mathematical Theory of Non-Uniform Gases. Cambridge at the University Press, 1964. [30] P.G. Richards, M.P. Hickey, and D.G. Torr. New sources for the hot oxygen geo-corona. Geop. Res. Lett, 21:657-660, 1994. Bibliography 246 Y . H . Yee, J.W. Meriwether, and P.B. Hays. Detection of a corona of fast oxygen atoms during solar maximum. J. Geophys. Res., 85:3396-3400, 1980. A . F . Nagy, J . K i m , and T .E . Cravens. Hot hydrogen and oxygen atoms in the upper atmospheres of Venus and Mars. Annal. Geop., 8:251-256, 1990. W . H . Ip. The fast atomic oxygen corona extent of Mars. Geop. Res. Lett., 17:2289-2292, 1990. V. I . Shematovich, D.V. Bisikalo, and J.C. Gerard. A kinetic model of the formation of the hot oxygen geocorona 1. quiet geomagnetic conditions. J. Geophys. Res., 99:23,217-23,228, 1994. M . B . McElroy, M . J . Prather, and J . M . Rodriguez. Escape of hydrogen from Venus. Science, 215:1614-1615, 1982. T .E . Donahue, J . H . Hoffman, R.R. Hodges, and A . J . Watson. Venus was wet: A measurement of the ratio of D to H. Science, 216:630-633, 1982. J . F . Kasting and J . B . Pollack. Loss of water from Venus. I. hydrodynamic escape of hydrogen. Icarus, 53:479-508, 1983. S.C. Liu and T . M . Donahue. The aeronomy of hydrogen in the atmosphere of the Earth. J. Atmos. Sci., 31:1118-1136, 1974. S.C. Liu and T . M . Donahue. Mesospheric hydrogen related to exospheric escape mechanisms. J. Atmos. Sci., 31:1466-1470, 1974. S.C. Liu and T . M . Donahue. Realistic model of hydrogen constituents in the lower atmosphere and the escape flux from the upper atmosphere. J. Atmos. Sci., 31:2238-2242, 1974. D . M . Hunten and D.F. Strobel. Production and escape of terrestrial hydrogen. J. Atmos. Sci., 31:305-317, 1974. J . L . Bertaux. Observed variations of the exospheric hydrogen density with the exospheric temperature. J. Geophys. Res., 80:639-642, 1975. J . W . Chamberlain. Charge exchange in a planetary corona; its effect on the distri-bution and escape of hydrogen. J. Geophys. Res., 82:1-9, 1977. R.R. Hodges, R.P. Rohrbaugh, and B . A . Tinsley. The effect of charge exchange source on the velocity and 'temperature' distributions and their anisotropics in the Earth's exosphere. J. Geophys. Res., 86:6917-6925, 1981. Bibliography 247 [45] B . Shizgal and M . J . Lindenfeld. A simple kinetic theory calculation of terrestrial atomic hydrogen escape induced by charge exchange collisions. J. Geophys. Res., 87:853-860, 1982. [46] D.R. Bates and M . R . C . McDowell. Escape of helium. J. Atmos. & Terr. Phys., 16:393-395, 1959. [47] 0 . Svendsen, M . H . Rees, and K . St amnes. Helium escape from the Earth's at-mosphere: the charge exchange mechanism revisited. Plan. Sp. Sci., 40:1639-1662, 1992. [48] T. Torgerson. Terrestrial helium degassing fluxes and the atmospheric helium bud-get: Implications with respect to the degassing processes of continental crust,. Chem. Geo/., 79:1-14, 1989. [49] 0 . Svendsen, M . H . Rees, and K . St amnes. Helium ion outflow from the Earth's atmosphere. EOS Trans. AGU, 75:490, 1994. [50] W . B . Maier. Reactions of H e + with N 2 and 0 2 in the upper atmosphere. Plan. Sp. Sci., 16:477-493, 1968. [51] T .E . Donahue. Evolution of water reservoirs on Mars from D / H ratios in the atmosphere and crust. Nature, 374:432-434, 1995. [52] T .E . Donahue. Water on Mars and Venus. Proc. of the Conference on Deep Earth and Planetary Volatiles, 1995. [53] S.W. Squyres and J.F. Kasting. Early Mars: How warm and how wet? Science, 265:744-749, 1994. [54] J .B. Pollack, J.F. Kasting, S.M. Richardson, and K . Poliakoff. The case for a wet, warm climate on early Mars. Icarus, 71:203-224, 1987. [55] T. Owen, J.P. Maillard, C. De Bergh, and B . J . Lutz. Deuterium on Mars: The abundance of DHO and the value of D / H . Science, 240:1767-1770, 1988. [56] R.O. Pepin. On the origin and early evolution of terrestrial planet atmospheres and meteoritic volatiles. Icarus, 92:2-79, 1991. [57] H . Y . McSween. What we have learned about Mars from SNC meteorites. Mete-oritics, 29:757-779, 1994. [58] D . M . Hunten. Atmospheric evolution of the terrestrial planets. Science, 259:915-920, 1993. [59] D . M . Hunten, R.O. Pepin, and J .C .G. Walker. Mass fractionation in hydrodynamic escape. Icarus, 69:532-549, 1987. Bibliography 248 [60] J .L. Fox. The production and escape of nitrogen atoms on Mars. J. Geophys. Res., 98:3297-3310, 1993. [61] M . B . McElroy and Y . L . Yung. Oxygen isotopes in the Martian atmosphere: Impli-cations for the evolution of volatiles. Plan. Sp. Sci., 24:1107-1113, 1976. [62] K . Zahnle, J .F. Kasting, and J.B. Pollack. Mass fractionation of noble gases in diffusion-limited hydrodynamic hydrogen escape. Icarus, 84:502-527, 1990. [63] R.R. Hodges. Isotopic fractionation of hydrogen in planetary exospheres due to ionosphere-exosphere coupling: Implications for Venus. J. Geophys. Res., 98:10833-10838, 1993. [64] R.O. Pepin. Evolution of the Martian atmosphere. Icarus, 111:289-304, 1994. [65] Jakosky, R.O. Pepin B . M . , R .E . Johnson, and J.L. Fox. Mars atmospheric loss and isotopic fractionation by solar-wind induced sputtering and photochemical escape. Icarus, 111:271-288, 1994. [66] X . He. Hydrogen Escape from the Earth's Atmosphere. PhD thesis, Boston Univer-sity, 1995. [67] F .F . Chen. Introduction to Plasma Physics and Controlled Fusion: Volume 1. Plenum Press, second edition, 1984. [68] D.A. McQuarrie. Statistical Mechanics. Harper Collins Publishers, 1976. [69] A . F . Nagy, T .E . Cravens, J . H . Yee, and A.I.F Stewart. Hot oxygen atoms in the upper atomosphere of Venus. Geop. Res. Lett., 8:629-632, 1981. [70] D.V. Bisikalo, V. I . Shematovich, and J.C. Gerard. A kinetic model of the formation of the hot oxygen geocorona 2. influence of 0 + ion precipitation. J. Geophys. Res., 100:3715-3720, 1995. [71] B . Shizgal. Hot hydrogen and deuterium in the exosphere of Venus. Adv. Sp. Res., 7:73-77, 1987. [72] A.S. Clarke and B.D. Shizgal. Relaxation dynamics of hot protons in a thermal bath of atomic hydrogen. Phys. Rev. E., 49:347-358, 1994. [73] R.R. Hodges and E .L . Breig. Ionosphere-exosphere coupling through charge ex-change and momentum transfer in hydrogen-proton collisions. J. Geophys. Res., 96:7697-7708, 1991. [74] R.R. Hodges and E .L . Breig. Charge transfer and momentum exchange in exo-spheric D-H+ and H-D+ collisions. J. Geophys. Res., 98:1581-1588, 1993. Bibliography 249 [75] D.L. Cooper, Y . H . Yee, and A. Dalgarno. Energy transfer in oxygen-hydrogen collisions. Plan. Sp. Sci., 32:825-830, 1984. [76] M . A . Gurwell and Y . H . Yung. Fractionation of hydrogen and deuterium on Venus due to collisional ejection. Plan. Sp. Sci., 41:91-104, 1993. [77] L . J . Maher and B . A . Tinsley. Atomic hydrogen escape rate due to charge exchange with hot plasmaspheric ions. J. Geophys. Res., 82:689-695, 1977. [78] Y . L . Yung, J.S. Wen, J.I. Moses, B . M . Landry, M . Allen, and K . J . Hsu. Hydrogen and deuterium loss from the terrestrial atmosphere: A quantitative assessment of nonthermal escape fluxes. J. Geophys. Res., 94:14,971-14,989, 1989. [79] M . B . McElroy, M . J . Prather, and J . M . Rodriguez. Loss of oxygen from Venus. Geop. Res. Lett., 9:649-651, 1982. [80] J .L. Fox. On the escape of oxygen and hydrogen from Mars. Geop. Res. Lett., 20:1847-1850, 1993. [81] M . B . McElroy. Mars: An evolving atmosphere. Science, 175:443-445, 1972. [82] M . G . Shepherd, J .C. McConnell, W . K . Tobiska, G.R. Gladstone, S. Chakrabati, and G. Schmidtke. Inference of atomic oxygen concentration from remote sensing of optical aurora. J. Geophys. Res., 100:17,415-17,428, 1995. [83] D . M . Cotton, G.R. Gladstone, and S. Chakrabarti. Sounding rocket observations of a hot atomic oxygen geocorona. J. Geophys. Res., 98:21,651-21,657, 1993. [84] R.P. Rohrbaugh and J.S. Nisbet. Effects of energetic oxygen atoms on neutral density models. J. Geophys. Res., 78:6768-6772, 1973. [85] B . Shizgaland M . J . Lindenfeld. Energy distribution function of translationally hot 0( 3p) atoms in the atmosphere of Earth. Plan. Sp. Sci., 27:1321-1332, 1979. [86] G .A. Schmitt, V . J . Abreu, and P.B. Hays. Nonthermal 0(x