Molecular-Kinetic Simulations of Escape from the Ex-planet and Exoplanets: Criterion for Transonic Flow
MMolecular-Kinetic Simulations of Escape from the Ex-planet and Exoplanets: Criterion for Transonic Flow
Robert E. Johnson , Alexey N. Volkov and Justin T. Erwin Engineering Physics, University of Virginia, Charlottesville, VA 22904-4745, USA Physics Department, New York University, NY 10003-6621, USA Department of Mechanical Engineering, University of Alabama, Tuscaloosa, AL 35487 Department of Planetary Science, University of Arizona, Tucson, AZ 85721, USA
Abstract:
The equations of gas dynamics are extensively used to describe atmospheric loss from solar system bodies and exoplanets even though the boundary conditions at infinity are not uniquely defined. Using molecular-kinetic simulations that correctly treat the transition from the continuum to the rarefied region, we confirm that the energy-limited escape approximation is valid when adiabatic expansion is the dominant cooling process. However, this does not imply that the outflow goes sonic. In fact in the sonic regime, the energy limited approximation can significantly under estimate the escape rate. Rather large escape rates and concomitant adiabatic cooling can produce atmospheres with subsonic flow that are highly extended. Since this affects the heating rate of the upper atmosphere and the interaction with external fields and plasmas, we give a criterion for estimating when the outflow goes transonic in the continuum region. This is applied to early terrestrial atmospheres, exoplanet atmospheres, and the atmosphere of the ex-planet, Pluto, all of which have large escape rates. The paper and its erratum, combined here, are published: ApJL 768, L4 (2013); ApJ, 779, L30 (2013).
Introduction
Rapid atmospheric escape is often described as a gas that goes sonic, sometimes called blow-off (Hunten 1982), a process that can account for certain isotope ratios on terrestrial planets. Transonic models have also been used to describe rapid escape from exoplanets (e.g., Murray-Clay et al. 2009) and from Pluto (e.g., Strobel 2008). However, we recently showed that this model for Pluto gave an incorrect upper atmospheric structure (Tucker et al. 2012; Erwin et al. 2013). In simulating rapid escape using continuum gas dynamics, the Jeans expressions at the exobase have been applied for the uncertain boundary conditions at infinity (e.g., Tian et al. 2008). More often, a sonic point is assumed to occur at some altitude, above which the density and temperature dependence can be simply characterized (Parker 1964a,b). The so-called energy-limited escape rate, extensively applied to exoplanet atmospheres (e.g., Lammer et al. 2009), is often assumed to imply that sonic boundary conditions are applicable (e.g., Erkaev et al. 2012). Here we use molecular-kinetic simulations to show that is not the case. We briefly review the continuum and molecular-kinetic models, and then present results of our simulations. These test the applicability of the energy-limited escape rate and our proposed criterion for determining whether sonic or kinetic upper boundary conditions are applicable. The results are applied to escape from Pluto, early terrestrial planet, and exoplanet atmospheres. odels
We describe escape from a one-dimensional (1D), steady state, single component atmosphere as illustrative, leaving out thermal transport by horizontal flow. For radial distance (cid:1870) , flow speed (cid:1873) , number density (cid:1866) , temperature (cid:1846) , pressure (cid:1868) (cid:3404) (cid:1866)(cid:1863)(cid:1846) , and escape rate
Φ (cid:3404) 4(cid:2024)(cid:1870) (cid:2870) (cid:1866)(cid:1873) (cid:3404) const, the momentum and energy equations are often used ignoring viscosity: (cid:1866) (cid:3031)(cid:3031)(cid:3045) (cid:4666)(cid:1865)(cid:1873) (cid:2870) /2 – (cid:1847)(cid:4667) (cid:3404) (cid:3398) (cid:3031)(cid:3043)(cid:3031)(cid:3045) (1) (cid:3031)(cid:3031)(cid:3045) (cid:4674)Φ(cid:4666)(cid:1865)(cid:1873) (cid:2870) /2 (cid:3397) (cid:1829) (cid:3043) (cid:1863)(cid:1846) (cid:3398) (cid:1847)(cid:4667) (cid:3398) 4(cid:2024)(cid:1870) (cid:2870) κ (cid:3031)(cid:3021)(cid:3031)(cid:3045) (cid:4675) (cid:3404) 4(cid:2024)(cid:1870) (cid:2870) (cid:1866) (cid:1869) (cid:3028) (cid:4666)(cid:1870)(cid:4667) (2) Here (cid:1863) is the Boltzmann constant, κ (cid:3404) κ (cid:4666)(cid:1846)(cid:4667) the thermal conductivity, (cid:1829) (cid:3043) the specific heat at constant pressure, (cid:1865) the molecular mass, (cid:1847) (cid:3404) (cid:1847)(cid:4666)(cid:1870)(cid:4667) (cid:3404) (cid:1833)(cid:1839)(cid:1865)/(cid:1870) the gravitational energy ( (cid:1833) the gravitational constant, (cid:1839) the planet’s mass), and (cid:1869) (cid:3028) (cid:4666)(cid:1870)(cid:4667) is the net heating rate per molecule produced by incident photons or plasma particles, in which we include radiative cooling. Knowing the density, (cid:1866) (cid:2868) , and temperature, (cid:1846) (cid:2868) , at a lower boundary, (cid:1870) (cid:3404) (cid:1870) (cid:2868) , a unique solution requires two other conditions, typically at the upper boundary. The gravitational energy is characterized by the Jeans parameter, λ ( (cid:1870) ) (cid:3404) (cid:1847)/(cid:1863)(cid:1846) and the rarefaction by the Knudsen number, (cid:1837)(cid:1866)(cid:4666)(cid:1870)(cid:4667) (cid:3404) (cid:1864) (cid:3030) /(cid:1834) , the ratio of the mean free path of gas molecules, (cid:1864) (cid:3030) , to the local scale height, (cid:1834) (cid:3404) (cid:3398)(cid:1866)/(cid:4666)(cid:1856)(cid:1866)/(cid:1856)(cid:1870)(cid:4667) . For an escaping gas at large distances from the source, where free molecular flow occurs, (cid:1834) → (cid:1870)/2 ; in the hydrostatic regime (cid:1834) → (cid:1870)/ λ (r). The Jeans expressions for the number, Φ (cid:3011) , and thermal, 〈(cid:1831)Φ〉 (cid:3011) , escape rates have been used as upper boundary conditions for Equations (1-2): Φ (cid:3011) (cid:3404) 4(cid:2024)(cid:1870) (cid:3051)(cid:2870) (cid:1866) (cid:3051) (cid:3495) (cid:3038)(cid:3021) (cid:3299) (cid:2870)(cid:3095)(cid:3040) (cid:4666)1 (cid:3397) (cid:2019) (cid:3051) (cid:4667)exp (cid:4666)(cid:3398)(cid:2019) (cid:3051) (cid:4667) (3a) 〈(cid:1831)Φ〉 (cid:3011) (cid:3404) (cid:1863)(cid:1846) (cid:3051) Φ (cid:3011) (cid:4672) (cid:2869)(cid:2869)(cid:2878)(cid:3090) (cid:3299) (cid:3397) (cid:1829) (cid:3043) (cid:3398) (cid:2871)(cid:2870) (cid:4673) (3b) The subscript ‘ x ’ indicates quantities evaluated at the nominal exobase, (cid:1870) (cid:3404) (cid:1870) (cid:3051) , where (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:3051) (cid:4667) (cid:3406) 1 , often assumed to be the upper boundary of the continuum region. When the upper atmosphere heating rate is large, the equations are more often solved through a sonic point, (cid:1870) (cid:3404) (cid:1870) ∗ , where (cid:1873) ∗ (cid:3404) (cid:1855) ( (cid:1855) (cid:3404) (cid:3493) γ (cid:1863)(cid:1846)/(cid:1865) the sound speed, γ (cid:3404) (cid:1829) (cid:3043) /(cid:1829) (cid:3023) , (cid:1829) (cid:3023) the specific heat at constant volume). For (cid:1870) ≫ (cid:1870) ∗ , (cid:1866) and (cid:1846) decay as power laws (Parker 1964a,b). Unfortunately those continuum solutions for which Jeans escape is applicable and those for which a sonic point is reached in the continuum region do not simply track from one to the other as the heating rate increases. Kinetic models can simulate both continuum (cid:4666)(cid:1864) (cid:3030) ≪ (cid:1834)(cid:4667) and non-continuum (transitional, (cid:1864) (cid:3030) ~(cid:1834) , and free molecular, (cid:1864) (cid:3030) ≫ (cid:1834)(cid:4667) gas flows and can, therefore, describe the change from Jeans-like to transonic escape. Since such simulations track particles in the potential well of the body (or bodies) of interest, escape is a natural outcome. We numerically implement a kinetic description of rarefied gas flow in an upper atmosphere based on the Boltzmann kinetic equation using the direct simulation Monte Carlo DSMC) method (Bird 1994). In this method, the gas flow is represented by a large set of representative atoms or molecules that are tracked subject to binary collisions and gravity (Volkov 2011a,b). Heating of the atmosphere is implemented by scaling the thermal velocities of the representative molecules according to the local energy deposition rate. The lower boundary of the simulation region, (cid:1870) (cid:3404) (cid:1870) (cid:2868) , is below the depth at which the UV/EUV or plasma energy deposition occurs and (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:2868) (cid:4667) ≪ 1 . Because the density drops rapidly with increasing (cid:1870) , but escape occurs at large (cid:1870) where the density is low, DSMC simulations starting at small (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:2868) ) can require an enormous number of particles to accurately describe escape. Therefore, we also use a hybrid continuum/kinetic model (Tucker et al. 2012; Erwin et al. 2013) in which the Equations (1-2) are solved at (cid:1837)(cid:1866)(cid:4666)(cid:1870)(cid:4667) (cid:3407) ~0.1 (cid:3398) 0.01 , where the gas is collisionally dominated, and the velocity and internal energy distributions are reasonably well represented by Maxwellians, and then iteratively couple it to a DSMC simulation in the rarefied region. Heating
Parker (1964a,b) used Equations (1-2) to describe escape when the dominant heat source is internal, as it is for expansion of the solar corona: i.e., (cid:1869) (cid:3028) (cid:4666)(cid:1870)(cid:4667) (cid:3404) 0 for (cid:1870) (cid:3408) (cid:1870) (cid:2868) . This model was subsequently applied to planetary atmospheres primarily heated at (cid:1870) (cid:3407)(cid:1870) (cid:2868) . For Jeans parameters at (cid:1870) (cid:3404) (cid:1870) (cid:2868) as large as λ ( (cid:1870) (cid:2868) )= λ (cid:2868) ~ 40 , such models were assumed to produce a transonic expansion, often referred to as ‘slow hydrodynamic escape’ (e.g., Strobel 2008). Although, rapid escape can occur for relatively large λ (cid:2868) and (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:2868) (cid:4667) ≪ 1 , for λ (cid:2868) (cid:3408) ~(cid:4666)(cid:1829) (cid:3043) (cid:3397) (cid:2011)/2(cid:4667) the gas does not go sonic in the collision-dominated region and the escape rate is a factor of a few larger than the Jeans rate (Volkov et al. 2011a,b). Ignoring the thermal conductivity in Equation (2), this corresponds to the enthalpy of fluid particles becoming sufficient to produce a transonic, isentropic expansion starting at r . For smaller λ (cid:2868) a non-equilibrium region, called a Knudsen layer, forms above (cid:1870) (cid:3404) (cid:1870) (cid:2868) and there is a steep transition with decreasing λ (cid:2868) to supersonic escape at λ (cid:2868) ~2.1 and ~2.8-3.5 for monatomic and diatomic gases respectively (Volkov & Johnson 2013). In an upper atmosphere heated by short wavelength radiation, or by incident plasma particles, escape is driven by the energy absorbed. Energy absorbed below (cid:1837)(cid:1866)(cid:4666)(cid:1870)(cid:4667) ~0.1 is typically converted to heat using an efficiency, ε , that depends on the radiation type and atmospheric composition. The heating rate is either directly calculated or a value of ε is estimated: often a constant (cid:4666)~0.15 (cid:3398) 0.4 ) up to the exobase where it goes to zero. Typically the gas-dynamic equations for exoplanet or early terrestrial atmospheres are then solved with Jeans-like or sonic upper boundary conditions. As discussed below, we used a hybrid continuum/DSMC model to describe escape from Pluto heated by the solar UV/EUV and a DSMC model to describe escape from an atmosphere in which, the heating is assumed to occur in a narrow layer. Energy Limited Escape
Because thermal conduction in the upper atmosphere is inefficient, adiabatic cooling by escape or horizontal transport often dominates (e.g., Erwin et al. 2013). For a globally averaged heating rate and adiabatic cooling, integration of Equation (2) gives a rough upper bound to the escape rate, Φ (cid:3006)(cid:3013) (e.g., Lammer et al. 2009): Φ (cid:3006)(cid:3013) (cid:3406) (cid:1843) (cid:3041)(cid:3032)(cid:3047) /(cid:4666)(cid:1847) (cid:3398) (cid:1829) (cid:3043) (cid:1863)(cid:1846) (cid:3398) (cid:1865)(cid:1873) (cid:2870) /2(cid:4667)| (cid:3045)(cid:3048)(cid:3045)(cid:2868) . (4) where the denominator is the difference in the average energy of a molecule between the lower, (cid:1870) (cid:2868) , and upper, (cid:1870) (cid:3048) , boundaries of the simulation. In Equation (4) (cid:1870) (cid:2868) is assumed to be below the heated region, above which adiabatic cooling dominates; (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3404)4(cid:2024) (cid:1516) (cid:1870) (cid:2870) (cid:1869) (cid:3028) (cid:4666)(cid:1870)(cid:4667)(cid:1856)(cid:1870) (cid:2998)(cid:3045) (cid:3116) is the integrated heating + radiative cooling rate. This expression is often referred to as the energy-limited rate, although Watson et al. (1981) discussed a related quantity. For (cid:1870) (cid:2868) deep in the gravitational well and λ (cid:2868) ≫ (cid:1829) (cid:3043) , Equation (4) is often approximated as Φ (cid:3006)(cid:3013) (cid:3406) (cid:1843) (cid:3041)(cid:3032)(cid:3047) /(cid:1847)(cid:4666)(cid:1870) (cid:2868) (cid:4667) . (5) Assuming a small fraction of (cid:1843) (cid:3041)(cid:3032)(cid:3047) is deposited for (cid:1837)(cid:1866)(cid:4666)(cid:1870)(cid:4667) (cid:3408) ~0.1 , so non-thermal escape processes (Johnson et al. 2008) can be ignored, we showed that UV/EUV heating of Pluto’s atmosphere resulted in an escape rate very close to that in Equation (5), but the gas did not go sonic below the exobase. Rather, a large expansion of the upper atmosphere occurred (Tucker et al. 2012; Erwin et al. 2013). Similarly, Tian et al. (2008) found that above a heating threshold, the atmosphere rapidly expanded and the escape rate increased with increasing EUV heating, consistent with energy-limited escape, even though the gas remained subsonic. Therefore, the energy limited escape rate is not contingent on the flow going sonic below the exobase. If the heating rate is increased to the point where the atmosphere does go sonic in the continuum region, energy limited escape can still apply if one accounts for the large (cid:1873) , enhanced radiative cooling, recombination in an ionized atmosphere, etc. (e.g., Murray-Clay et al. 2009). Criteria for Transonic Solutions
Since the isentropic approximation and energy-limited escape are applicable to both sub-sonic and transonic rapid outflows, we use the Mach number, (cid:1839)(cid:1853) (cid:3404)(cid:1873)/(cid:3493)γ(cid:1863)(cid:1846)/(cid:1865) , and rewrite Equation (5): (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3406) 4(cid:2024)(cid:1870) (cid:2870) (cid:1866) (cid:1839)(cid:1853)(cid:3495) (cid:3082)(cid:3090) (cid:3022)(cid:4666)(cid:3045)(cid:4667)(cid:3040) (cid:1847)(cid:4666)(cid:1870) (cid:2868) (cid:4667) (6) Because the boundary conditions for sub-sonic and transonic solutions differ, it is important to be able to estimate the minimum value of (cid:1843) (cid:3041)(cid:3032)(cid:3047) required to apply sonic boundary conditions, which we will call (cid:1843) (cid:3030) . Assuming the sonic point, (cid:1870) (cid:3404) (cid:1870) ∗ , occurs in the continuum region, the flow can be effectively approximated by the isentropic model: (cid:2870) (cid:4666)(cid:1870) ∗ (cid:4667) (cid:3404) (cid:1847)(cid:4666)(cid:1870) ∗ (cid:4667) (cid:3397) (cid:4666) γ (cid:3398) 1(cid:4667) (cid:3045) ∗ (cid:3044)(cid:4666)(cid:3045) ∗ (cid:4667)(cid:3048)(cid:4666)(cid:3045) ∗ (cid:4667) . (7) (e.g., Murray-Clay et al. 2009). For (cid:1869)(cid:4666)(cid:1870) ∗ (cid:4667) (cid:3406) 0 Equation (7) reduces to λ ( (cid:1870) ∗ )= λ ∗ (cid:3406) 2(cid:2011) and transonic escape occurs when (cid:1839)(cid:1853) (cid:3408) 1 in Equation (6) giving: (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3408) (cid:1843) (cid:3030) (cid:3406) 4(cid:2024)(cid:1870) ∗(cid:2870) (cid:1866) ∗ (cid:3495) (cid:3022)(cid:4666)(cid:3045) ∗ (cid:4667)(cid:2870)(cid:3040) (cid:1847)(cid:4666)(cid:1870) (cid:2868) (cid:4667), (8) ith (cid:1866) ∗ (cid:3404) (cid:1866)(cid:4666)(cid:1870) ∗ (cid:4667). Below we estimate (cid:1843) (cid:3030) and test it against molecular kinetic simulations. We first consider a narrow heating layer located at (cid:1870) (cid:3404) (cid:1870) (cid:3028) (cid:3408) (cid:1870) (cid:2868) with (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:3028) (cid:4667) ≪ 1 , as in Watson et al. (1981) and McNutt (1989), and give an approximate analytic solution to Equations (1-2) in the Appendix. For this case, we performed DSMC simulations for a monatomic gas of hard spheres with λ (cid:2868) (cid:3404) 10 , (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:2868) (cid:4667) (cid:3404) 10 (cid:2879)(cid:2871) , (cid:1870) (cid:3028) /(cid:1870) (cid:2868) (cid:3404) 1.1 . Although the results can be roughly scaled, the simulation results shown in Figures 1 and 2 are for an N atmosphere with hard sphere collision diameter (cid:2879)(cid:2869)(cid:2868) (cid:1865) , planet mass (cid:1839) (cid:3404)4.45(cid:1876)10 (cid:2870)(cid:2869) (cid:1863)(cid:1859) , (cid:1870) (cid:2868) (cid:3404) 10 (cid:2874) (cid:1865) , (cid:1866) (cid:2868) (cid:3404) 0.993(cid:1876)10 (cid:2869)(cid:2874) /(cid:1865) (cid:2871) , and (cid:1846) (cid:2868) (cid:3404) 100(cid:1837) . For (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3404) 0, we showed earlier that these conditions correspond to enhanced Jeans-like escape with Φ (cid:3406) 1.6Φ (cid:3011) (Volkov et al. 2011a,b) which for our simulation parameters is ~0.7(cid:1876)10 (cid:2879)(cid:2872) Φ (cid:2868) , where Φ (cid:2868) is the upward flow across the lower boundary of the simulation, (cid:1870) (cid:2868) : Φ (cid:2868) (cid:3404) 4(cid:2024)(cid:1870) (cid:2868)(cid:2870) (cid:1866) (cid:2868) (cid:4666)(cid:1863)(cid:1846) (cid:2868) /2(cid:2024)(cid:1865)(cid:4667) (cid:2868).(cid:2873) . Increasing (cid:1843) (cid:3041)(cid:3032)(cid:3047) to find when the outflow goes sonic in the continuum regime, it is seen in Figure 1 that for (cid:1843) (cid:3041)(cid:3032)(cid:3047) ≫ 0 a non-equilibrium layer forms near (cid:1870) (cid:3028) in which the parallel and perpendicular components of temperature differ. For small (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:3028) (cid:4667) the flow properties in this layer are analogous to those in the Knudsen layer discussed above for heating below (cid:1870) (cid:2868) . Near (cid:1870) (cid:3028) the density and temperature change dramatically but the pressure from (cid:1870) (cid:2868) to ~(cid:1870) (cid:3028) can be estimated from the hydrostatic approximation. When the outflow goes sonic, a pressure drop occurs from ~ (cid:1870) (cid:3028) to (cid:1870) ∗ : (cid:1842) (cid:3030) (cid:3404) (cid:1868) ∗ /(cid:1868) (cid:3028) (Volkov and Johnson 2013). Since the transition layer is narrow, (cid:1870) ∗ (cid:3406) (cid:1870) (cid:3028) , we rewrite Equation (8) using values at (cid:1870) (cid:2868) : (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3408) (cid:1843) (cid:3030) (cid:3406) 〈(cid:1831)Φ〉 (cid:2868) (cid:4670)(cid:2011)(cid:3493)(cid:2024)(cid:2019) (cid:2868) (cid:4672) (cid:3045) (cid:3276) (cid:3045) (cid:3116) (cid:4673) (cid:3121)(cid:3118) (cid:1842) (cid:3030) (cid:3043) (cid:3276) (cid:3043) (cid:3116) (cid:4671) (9) Here (cid:3407) (cid:1831)Φ (cid:3408) (cid:2868) (cid:3404) (cid:4666)2(cid:1863)(cid:1846) (cid:2868) (cid:4667)Φ (cid:2868) (cid:3404) 4(cid:2024)(cid:1870) (cid:2868)(cid:2870) (cid:1866) (cid:2868) (cid:1863)(cid:1846) (cid:2868) (cid:4666)2(cid:1863)(cid:1846) (cid:2868) /(cid:2024)(cid:1865)(cid:4667) (cid:2868).(cid:2873) is the upward Maxwellian energy flux of molecules leaving the source at (cid:1870) (cid:3404) (cid:1870) (cid:2868) . When the effects of thermal conduction are small relative to adiabatic cooling, then Equation (9) can be derived from the analytic expression in the Appendix. The barometric equation gives (cid:1868) (cid:3028) /(cid:1868) (cid:2868) (cid:3406) exp (cid:4670) λ (cid:2868) (cid:4666)(cid:1870) (cid:2868) /(cid:1870) (cid:3028) (cid:3398) 1(cid:4667)(cid:4671) with (cid:1842) (cid:3030) depending on the number of degrees of freedom: ~0.4 for a monatomic gas. From Figure 1, the transition to supersonic flow occurs for (cid:3041)(cid:3032)(cid:3047) /(cid:1843) (cid:3030) (cid:3407) 0.67 with Q (cid:3030) ~4.5(cid:1876)10 (cid:2869)(cid:2868) (cid:1849) in Equation 9 for our simulation parameters. Due to the above approximations, Equation (9) overestimates (cid:1843) (cid:3030) by about a factor of two. It is also seen in Figure 1(c) that (cid:1866)(cid:4666)(cid:1870)(cid:4667) for the subsonic solution increases slowly for (cid:1870) ≫ (cid:1870) (cid:3028) , resulting in a significantly expanded atmosphere, but (cid:1866)(cid:4666)(cid:1870)(cid:4667) for the transonic solution in Figure 1(e) roughly decreases as a power law consistent with transonic escape. Figure 1.
DSMC simulations: a monatomic gas of hard spheres at (cid:2019) (cid:2868) (cid:3404) 10 , (cid:1837)(cid:1866)(cid:4666)(cid:1870) (cid:2868) (cid:4667) (cid:3404) 10 (cid:2879)(cid:2871) heated at (cid:1870) (cid:3028) /(cid:1870) (cid:2868) (cid:3404) 1.1 ,with (cid:1870) (cid:2868) (cid:3404) (cid:1844) (cid:2868) . Scaled density, (cid:1866)/(cid:1866) (cid:2868) (black), parallel, (cid:1846) || /(cid:1846) (cid:2868) (green), and perpendicular, (cid:1846) (cid:2884) /(cid:1846) (cid:2868) (blue), temperatures, and local Mach (red) and Knudsen (magenta) numbers for (cid:1843) (cid:3041)(cid:3032)(cid:3047) /(cid:1843) (cid:3030) (cid:3404) 0 (a,b; Jean-like escape), = 0.46 (c,d; subsonic; Φ/Φ (cid:3006)(cid:3013) (cid:3404) 0.87 ), and = 0.67 (e,f; transonic at (cid:1870) ∗ /(cid:1870) (cid:2868) (cid:3404) 1.13 ; Φ/Φ (cid:3006)(cid:3013) (cid:3404)0.65 ). Lines: Heated and sonic surfaces. s important, the escape rate increases dramatically above (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3404) 0 and becomes close to Φ (cid:3006)(cid:3013) in Equation (5) in the subsonic regime as seen in Figure 2. Although Φ does not change significantly in the transition to transonic escape, a steep increase is seen in the thermal + flow energy removed, (cid:3407) (cid:1831) Φ (cid:3408) . It is also seen that additional heating primarily increases the average kinetic energy of escaping molecules so that Equation (5) becomes a poor approximation if (cid:1843) (cid:3041)(cid:3032)(cid:3047) ≫ (cid:1843) (cid:3030) . For very large Q (cid:3041)(cid:3032)(cid:3047) in Figure 2 the escape rate is limited to a fraction of Φ (cid:2868) and (cid:3407) (cid:1831)Φ (cid:3408) approaches Q (cid:3041)(cid:3032)(cid:3047) . For small Q (cid:3041)(cid:3032)(cid:3047) downward thermal conduction affects the escape rate (Erwin et al. 2013), which is eventually driven by the temperature at (cid:1870) (cid:2868) . Figure 2.
Number,
Φ (cid:4670)green circles , red squares] and energy, 〈(cid:1831)Φ〉 [blue triangles] escape rates vs. (cid:1843) (cid:3041)(cid:3032)(cid:3047) /(cid:1843) (cid:3004) calculated as for Figure 1: Scaled to Φ (cid:3006)(cid:3013) from Equation (5) and also to Φ (cid:2868) , the upward flow across the lower boundary of the simulation, (cid:1870) (cid:2868) , Φ (cid:2868) (cid:3404) 4(cid:2024)(cid:1870) (cid:2868)(cid:2870) (cid:1866) (cid:2868) (cid:4666)(cid:1863)(cid:1846) (cid:2868) /2(cid:2024)(cid:1865)(cid:4667) (cid:2868).(cid:2873) (in Volkov et al. (2011a,b) the symbol Φ (cid:2868),(cid:2868) is used). (cid:3407) (cid:1831) Φ (cid:3408) is scaled to the energy flux of molecules across (cid:1870) (cid:2868) , (cid:3407) (cid:1831)Φ (cid:3408) (cid:2868) (cid:3404)(cid:4666)2(cid:1863)(cid:1846) (cid:2868) (cid:4667)Φ (cid:2868) (cid:3404) 4(cid:2024)(cid:1870) (cid:2868)(cid:2870) (cid:1866) (cid:2868) (cid:1863)(cid:1846) (cid:2868) (cid:4666)2(cid:1863)(cid:1846) (cid:2868) /(cid:2024)(cid:1865)(cid:4667) (cid:2868).(cid:2873) . Rectangle indicates transition from subsonic to supersonic flow below the exobase. For these simulations Φ (cid:2868) (cid:3404) 8.6(cid:1876)10 (cid:2871)(cid:2868) (cid:1871) (cid:2879)(cid:2869) ; (cid:3407) (cid:1831)Φ (cid:3408) (cid:2868) (cid:3404) 2.4(cid:1876)10 (cid:2869)(cid:2868) (cid:1849); (cid:1843) (cid:3030) ~4.5(cid:1876)10 (cid:2869)(cid:2868) (cid:1849) from Equation 9. We now consider more realistic heating profiles. In order to apply sonic boundary conditions in Equations (1-2), (cid:1837)(cid:1866)(cid:4666)(cid:1870) ∗ (cid:4667) must be in the continuum region below some maximum, (cid:1837)(cid:1866) (cid:3040) : i.e., (cid:1837)(cid:1866)(cid:4666)(cid:1870) ∗ (cid:4667) (cid:3404) (cid:1864) (cid:3030)∗ /(cid:1834) ∗ (cid:3407) (cid:1837)(cid:1866) (cid:3040) . Using (cid:1834) ∗ ~ (cid:1870) ∗ /(cid:2019) ∗ and (cid:1864) (cid:3030)∗ (cid:3404) 1/(cid:4666)(cid:1855) (cid:3030) (cid:2026) (cid:3030) (cid:1866) ∗ (cid:4667) , where (cid:1855) (cid:3030) is determined by the energy dependence of the total collision cross section (cid:2026) (cid:3030) (e.g., (cid:1855) (cid:3030) (cid:3404) √2 , (cid:2026) (cid:3030) (cid:3404) (cid:2024)(cid:1856) (cid:2870) for the gas of hard sphere molecules of diameter d ), then (cid:3030) (cid:1870) ∗ (cid:2026) (cid:3030) (cid:1866) ∗ (cid:4667) (cid:3407) (cid:1837)(cid:1866) (cid:3040) . From Equation (8), we estimate (cid:1843) (cid:3030) requiring as that (cid:1870) ∗ occurs in the continuum domain: (cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3408) (cid:1843) (cid:3030) (cid:3406) 4(cid:2024)(cid:1870) ∗ (cid:3082)(cid:3030) (cid:3278) (cid:3097) (cid:3278) (cid:3012)(cid:3041) (cid:3288) (cid:3495) (cid:2870)(cid:3022)(cid:4666)(cid:3045) ∗ (cid:4667)(cid:3040) (cid:1847)(cid:4666)(cid:1870) (cid:2868) (cid:4667), (10) here (cid:1870) (cid:2868) (cid:3407) (cid:1870) ∗ (cid:3407) (cid:1870) (cid:2934) . When there is a sharp change in the gas properties, as for the heated layers discussed above, then (cid:1837)(cid:1866) (cid:3040) (cid:3407) ~ 0.1 (Volkov and Johnson 2013). However, if the heat is primarily absorbed over a board range of (cid:1870) below (cid:1870) (cid:3051) , then (cid:1837)(cid:1866) (cid:3040) ~ 1 is sufficient. It is seen in Equation 10 that (cid:1843) (cid:3030) does not explicitly depend on (cid:1846) (cid:2868) , consistent with simulations when Φ is large. Because (cid:1843) (cid:3030) depends on the sonic point only via (cid:4666)(cid:1870) ∗ (cid:4667) (cid:2869)/(cid:2870) , a rough lower bound can be obtained by replacing (cid:1870) ∗ with the mean energy absorption depth, (cid:1870) (cid:3028) estimated from the absorption cross section, (cid:2026) (cid:3028) . More accurately, at threshold the sonic point approaches (cid:1870) (cid:3051) so that (cid:1870) ∗ ~(cid:1870) (cid:3028) (cid:4670)1 (cid:3397) (cid:4666)(cid:2026) (cid:3028) /(cid:1855) (cid:3030) (cid:2026) (cid:3030) (cid:4667) (cid:2019) (cid:3028)(cid:3049)(cid:3032) (cid:4671) where (cid:2019) (cid:3028)(cid:3049)(cid:3032) ~(cid:4666)(cid:2019) (cid:3028) (cid:3397)2(cid:2011)(cid:4667)/2 slightly increasing (cid:1843) (cid:3030) . For a close-in exoplanet, tidal heating can be included in (cid:1847)(cid:4666)(cid:1870)(cid:4667) and ion escape can dominate so that (cid:2026) (cid:3030) becomes large due to ion-neutral or ion-ion collisions reducing (cid:1843) (cid:3030) . For solar minimum, medium, and maximum conditions we simulated Pluto’s upper atmosphere at the New Horizons encounter using our hybrid fluid/kinetic model ignoring the interaction with the solar wind and Charon, as well as non-thermal escape (Erwin et al. 2013). For all three cases ( (cid:1843) (cid:3041)(cid:3032)(cid:3047) = 0.38, 0.78, ~1.6x10 W) the atmosphere became highly extended, but the flow remained subsonic contrary to all earlier models (e.g., Strobel 2008) with the escape rate close to the energy-limited estimate in Equation (4). It is seen in Figure 3 that at solar medium (cid:1870) (cid:3051) is more than twice that obtained when sonic boundary conditions are applied and, although the escape rate is large, using the Jeans boundary conditions results a much better approximation to the upper atmosphere structure.
Figure 3.
Fluid simulations of Pluto’s N atmosphere at 32Au at solar-medium; upper bc: (solid) fluid/ DSMC coupled at (cid:1870) (cid:3047) , (cid:1837)(cid:1866) (cid:3404) 0.1; Φ (cid:3404) 2.6x10 (cid:2870)(cid:2875) ; (dashed) transonic assumption: (cid:1870) ∗ (cid:3408) (cid:1870) (cid:3051) , Φ (cid:3404) 2.5x10 (cid:2870)(cid:2875) /s (Strobel 2008); (dotted) Jeans bc from Equations (3a,b) at (cid:1870) (cid:3051) ~(cid:4666)7 (cid:3398) 8(cid:4667)(cid:1870) (cid:3043) ; Φ (cid:3404) 2.6x10 (cid:2870)(cid:2875) /(cid:1871); (cid:1870) (cid:3043) (cid:3404) 1153 km
Pluto’s radius; (cid:1870) (cid:2868) (cid:3404) 1.25(cid:1870) (cid:3043) roughly the visible extinction radius; parameters in VHS-LB model (Erwin et al. 2013): (cid:1865) (cid:3404) 28 amu, (cid:1866) (cid:2868) (cid:3404) 4x10 (cid:2869)(cid:2870) cm (cid:2879)(cid:2871) , (cid:1846) (cid:2868) (cid:3404) 88.2 K; (cid:1847)(cid:4666)(cid:1870) (cid:2868) (cid:4667) (cid:3404) 2.8x10 (cid:2879)(cid:2869)(cid:2871) ergs ; λ (cid:2868) (cid:3404) 23; γ (cid:3404)7/5 ; σ (cid:2913) ~9x10 (cid:2879)(cid:2869)(cid:2873) cm (cid:2870) ; Kn(cid:4666)(cid:1870) (cid:2868) (cid:4667) ~ 10 (cid:2879)(cid:2874) . For UV/EUV absorption at (cid:1870) (cid:3028) ~1.5 times Pluto’s radius, (cid:1870) (cid:3043) , using (cid:1837)(cid:1866) (cid:3040) ~1 to get a very rough lower bound and (cid:1870) ∗ ~(cid:1870) (cid:3028) ~(cid:1870) (cid:2868) , Equation (9) gives (cid:1843) (cid:3030) ≫ ~10 (cid:1876)10 (cid:2876) (cid:1849) . This is R ad i u s (r p ) r x r (cid:31) r x r x r t Fluid-DSMCFluid-JeansFluid-Sonic Density (m − ) ell above the largest heating rate (~1.6 x W ) in Erwin et al. (2013). Therefore, Pluto’s atmosphere at the New Horizon encounter will be highly extended with an escape rate very close to the energy-limited rate, but the flow in the continuum region will be sub-sonic. In order to further test Equation (10), we performed DSMC simulations with a distributed heating model in which, for simplicity, we used Beer’s law along the radial direction (Murray-Clay et al. 2009): (cid:1869) (cid:3028) (cid:4666)(cid:1870)(cid:4667) (cid:3404) ε (cid:2026) (cid:3028) (cid:1866)(cid:4666)(cid:1870)(cid:4667)(cid:1832) (cid:3022)(cid:3023)/(cid:3006)(cid:3022)(cid:3023) exp (cid:4666)(cid:3398)(cid:2026) (cid:3028) (cid:1840)(cid:4666)(cid:1870)(cid:4667)(cid:4667) , where (cid:2026) (cid:3028) , ε , and (cid:1832) (cid:3022)(cid:3023)/(cid:3006)(cid:3022)(cid:3023) are the absorption cross section, the heating efficiency and solar energy flux at the upper boundary of our domain, (cid:1870) (cid:3404) (cid:1870) (cid:3048) , and (cid:1840)(cid:4666)(cid:1870)(cid:4667) (cid:3404) (cid:1516) (cid:1866)(cid:4666)(cid:1870)̅(cid:4667)(cid:1856)(cid:1870)̅ (cid:3045) (cid:3296) (cid:3045) , with constant (cid:2013) out to (cid:1870) (cid:3048) . The transition to a transonic solution with (cid:1870) ∗ in the continuum region of the atmosphere was found to occur at (cid:3041)(cid:3032)(cid:3047) /(cid:1843) (cid:3030) (cid:3407) 1.6 , in agreement with our criterion (cid:1843) (cid:3030) calculated from Equation (10). The threshold for transonic flow can be related to the absorbed energy converted to heat, (cid:1843) (cid:3028) , by accounting for the non-adiabatic cooling processes , (cid:1843) (cid:3030)(cid:3042)(cid:3042)(cid:3039) : (cid:1843) (cid:3041)(cid:3032)(cid:3047) ~(cid:1843) (cid:3028) (cid:3398)(cid:1843) (cid:3030)(cid:3042)(cid:3042)(cid:3039) (Erwin et al. 2013, Figure 4b). If (cid:1843) (cid:3030)(cid:3042)(cid:3042)(cid:3039) is small, writing (cid:1843) (cid:3028) ~ επ (cid:1870) (cid:3028)(cid:2870) (cid:1832) (cid:3022)(cid:3023)/(cid:3006)(cid:3022)(cid:3023) , then Φ (cid:3006)(cid:3013) in Equation (5) can be roughly scaled by the luminosity if the solutions are subsonic or do not significantly exceed (cid:1843) (cid:3030) . Because the atmospheric expansion affects (cid:1870) (cid:3028) , iterative solutions can improve estimates of (cid:1843) (cid:3028) . Ignoring this, we note that for an early earth-like upper atmosphere dominated by N or O, a solar EUV flux more than 100 times the present would be required for the escaping gas to go sonic in the continuum regime using data from Tian et al. (2008) in Equation (10). It is, therefore, unlikely that escape from such an atmosphere on a super-earth in the habitable zone would have a sonic point in the continuum regime. Lammer et al. (2013) calculated the escape rate from a hydrogen atom thermosphere due to XUV radiation on super-earths observed orbiting close to their star. They used Parker’s upper boundary conditions to solve the continuum equations and then decided that blow-off occurred if (cid:2019) (cid:3051) (cid:3407) 3/2 (Öpik 1963). Rather than solving these equations and then deciding whether the sonic or Jeans conditions should have been used, Equation (10) can be used to estimate whether a sonic point might occur in the continuum regime for a given (cid:1843) (cid:3041)(cid:3032)(cid:3047) . For example, using (cid:1837)(cid:1866) (cid:3040) ~1 and (cid:1870) ∗ ~(cid:1870) (cid:2868) we find that (cid:1843) (cid:3030) ~ 0.6 (cid:1876)10 (cid:2869)(cid:2871) W and (cid:2869)(cid:2871)
W for Kepler11b and 11c respectively, using data in Table 1 of Lammer et al. (2013). These values can be compared to their heating rates obtained using ε (cid:3404) 0.15: (cid:1843) (cid:3041)(cid:3032)(cid:3047) ~1.2 and x (cid:2869)(cid:2871) W respectively. For this ε , Kepler11b has a sonic point (cid:4666)(cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3408) (cid:1843) (cid:3030) (cid:4667) in the continuum region, whereas Kepler11c does not (cid:4666)(cid:1843) (cid:3041)(cid:3032)(cid:3047) (cid:3407) (cid:1843) (cid:3030) ). Therefore, Φ (cid:3006)(cid:3013) in Equation (5) is most applicable to Kepler11c which requires kinetic, not sonic, boundary conditions to obtain for an accurate description of its upper atmosphere. Of course, increasing ε or reducing the gravitational energy due to the tides can change this. Summary
We have used results from DSMC simulations to show that the oft-used energy limited rate in Equation (5) for an isentropic expansion of a heated upper atmosphere is most reasonable for a subsonic expansion with a large escape rate. The accuracy depends on how well one estimates (cid:1843) (cid:31) (cid:3032)(cid:3047) . Conversely, agreement with the energy limited escape rate does not imply that sonic boundary conditions are applicable to continuum models of hermal escape. In fact the simple approximation in Equation (5) becomes less good with increasing (cid:1843) (cid:3041)(cid:3032)(cid:3047) above (cid:1843) (cid:3030) as seen in Figure 2. Although the size of the escape rate might not be strongly dependent on whether sonic, Jeans, or kinetic boundary conditions are used, the upper atmosphere can be significantly affected as seen in Figure 3. Therefore, past applications Parker’s (1964a,b) model have led to incorrect descriptions of the upper atmosphere when rapid escape occurs (Tucker and Johnson 2009; Tucker et al. 2012; Erwin et al 2013). Since the upper atmosphere structure affects the interaction of the escaping gas with the ambient plasma and with neighboring bodies, we have given an expression in Equation (10) to estimate when sonic boundary conditions are likely to be applicable in calculating the escape from and the expansion of the upper atmosphere of a planetary body.
Acknowledgement
We acknowledge support from NASA’s Planetary Atmospheres Program
Appendix