Comparative planetary nitrogen atmospheres: Density and thermal structures of Pluto and Triton

Icarus 291 (2017) 55–64

Contents lists available at ScienceDirect

Icarus journal homepage: www.elsevier.com/locate/icarus

Comparative planetary nitrogen atmospheres: Density and thermal structures of Pluto and Triton Darrell F. Strobel a,∗, Xun Zhu b a b

Departments of Earth & Planetary Sciences and Physics & Astronomy, Johns Hopkins University, 3400 N. Charles Street, Baltimore, MD 21218, United States The Johns Hopkins University Applied Physics Laboratory, 11100 Johns Hopkins Road, Laurel, MD 20723, United States

a r t i c l e

i n f o

Article history: Received 7 October 2016 Revised 2 March 2017 Accepted 10 March 2017 Available online 18 March 2017 Keywords: Atmospheres Structure Pluto Atmosphere Satellites Atmospheres Triton Radiative transfer

a b s t r a c t Both atmospheres of Pluto and Neptune’s largest satellite Triton have cold surfaces with surface gravitational accelerations and atmospheric surface pressures of comparable magnitude. To study their atmospheres we have updated Zhu et al. (2014) model for Pluto’s atmosphere by adopting Voigt line profiles in the radiation module with the latest spectral database and extended the model to Triton’s atmosphere by including additional parameterized heating due to the magnetospheric electron transport and energy deposition. The CH4 mixing ratio profiles play central roles in differentiating the atmospheres of Pluto and Triton. On Pluto the surface CH4 mole fraction is in the range of 0.3-0.8%, sufficiently high to ensure that it is well mixed in the lower atmosphere and not subject to photochemical destruction. Near the exobase CH4 attains comparable density to N2 due to gravitational diffusive separation and escapes at 500 times the N2 rate (= 1 × 1023 N2 s−1 ). In Triton’s atmosphere, the surface CH4 mole fraction is on the order of 0.015%, sufficiently low to ensure that it is photochemically destroyed irreversibly in the lower atmosphere and that N2 remains the major species, even at the exobase. With solar EUV power only, Triton’s upper thermosphere is too cold and magnetospheric heating, approximately comparable to the solar EUV power, is needed to bring the N2 tangential column number density in the 50 0–80 0 km range up to values derived from the Voyager 2 UVS observations (Broadfoot et al., 1989). Due to their cold exobase temperatures relative to the gravitational potential energy wells that N2 resides in, atmospheric escape from Triton and Pluto is not dominated by N2 Jeans escape but by CH4 from Pluto and H, C, N and H2 from Triton. The atmospheric thermal structure near the exobase is sensitive to the atmospheric escape rate only when it is significantly greater than 2 × 1027 amu s−1 , above which enhanced Jeans escape and larger radial velocity adiabatically cools the atmosphere to a lower temperature. Finally we suggest that Pluto’s thermosphere is a cold ∼ 70 K due to ablation of H2 O molecules from the influx of dust grains detected by New Horizons Student Dust Counter. © 2017 Elsevier Inc. All rights reserved.

1. Introduction In the outer solar system, Pluto and Triton are widely regarded as the largest end-members of Kuiper-Belt objects and as “twins” with thin buffered N2 atmospheres controlled by interactions with surface ice, primarily N2 frost, by the sublimation and condensation processes. As a result, the atmospheric surface pressure is tenuous and is close to nitrogen vapor pressure equilibrium with the cold surface temperature. During the Voyager flyby in 1989, Triton’s N2 surface pressure was 14 ± 1 μbar and consistent with a surface temperature ∼37.8 K (Gurrola, 1995), whereas in 2009 the N2 surface pressure had risen to ∼40 μbar and surface temperature ∼39 K (Lellouch et al., 2010). From the New Horizons space∗

Corresponding author. E-mail address: [email protected] (D.F. Strobel).

http://dx.doi.org/10.1016/j.icarus.2017.03.013 0019-1035/© 2017 Elsevier Inc. All rights reserved.

craft flyby and the radio occultations performed by the Radio Experiment (REX) instrument, Gladstone et al. (2016) reported a surface pressure of ∼10-12 μbar and surface temperature of 37 ± 3 K on Pluto. Although they have small surface pressures, these icy planetary bodies with small surface gravities have extended atmospheres with exobase altitudes for Pluto ∼1700 and Triton ∼850 km. Their atmospheric radial density and thermal structures can be well described by fluid models in a great spatial domain and may also be significantly influenced by the atmospheric escape processes (Strobel, 2002; Zhu et al., 2014). Zhu et al. (2014) developed a one-dimensional (1D) radiativedynamical model for Pluto’s vertical density and thermal structure that includes both the atmospheric thermal energy sources associated with the radiative heating and cooling rates and dynamical energy sink of the adiabatic cooling associated with the radial motion induced by the atmospheric escape processes. The

56

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

model was significantly updated to simulate Pluto’s density and thermal structure as observed by New Horizons (Gladstone et al., 2016). In this paper, we provide a description of the updated radiation module of the Pluto model that adopts the more accurate Voigt line profiles and uses the latest spectral database. By including an additional heating term by the magnetospheric electron transport and energy deposition for Triton’s atmosphere the updated 1D radiative-dynamical model is used to investigate the vertical density and thermal structure and the associated escape rate for Triton’s atmosphere. Section 2 describes the improvements in our updated model for the icy planetary bodies. For readers only interested in results, they can skip Section 2 and go directly to Section 3 that presents and discusses the vertical density and thermal structure of Triton’s atmosphere for several reference cases. In Section 4 we compare Triton’s atmosphere with Pluto’s and Titan’s. Section 5 summarizes the paper. 2. Numerical model of the vertical density and thermal structure for an icy planetary body’s atmosphere

2.2. Improved spectral integration of the k-distribution coefficients and the escape function based on a Voigt line profile

2.1. Radiative-dynamical atmosphere model including escape processes Zhu et al. (2014) developed a 1D radiative-dynamical model that includes the atmospheric escape processes for an icy planetary body such as Pluto’s. The total energy balance equation for an expanding atmosphere under a gravitational potential is given by

ρ cp

    Fρ ∂ 1 2 GM ∂T 1 ∂ ∂T = 2 r2 κ + Rnet − 2 u + c pT − , (1) ∂t ∂r r r ∂r r ∂r 2

where T is the temperature of an ideal gas with cp being the specific heat at constant pressure, ρ is the air density, r is radius, κ = κ0 T s is the thermal conductivity with κ0 = 5.63 × 10−5 J m−1 s−1 K−(1+s ) and s = 1.12 (Hubbard et al., 1990). Rnet is the net radiative heating shown in Strobel et al. (1996) but updated by using the latest molecular kinetics (Zalucha et al., 2011) and including the additional heating of solar far ultraviolet and extreme ultraviolet (FUV-EUV) absorption by N2 and CH4 constituents. In this work, we have also updated the k-coefficients and the calculations of the escape functions by use of the Voigt line profiles in the radiation module based on the latest spectral database (Rothman et al., 2013). The three terms in the brackets of Eq. (1) are the total energy for an ideal gas with u being the bulk radial velocity of the gas, G the Newton’s gravitational constant and M the mass of either Pluto or Triton. Finally, Fρ is the mass escape rate that is an integration constant of the steady state continuity equation to relate air density and radial velocity (Zhu et al., 2014):

r 2 ρ u = Fρ .

(2)

Eq. (1) needs to be solved with the steady-state radial momentum equation

1 ∂p −1 = 2 p ∂r r



Fρ GM + RT p

 ∂u , ∂r

(3)

where pressure p, T, ρ are related through the ideal gas law

p = ρ RT .

escape rate Fρ determined by matching the altitude of the upper boundary to the exobase (Zhu et al., 2014). We have modified the numerical model by flexible parameter settings so that it can be used to simulate the vertical thermal structure for either Pluto or Triton’s atmosphere. By setting the parameters such as planetary orbit, mass, and radius to those for Pluto, the Pluto-Triton model is used to explain the vertical density and thermal structure for Pluto’s atmosphere as observed by the New Horizons (Gladstone et al., 2016). In addition to the two basic parameters of the planetary mass and radius that are distinctively different between Pluto and Triton, there is an additional heating term by the magnetospheric electron transport and energy deposition for Triton’s atmosphere that has also been included (Strobel et al., 1990b; Stevens et al., 1992) in the updated Pluto-Triton model.

(4)

Eqs. (1)-(4) are integrated to a steady state under constant radiative forcing and a set of fixed boundary conditions. The lower boundary condition is the fixed surface temperature and pressure. The physically self-consistent upper boundary condition is a parameterized escape rate at the exobase derived from the direct simulation Monte Carlo kinetic particle model, i.e., the type-I boundary condition in Zhu et al. (2014). The four unknowns (T, ρ , u, p) are solved as functions of r by an iteration procedure with the

In (Strobel et al., 1996), we assumed Doppler broadening in calculating the k-distribution coefficients and escape function because the Voigt parameter, d ≡ 2α L /α D with α L and α D being the half widths of Lorentzian and Doppler line profiles, respectively, was at most 0.001 for an air pressure less than 10 μbar so that the absorption is mostly contributed from the core regions. Currently, we have adopted the latest HITRAN2012 spectral database (Rothman et al., 2013) that contains many more lines in hot bands of weaker line strengths. It is found that line parameters for interpolating the line strength into low temperature regime have also been improved in HITRAN2012. Hot bands of weaker line strengths could make contributions in absorption comparable to those from the wing regions of absorption lines with the strong line strengths. As a result, we recalculate the k-distribution coefficients by adopting a Voigt line profile. The Voigt line profile calculation is based on an algorithm by Fomichev and Shved (1985). The derived lookup tables of the k-distribution coefficients are now functions of both T and p. Furthermore, we have also assumed a Voigt line profile in calculating the escape function of CO and HCN rotational lines and included 72 and 141 rotational lines shown in HITRAN12 (Rothman et al., 2013), respectively. The treatment of numerical integration for each rotational line in calculating the escape function is described as follows. The atmospheric cool-to-space cooling rate of an absorption line is proportional to its escape function that is defined as (Andrews et al., 1987; Zhu, 2004)

(b, p, p ) =



∞

−∞

fν ( p)E2 [b×gν ( p, p )]dν,

(5)

where b is the product of line strength and the absorber amount, fν is the line profile, E2 (x) is an exponential integral, gν is an integrated or an averaged line profile between pressures p and p . It remains nearly a Doppler profile if α D > >α L . It becomes Voigt profile when α D ≤ α L (Goody and Yung, 1989). For CO rotational lines, α D ∼ α L mainly because α D decreases linearly with center wavenumber ν 0 whereas α L is almost independent of ν 0 (Goody and Yung, 1989). Thus, fν (p) is a Voigt line profile with significant cooling rate contributions from its Lorentzian wings. Note that if we let the normalized fν and gν approach a Dirac delta function in Eq. (5), then the escape function  approaches 0 for any finite value of b. In reality, the composite half width of a line profile will always be limited by a finite α D even if under a special case of α L → 0 as p → 0. Here, we approximate gν (p, p ) by a Voigt line profile with averaged α¯ D and α¯ L between p and p .

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

We then get

   b, p, p =



α D (α L ∼α D or d∼1) may lead to an increase in cooling rate by, ∞

−∞



fν [αD (T ( p) ), αL ( p, T ( p) )] ×

E2 b × f ν or

(b, αD , αL , α¯ D , α¯ L ) =



    αD T , αL p, T dν,



∞ −∞

(6)

fν (αD , αL )E2 [b fν (α¯ D , α¯ L )]dν.

(7)

An enclosed analytic expression can be derived if we assume a simple Voigt line profile with a rectangular core plus a set of Lorentzian wings (Fels, 1979):

fν (αD , αL ) =

|ν| ≤ ν0 , |ν| > ν0

C,

αL /(π ν 2 ),

(8)

where the tuned profile parameters are given by

1 αL + 1.5αD , C = − , π 2 ν0 π ν02 4α¯ 1 α¯ L ν¯ 0 = L + 1.5α¯ D , C¯ = − . π 2ν¯ 0 π ν¯ 02

ν0 =

57

4αL

(9)

Substituting Eqs. (8) and (9) into Eq. (7), we obtain

 αL ¯ (m, αD , αL , α¯ D , α¯ L ) = 2 ν˜0C E2 (mC ) + √ erf(x0 ) , 2 πτ

(10)

√ where ν˜0 = (ν0 + ν¯ 0 )/2, x0 = τ /ν˜0 , τ = ηbα¯ L /π with η = 2.1 being the diffusivity factor for cooling rate calculations (Apruzese, 1980), and finally, erf(x ) is the error function that can be efficiently evaluated (Lether, 1990). Eq. (10) can be evaluated adequately with peak fractional errors ranging from 5% to 30% for the Voigt parameter d ranging from 0.1 to 10 when the cooling rate becomes significant. A more accurate but less efficient method to evaluate Eq. (7) is a direct numerical integration. Setting ν1 = 2(αD + αL + α¯ D + α¯ L ), we can evaluate the escape function numerically in two frequency intervals (0, ν 1 ) and (ν 1 , ∞) with 8 and 12 points, respectively,

(b, αD , αL , α¯ D , α¯ L ) = 2

 ν1

fν E2 [b f¯ν ]dν

0



+2

1/ν1

0

(1/ν 2 ) fν E2 [b f¯ν ]d (1/ν ),

(11)

or

(b, αD , αL , α¯ D , α¯ L ) =

8 ν1 

4

i=1

12 1  fi E2 [b f¯i ] + f j E2 [b f¯j ]/ν 2j . 6 ν1 j=1

(12) A slightly better approximation of numerical integration is by the Gaussian quadratures where both the weights (wj ) and abscissas (xj ) are optimally chosen (Press et al., 1992). We may select 4 and 6 points in the frequency intervals (0, ν 1 ) and (ν 1 , ∞), respectively, to obtain

(b, αD , αL , α¯ D , α¯ L ) = 2

4  i=1

wi fi E2 [b f¯i ] + 2

6 

w j f j E2 [b f¯j ]/ν 2j .

j=1

(13) In Fig. 1a, we show the plots of escape function  versus absorber amount b under three different Voigt parameters of d = 0.1, 1, and 10. Panels (b)-(d) are the corresponding fractional errors when approximations (10), (12) and (13) are used to evaluate  . The ratios of the central processing times for computing  in Panels (b)-(d) are approximately 1:5:3, respectively. Fig. 1 suggests that, for a given α D , inclusion of the effect α L as it gets close to

say, from ∼5% to a factor of ∼3, all depending on whether the absorption of a particular line is in a weak or a strong regime. It is also noted that fractional errors of the simple analytic expression (10) become significantly great only for a nearly pure Doppler line (d << 1) with strong absorption when the rectangular core given by (8) is no longer a good approximation and the cooling rate as scaled by the escape function is quite small. Eq. (10) is routinely used in calculating Rnet in the current model except at the final stage of iteration where Eq. (12) is used for calculating the cooling rate by CO and HCN rotational lines. 2.3. Parameterization of diffuse scattered solar radiation from Pluto’s surface At the suggestion of Leslie Young, we have included backscattered, diffuse solar radiation off Pluto’s surface to boost the CH4 near-IR bands heating in the lower atmosphere of Pluto. Our parameterization follows Strobel (1978) with a Lambert surface having an assumed globally averaged albedo of 0.8 (Buratti et al., 2016), independent of wavelength for our purposes. As discussed below this increases the CH4 near-IR bands overall heating rate by ∼ 30%, mostly in the 2.3 and 1.6 μm bands. 3. The vertical density and thermal structure of Triton’s atmosphere One of the main focuses in this paper is to examine the vertical density and thermal structure of Triton’s atmosphere. Lellouch et al. (2010) detected CO in Triton’s atmosphere from its (2-0) band at 2.35 μm, with a column density of ∼8 × 1018 cm−2 . For their estimated 40 μbar surface pressure, the corresponding CO mixing ratio would be ∼6 × 10−4 . During the Voyager 2 Triton flyby ultraviolet solar occultations revealed that CH4 had a smaller scale height (∼9 km) than N2 did (∼20 km) as a consequence of photochemistry (Strobel et al., 1990a) and when extrapolated to the surface yielded an averaged surface density of ∼4.7 × 1011 cm−3 and surface mixing ratio of 0.0 0 015 (cf. Herbert and Sandel, 1991; Strobel et al., 1990a; Strobel and Summers, 1995; Krasnopolsky and Cruikshank, 1995). In July 2009, Lellouch et al. (2010) inferred that the CH4 surface partial pressure was about 4 times larger than during the Voyager encounter with a surface CH4 density ∼1.9 × 1012 cm−3 from which we would infer a somewhat larger CH4 scale height (∼15 km). Due to the higher surface pressure the surface CH4 mixing ratio (0.0 0 024) was only slightly higher than during the Voyager 2 flyby. Since the same radiative module has been well tested while simulating Pluto’s atmospheric vertical structure, our focus will be on the effect of the heating rate by the magnetospheric electron transport and energy deposition. We first run the Pluto-Triton 1D radiative-dynamical model for three cases (Triton-1, Triton-2 and Triton-3) for Triton’s atmosphere with the basic model parameters shown in Table 1. The Triton-1 case includes only solar heating in the EUV, FUV, and near-IR (no magnetospheric heating) with the CH4 mixing ratio profile inferred from the Voyager 2 UV occultations. The Triton-2 case investigates what the required CH4 constant mole fraction would have to be if solar heating alone (no magnetospheric heating) could explain the N2 density and thermal structure inferred from the UV occultations, while Triton-3 adds an intermediate amount of magnetospheric heating by the magnetospheric electron transport and energy deposition to Triton-1 case. The latter heating rate has been parameterized in Table 1 by a factor hmag with a vertical profile following a model calculation by Stevens et al. (1992). The Triton-1 & 2 cases shown in Table 1 correspond to vanishing (hmag = 0), and intermediate (hmag = 0.73) magnetospheric heating for Triton-3. The criteria for

58

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

Fig. 1. Escape function ( ) and fractional errors based on approximations given by Eqs. (10) labeled Fels, (12) labeled Sum1 and (13) labeled Sum2, respectively. Table 1 Basic parameters for three cases of model runs for Triton atmosphere.

Triton - 1 Triton - 2 Triton - 3

Ts (K)

ps (μbar)

dS−P (AU)

g0 (m s−2 )

r(CH4 )

r(CO)

hmag

zexo (km)

Fn (s−1 )

37.8 37.8 37.8

16.0 16.0 16.0

30.3 30.3 30.3

0.779 0.779 0.779

1.5 × 10−4 0.9 × 10−4 1.5 × 10−4

6 × 10−4 6 × 10−4 6 × 10−4

0 0 0.73

696 858 884

5.28 ×1016 1.85 ×1019 5.21 ×1019

Triton-2 case CH4 constant mole fraction value and the value of hmag for Triton-3 were a best fit to the measured tangential N2 column density reported in Stevens et al. (1992). It should be noted that the Triton-2 constant CH4 mole fraction case is invalidated by the Voyager 2 solar and stellar occultation data (Strobel et al., 1990a) and also by physics and chemistry unless Triton’s lower atmosphere were highly turbulent and vigorously mixed, which would be improbable as it is very stably stratified. This is easily understood by noting that the integrated CH4 photochemical loss rate at 30 AU is ∼2 × 108 cm−2 s−1 (Strobel et al., 1990a). The ability to replace dissociated CH4 depends on the rate that CH4 can be diffusively transported upward. There is an upper limit known as the Hunten limiting flux (Hunten, 1973) ϕ l = binary collision coefficient (∼2 × 1018 cm−1 s−1 ) times CH4 mixing ratio times (1 − 16/28) divided

by N2 scale height (∼20 km). For 0.009% CH4 , the Triton-2 value, the maximum upward diffusive flux (ϕ l ∼4.5 × 107 cm−2 s−1 ) fails by almost an order of magnitude. But if the CH4 mixing ratio were increased to Pluto’s value ∼ 0.6–0.8%, then ϕ l ∼ 3 × 109 cm−2 s−1 and almost a factor of 10 larger than the integrated CH4 photochemical loss rate, so that diffusion can overwhelm chemical loss. Fig. 2 shows three cases of temperature and density profiles calculated with the input quantities given in Table 1. The exobase altitudes (zexo ) correspond to the altitude where the atmospheric Knudsen number is 0.7, which is the equivalent of the classic definition of the exobase where √ the mean free path is equal to the scale height divided by 2 (Strobel, 2002), and are listed in Table 1. We first note that the vertical temperature gradient in the lower thermosphere (<100 km) is substantially reduced in comparison with Pluto’s lower atmosphere (Zhu et al., 2014;

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

59

Fig. 2. Model (a) temperature and (b) density profiles of Triton’s atmosphere for the model parameters given in Table 1.

Fig. 3. Solar near infrared and FUV-EUV-Magnetosphere heating rate profiles for Triton’s atmosphere.

Gladstone et al., 2016). This is mainly due to the much lower value of r(CH4 ) in Triton’s atmosphere at the surface (Table 1) and to the fact that its density decreases more rapidly than N2 ’s density (Table 1), whereas in Pluto’s atmosphere r(CH4 ) ∼ 0.0 06–0.0 08 and well mixed in the lower atmosphere (Gladstone et al., 2016; Lellouch et al., 2017). As a consequence Triton’s radiative heating in the 1.7 μm, 2.3 μm and 3.3 μm CH4 bands is correspondingly reduced in comparison to Pluto’s lower atmospheric radiative heating (Zhu et al., 2014). The much lower temperatures below ∼200 km in all Triton cases (Fig. 2a) and the absence of a peak temperature ∼110 K just 30 km above the surface as on Pluto render Triton’s atmosphere more compact and compressed. We show in Fig. 3 the vertical profiles of the solar near-IR and FUV-EUV heating rates for Triton’s atmosphere’s three cases. The FUV-EUV heating rate for Triton-3 includes the heating by the magnetospheric electron transport and energy deposition with a secondary peak in the lower atmosphere around 60–130 km (Stevens et al., 1992) due to ionospheric produced N atoms diffusing downward to recombine in termolecular reactions which release chemical and vibrational energy eventually converted to heat. The large peak

Fig. 4. Tangential N2 column number density profiles for Triton’s atmosphere based on density profiles shown in Fig. 2b. The measured values are marked by square boxes.

heating at ∼15 km is due to FUV, mostly solar Lyman-α , absorption by CH4 . In Fig. 4, we plot the Voyager 2 Ultraviolet Spectrometer (UVS) tangential N2 column number density as a function of altitude for Triton’s atmosphere (Broadfoot et al., 1989; Stevens et al., 1992), for comparison with our model calculations. The figure suggests that our best fit model requires the magnitude of the magnetospheric heating rate to be reduced by ∼30% in comparison with the estimation made in (Stevens et al., 1992), as they did not include near-IR CH4 band heating. But Fig. 4 does confirm their conclusion that a solar heating only model is incompatible with the UVS data for N2 above 450 km (Triton-1) and for CH4 below 100 km (Triton-2, Strobel et al., 1990a). Furthermore, without the contribution from the magnetospheric heating component the model produces a much colder Triton’s atmosphere with its tangential column number density in the 50 0–80 0 km altitude range being a factor of 10 smaller than derived from the measurements.

60

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

(a)

(b)

(c)

Fig. 5. Comparison of the density, pressure, and temperature profiles of Pluto, Titan, and Triton based on New Horizons data, Huygens Probe data, and Voyager 2 data, respectively. For Titan, the profiles start at an altitude of 400 km (right y axis), whereas for Pluto and Triton the left axis is the true altitude, for (a) temperature, (b) pressure, and (c) densities. Note that red temperature profile (5a) has been slightly updated from (Gladstone et al., 2016) to represent value and location of peak temperature and the shallow temperature minimum between 30 0–50 0 km as inferred from ALMA data (Lellouch et al., 2017).

4. Comparison with the other N2 atmospheres of Pluto and Titan In Fig. 5 the density, pressure, and temperature profiles of the atmospheres of Pluto, Titan, and Triton are compared. For Titan, the Huygens atmospheric structure instrument (HASI) profiles extracted from the Huygens Probe entry deceleration data (Fulchignoni et al. 2005) were used starting at altitude of 400 km, where the pressure was ∼14 μbar and comparable to the surfaces pressures on Pluto and Triton. For Pluto, the New Horizons profiles were obtained from the radio and solar occultation data reported in Gladstone et al. (2016) with a slightly updated temperature profile, while for Triton our theoretical Triton-3 model that was constrained by the Voyager 2 solar UV occultation data was adopted. Note that two altitude scales are used in the figures: one, for absolute altitude above the surface for Pluto and Triton plots and relative altitude for Titan corresponding to z = 0 at 400 km and two, for absolute altitude for Titan on right y axis. There is an ordering in upper atmosphere isothermal temperatures of 170:90:70 K for Titan, Triton, Pluto (Fig. 5a). When combined with the corresponding ratios of surface gravities, 80:47:34 cms−2 , at altitudes of 80 0:50 0:40 0 km, the base of their isother-

mal regions, they have comparable ratios of T/g and hence, as nitrogen atmospheres, the same scale height. This suggests that the major energy balance between the input solar radiative energy and its redistribution within the thermosphere through diffusive processes occurs over the same altitude range measured in terms of the scale height. The pressure profiles of Pluto and Titan (Fig. 5b) are very similar, whereas Triton with its colder lower atmosphere has a similar shape due to its comparable ratio of T/g, but two orders of magnitude lower pressures. Fig. 5c shows that the most distinctive difference among these N2 atmospheres is the CH4 density profile of Triton. As noted previously the photochemical destruction of CH4 below 100 km results in a very cold lower atmosphere and overall a more compact atmosphere. The N2 and CH4 density profiles of Titan and Pluto are remarkably similar over their entire atmospheres below ∼14 μbar, except at very high altitudes where CH4 is going into gravitational diffusive equilibrium on Pluto (Gladstone et al., 2016), but for reasons not well understood CH4 remains fairly well mixed to high altitudes, as if it were escaping rapidly from Titan (Cui et al., 2012; Strobel and Cui, 2014). Escape from these atmospheres is dominated by minor species, not N2 . From Titan, the largest escape rate is 1 × 1028 H2 s−1 , followed by (1.1–1.7) × 1027 H s−1 (Strobel and Cui, 2014), and

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

a highly variable and uncertain ∼3 × 1027 CH4 s−1 (Cui et al., 2012). On Triton, escape is dominated also by minor species: 5 × 1025 H s−1 , followed by 2 × 1025 H2 s−1 , 6 × 1024 N s−1 , (Strobel and Summers, 1995) to be compared with 2.4 × 1025 H s−1 , 4.5 × 1025 H2 s−1 , 7.7 × 1024 N s−1 , and 1.1 × 1024 C s−1 , (Krasnopolsky and Cruikshank, 1995). On Pluto, CH4 leads escape with a rate of 5 × 1025 CH4 s−1 , whereas the heavier N2 escapes at only 1 × 1023 N2 s−1 , (Gladstone et al., 2016). In the Zhu et al. (2014) model, the atmosphere is assumed to be fully compressible with no restrictions on the Mach number of the radial flow. The model also assumes that atmospheric escape is dominated by the major species. Their predicted escape rates for Pluto were large, but the radial velocity remained in the low Mach number regime (<0.05 at the exobase). Thus Pluto’s atmosphere was predicted to be hydrostatic. With the escape rates derived from New Horizons data which are lower than predicted, one can estimate the adiabatic cooling rate by rewriting the last term in (1) in the hydrostatic limit (u → 0) as,

Qad (r ) = kb

N N  Fρ i dT  c pi Fρ i + g( r ) → 1.1 × 10−13 (1700/r )4 , 2 dr mi r r2 i

i

(14) where kb is Boltzmann’s constant and cpi is the heat capacity in units of kb , e.g., 3.5 for N2 and 4 for CH4 . Above z = 500 km, the expression on the RHS of Eq. (14) is accurate to better than 5%. At Pluto’s exobase, the adiabatic cooling rate from Eq. (14) is ∼1 × 10−14 erg cm−3 s−1 with the Gladstone et al. (2016) escape rates. Based on radiative heating and cooling rates at the exobase, adiabatic cooling becomes important when escape rates exceed ∼1 × 1026 CH4 s−1 or more generally ∼2 × 1027 amu s−1 . Pluto’s upper atmosphere is a cold 70 K as shown in Fig. 5a. It is isothermal because the escape rate is low and adiabatic cooling is inconsequential. Gladstone et al. (2016) suggested that the combination of non-LTE cooling by C2 H2 13.7 μm ν 5 band emission and HCN LTE rotational line emission (if HCN were supersaturated) might explain the inferred temperature profile. From HCN observations with ALMA, Lellouch et al. (2017) showed that even though HCN was supersaturated by some seven orders of magnitude in the upper atmosphere it still was not sufficient to maintain a 70 K atmosphere, but instead led to the Zhu et al. (2014) solution of a high N2 escape rate (>1 × 1027 s−1 ) and large adiabatic cooling. We suggest that the additional “unknown” cooling agent is most likely H2 O molecules in spite of its very low vapor pressure at Pluto atmospheric temperatures. They are probably a volatile component of the dust grains intercepted by Pluto and detected by New Horizons’ Student Dust Counter (Poppe, 2015; Horanyi et al., 2015). They estimate the mass influx of dust grains to be about 50 kg per day. But only a small fraction of the grains’ mass must be H2 O molecules as discussed below. LTE H2 O rotational line emission is modeled with the HITRAN2012 database. The adopted H2 O mixing ratio is tweaked to yield a temperature profile acceptably close to the inferred temperature profile in Gladstone et al. (2016), as shown in Fig. 5a. In Fig. 6a the model heating rates for the nonLTE CH4 3.3, 2.3, 1.7, and 7.6 μm vibrational bands, the latter where it is absorbing blackbody radiation from the 100–110 K region, FUV, and N2 EUV are shown. The CH4 bands dominate heating up to an altitude of ∼ 360 km, where CH4 absorption of solar FUV radiation becomes more important and at ∼ 700 km, N2 absorption of solar EUV radiation is equally important. As in Zhu et al. (2014), heating efficiencies of 0.5 for CH4 and 0.25 for N2 were adopted. The corresponding cooling rates are shown in Fig. 6b for CO LTE rotation line emission with the CO constant mixing ratio r(CO) = 0.0 0 05 (Lellouch et al., 2017), CH4 7.6 μm vibrational band, C2 H2 ν 5 band emission, and H2 O LTE rotational line emission. The C2 H2 number density profile was constructed and constrained by

61

its line of sight column density profile in Gladstone et al., (2016). The C2 H2 collisional deactivation by N2 coefficient adopted was 4.4 × 10−14 exp(−(T−240)/105) cm3 s−1 (normalized to Hager et al., 1981 with CH4 ν 4 band temperature dependence). CO cooling plays a minor role except near the surface and is overtaken by C2 H2 13.7 μm band cooling up to ∼ 100 km, where the decreasing temperature from 110 to 70 K shifts the peak of the Planck function to ∼ 40 μm and C2 H2 band emission occurs from the tail of function. Above 50 km, H2 O rotational line emission dominates cooling. To demonstrate the potential importance of H2 O we present in Fig. 7 temperature profiles calculated with and without H2 O along with a comparison of the H2 O and HCN mixing ratio profiles for the “with H2 O” case. Without H2 O or an alternative cooling agent, one obtains the Zhu et al. (2014) solution of a high N2 escape rate (>1 × 1027 s−1 ) and large adiabatic cooling. To make the case for H2 O molecules playing an important role as a cooling agent in the atmosphere, we adopted the method of Strobel (2009) to solve the continuity and diffusion equations numerically for N2 , CH4 , and H2 O with a fourth order Runge–Kutta algorithm and all boundary conditions specified at the surface. The ablation source rate of H2 O was modeled as proportional to 2 exp[−( (r − 1690 )/80) ] based on Horanyi et al. (2015) with an integrated source rate of 3.4 × 1019 H2 O s−1 . In the calculation integrated source rate is much more important than the radial distance of peak ablation and shape of ablation profile. The boundary conditions for H2 O were a downward flow rate to the surface of 3.4 × 1019 H2 O s−1 in lieu of explicit condensation processes and an integrated enhanced escape rate at the exobase of ∼ 1 × 1014 H2 O s−1 , governed by the Hunten limiting escape rate (Hunten, 1973; Chamberlain and Hunten, 1987, p. 370). For the illustrative calculation, the eddy diffusion coefficient was a constant 104 cm2 s−1 . In Fig. 8, the 1D calculation of the H2 O density profile with an ablation source is compared with the required H2 O density profile to cool the atmosphere. Also shown for reference are the corresponding N2 and CH4 profiles. In practice to construct this H2 O profile, also shown in Fig. 7, to precisely yield the inferred temperature profile proved difficult. One should interpret the H2 O profile structure with a significant degree of smoothing. If the dust grain mass flux of 50 kg day−1 were all H2 O molecules, this would correspond to 1.94 × 1022 H2 O s−1 to be compared with only 3.4 × 1019 H2 O s−1 needed in Fig. 8. Dust grains as carriers of H2 O would appear to be adequate to supply the required H2 O, but, of course, subject to a rigorous calculation of how supersaturated H2 O could be above z = 50 km. Thus the dust grains could be mostly silicates with only a fraction of their mass comprising of volatiles such as H2 O. Note in Fig. 5a that the calculated temperature profile has a peak value of ∼ 106 K at 37 km. The altitude of the peak is strongly impacted by the adopted H2 O density and cooling profiles. The REX occultation temperature profile in Gladstone et al. (2016) had a peak value of ∼ 111 K at ∼ 32 km and the stellar occultation of Sicardy et al. (2016) has a peak temperature ∼ 110 K at 28 km. Additional heating may still be needed just above the surface and one possibility is the lowest, bright haze layer hugging the surface (cf. Fig. 4 in Gladstone et al., 2016) if it has an absorption optical depth of ∼ 4 × 10−5 over bulk of the solar constant (equal to 1260 erg cm−2 s−1 during Pluto encounter). Increasing only the surface pressure from 11 to 13 μbar and holding all other input quantities the same, even with the surface CH4 mole fraction at 0.8%, did not alter the location of the peak temperature. 5. Discussion and conclusions We updated the 1D radiative-dynamical Pluto model for Pluto’s atmosphere (Zhu et al., 2014) by adopting Voigt line profiles in calculating the k-distribution coefficients and line escape func-

62

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

(a)

(b)

Fig. 6. (a) Heating rates of Pluto’s atmosphere due to the non-LTE CH4 3.3, 2.3, 1.7, and 7.6 μm vibrational bands, FUV band, and N2 EUV band. (b) Cooling rates of Pluto’s atmosphere due to CO rotational lines, C2 H2 13.7 μm, CH4 7.6 μm vibrational bands, HCN rotational lines, adiabatic cooling by CH4 escaping, and H2 O rotational line cooling needed to account for the New Horizons occultation inferred temperature profile.

(a)

(b)

Fig. 7. (a) Temperature profiles of Pluto’s atmosphere with and without cooling rate by H2 O rotational lines. (b) Mixing ratio profiles for HCN and H2 O used in the Pluto’s model atmosphere.

tions with the latest spectral database (Rothman et al., 2013). With the inclusion of parameterized heating due to the magnetospheric electron transport and energy deposition, we applied the model to Triton’s atmosphere. The vertical density and thermal structure and the associated escape rate of Triton’s atmosphere are sensitive to its CH4 mixing ratio profile that determines where solar near-IR and FUV radiation is absorbed. In the Voyager epoch this solar power was deposited in the lower atmosphere. But with solar power only, the calculated upper thermosphere temperature is still too cold. Magnetospheric heating, approximately comparable to the solar EUV power, is needed to bring the N2 tangential column number density in the 50 0–80 0 km range up to values derived from the Voyager 2 UVS observations (Broadfoot et al., 1989). A stellar occultation light curve recorded with Fine Guidance Sensor (FGS) #3 aboard the Hubble Space Telescope in November 1997 by Elliot et al. (1998) yielded an isothermal region at T ∼ 50 K over the altitude region 25–50 km in Triton’s atmosphere. Elliot et al. (20 0 0) had difficultly explaining this occultation temperature profile. From Figs. 2 and 5a, one might initially conclude

that Triton has a “significant” positive temperature gradient there, but certainly not of the magnitude in Pluto’s lower atmosphere. But actually for the Triton best case (Triton-3) at 25 km, T = 51 K and at 50 km, 53.5 K, effectively isothermal. At 80 km, the calculated value T = 57 K is quite close to the inferred value ∼ 60 K. Note that initial conditions for the inversion were applied at 90 km. Our updated model for CH4 near-IR heating resulted in a more “isothermal region” than calculated in Elliot et al. (20 0 0). Due to their cold exobase temperatures relative to the gravitational potential energy well that N2 resides in, atmospheric escape from Triton and Pluto is not dominated by N2 Jeans escape but by minor species. The atmospheric thermal structure near the exobase is sensitive to the atmospheric escape rate only when it is significantly greater than 2 × 1027 amu s−1 above which enhanced Jeans escape and larger radial velocity adiabatically cools the atmosphere to a lower temperature. Our choice of H2 O as the cooling agent needed to maintain Pluto’s thermosphere at ∼ 70 K was made after rejecting other molecules. If other nitriles cool by LTE rotational line emission and

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

63

References

Fig. 8. Number density profiles for N2 , CH4 and H2 O derived from a onedimensional photochemical model plus a required H2 O density profile that produces a temperature profile needed to account for the New Horizons occultation inferred profile.

have a dipole moment comparable to HCN (= 2.99 Debye) that is typical of many nitriles that range from HNC (= 3.05 D) to CH3 CN (= 3.92 D), then their gas density would have to exceed the ALMA HCN density by ∼ a factor of 20 at high altitudes and ∼ a factor of 200 below 300 km, while not condensing out into haze particles. The required H2 O column density is smaller than ALMA HCN density above z = 300 km (cf. Fig. 7b). By way of comparison, the temperature profile calculated with the ALMA HCN density profile yields a mesopause at z ∼ 530 km with T ∼ 80 K, an approximately isothermal upper thermosphere with peak T ∼ 108 K at z = 20 0 0 km (Lellouch et al., 2017; cf. their Fig. 13), decreasing to 10 0 K at zexo ∼ 450 0 km with an N2 escape rate of ∼ 2.3 × 1027 N2 s−1 . In Fig. 7 we compare the temperature profiles calculated with and without H2 O rotational line emission. Without H2 O emission, the exobase temperature increases to 99 K at zexo ∼ 4670 km and the N2 escape rate is ∼ 2.6 × 1027 N2 s−1 . Zhu et al. (2014) predicted a slightly higher escape rate of ∼ 3.5 × 1027 N2 s−1 with exobase at zexo ∼ 7050 km where T = 69 K. The N2 Jeans λ parameter was ∼ 5 in both cases. The key difference is the newer calculations had C2 H2 13.7 μm ν 5 band emission which leads to a colder, more compact atmosphere with lower exobase but compensated by higher exobase temperature to provide comparable escape rates and adiabatic cooling in lieu of an unknown cooling agent. Finally in terms of comparative N2 atmospheres, Pluto’s atmosphere is a much closer twin to upper atmosphere of Titan above 400 km in terms of the N2 , CH4 density and pressure profiles than to Triton’s atmosphere.

Acknowledgments This research was supported by the New Horizons Mission through SWRI Contract No. 277043Q and NASA Grant NNX10AB84G from the NASA Planetary Atmosphere Program. DFS acknowledges beneficial discussions with the New Horizons Atmospheric Team in many telecoms. In particular, we thank Leslie Young for the suggestion of including backscattered, diffuse solar radiation off Pluto’s surface to boost the CH4 near-IR bands heating.

Andrews, D.G., Holton, J.R., Leovy, C.B., 1987. Middle Atmosphere Dynamics. Academic Press, New York, p. 489 pp. Apruzese, J.P., 1980. The diffusivity factor reexamined. J. Quant Spectrosc. Radiat. Transfer 24, 461–470. Broadfoot, A.L., Atreya, S.K., Bertaus, J.L., Blamont, J.E., Dessler, A.J., Donahue, T.M., Forrester, W.T., Hall, D.T., Herbert, F., Holberg, J.B., Hunten, D.M., Krasnopolsky, V.A., Linick, S., Lunine, J.I., McConnell, J.C., Moos, H.W., Sandel, B.R., Schnerder, N.M., Shemansky, D.E., Smith, G.R., Strobel, D.F., Yelle, R.V., 1989. Ultraviolet spectrometer observations of Neptune and Triton. Science 246, 1459–1466. Buratti, B.J., Hofgartner, J.D., Hicks, M.D., et al., 2016. Global albedos of Pluto and Charon from LORRI New Horizons observations. Icarus doi:10.1016/j.icarus.2016. 11.012. Chamberlain, J.W., Hunten, D.M., 1987. Theory of Planetary Atmospheres: An Introduction to Their Physics and Chemistry. Academic Press, p. 481 pp. Cui, J., Yelle, R.V., Strobel, D.F., M¨uller-Wodarg, I.C.F., Snowden, D.S., Koskinen, T.T., Galand, M., 2012. The CH4 structure in Titan’s upper atmosphere Revisited. J. Geophys. Res. E11. doi:10.1029/2012JE004222. Elliot, J.L., Hammel, H.B., Wasserman, L.H., Franz, O.G., McDonald, S.W., Person, M.J., Olkin, C.B., Dunham, E.W., Spencer, J.R., Stansberry, J.A., Buie, M.W., Pasachoff, J.M., Babcock, B.A., McConnochie, T.H., 1998. Global warming on Triton. Nature 393, 765–767. Elliot, J.L., Strobel, D.F., Zhu, X., Stansberry, J.A., Wasserman, L.H., Franz, O.G., 20 0 0. The thermal structure of Triton’s middle atmosphere. Icarus 143, 425–428. Fels, S.B., 1979. Simple strategies for inclusion of Voigt effects in infrared cooling rate calculations. Appl. Opt. 18, 2634–2637. Fomichev, V.I., Shved, G.M., 1985. Parameterization of the radiative flux divergence in the 9.6 μm O3 band. J. Atmos. Terr. Phys. 47, 1037–1049. Fulchignoni, M., Ferri, F., Angrilli, F., 2005. In situ measurements of the physical characteristics of Titan’s environment. Nature 438 (Dec.), 785–791. doi:10.1038/ nature04314. Gladstone, G.R., Stern, S.A., Ennico, K., Olkin, C.B., Weaver, H.A., et al., 2016. The atmosphere of Pluto as observed by New Horizons. Science 351 (6279). doi:10. 1126/science.aad8866. Goody, R.M., Yung, Y.L., 1989. Atmospheric Radiation: Theoretical Basis. Oxford Univ. Press, New York and Oxford, p. 519. Gurrola, E.M., 1995. Interpretation of Radar Data from the Icy Galilean Satellites and Triton Ph.D. thesis. Stanford University. Hager, W.Krieger, Pfab, J., 1981. Collisional Deactivation of Laser-excited Acetylene by H2, HBr, N2 and CO. J. Chem. Soc. Faraday Trans. 2 (77), 469–476. Herbert, F., Sandel, B.R., 1991. CH4 and haze in Triton’s lower atmosphere. J. Geophys. Res. 96, 19241–19252. Horanyi, M., Poppe, A., Sternovsky, Z, 2015. Dust ablation in Pluto’s atmosphere. AGU Abstract, Fall Meeting, San Francisco. Hubbard, W.B., Yelle, R.V., Lunine, J.I., 1990. Nonisothermal Pluto atmosphere models. Icarus 84, 1–11. Hunten, D.M., 1973. The escape of light gases from planetary atmospheres. J. Atmos. Sci. 30, 1481–1494. Krasnopolsky, V.A., Cruikshank, D.P., 1995. Photochemistry of Triton’s atmosphere and ionosphere. J. Geophys. Res. 100, 21271–21286. Lellouch, E., de Bergh, C., Sicardy, B., Ferron, S., Käufl, H.-U., 2010. Detection of CO in Triton’s atmosphere and the nature of surface–atmosphere interactions. Astron. Astrophys. 512, L8. http://dx.doi.org/10.1051/0 0 04-6361/201014339. Lellouch, E., Gurwell, M., Butler, B., Fouchet, T., Lavvas, P., Strobel, D.F., Sicardy, B., ´ Moullet, A., Moreno, R., Bockelee-Morvan, D., Biver, N., Young, L., Lis, D., Stansberry, J., Stern, A., Weaver, H., Young, E., Zhu, X., Boissier, J., 2017. Detection of CO and HCN in Pluto’s atmosphere with ALMA. Icarus 286, 289–307. doi:10.1016/j.icarus.2016.10.013. Lether, F.G., 1990. An elementary approximation for exp(x2)erfc(x). J. Quant. Spectrosc. Radiat. Transfer 43, 511–513. Press, W.H., Teukolsky, S.A., Vetterling, W.T., Flannery, B.P., 1992. Numerical Recipes in Fortran. The Art of Scientific Computing. Cambridge Univ. Press, p. 963. Poppe, A.R., 2015. Interplanetary dust influx to the Pluto–Charon system. Icarus 246, 352–359. doi:10.1016/j.icarus.2013.12.029. Rothman, L.S., Gorden, I.E. and collaborators, 2013. The HITRAN2012 molecular spectroscopic database. J. Quant. Spectrosc. Radiat. Transfer 130, 4–50. Sicardy, B.and 67 colleagues, 2016. Pluto’s atmosphere from the 2015 June 29 ground-based stellar occultation at the time of the New Horizons flyby, 2016. Astrophys. J. 819, L38. doi:10.3847/2041-8205/819/2/L38. Stevens, M.H., Strobel, D.F., Summers, M.E., Yelle, R.V., 1992. On the thermal structure of Triton’s thermosphere. Geophys. Res. Lett. 19, 669–672. Strobel, D.F., 1978. Parameterization of the Atmospheric Heating Rate from 15 to 120 km due to O2 and O3 Absorption of Solar Radiation. J. Geophys. Res. 83, 6225–6230. Strobel, D.F., 2002. Aeronomic systems on Planets, Moons, and Comets. In: Mendillo, M., Nagy, A., Waite, J.H. (Eds.), Atmospheres in the Solar System: Comparative Aeronomy. American Geophysical Union, Washington DC, pp. 7–22. Strobel, D.F., 2009. Titan’s hydrodynamically escaping atmosphere: Escape Rates and the Structure of the Exobase Region. Icarus 202, 632–641. doi:10.1016/j.icarus. 20 09.03.0 07.20 09. Strobel, D.F., Summers, M.E., Herbert, F., Sandel, B.R., 1990a. The Photochemistry of Methane in the Atmosphere of Triton. Geophys. Res. Lett. 17, 1729–1732. Strobel, D.F., Cheng, A.F., Summers, M.E., Strickland, D.J., 1990b. Magnetospheric interaction with Triton’s ionosphere. Geophys. Res. Lett. 17, 1661–1664.

64

D.F. Strobel, X. Zhu / Icarus 291 (2017) 55–64

Strobel, D.F., Summers, M.E., 1995. Triton’s upper atmosphere and ionosphere. In: Cruikshank, D.P. (Ed.), Neptune and Triton. University of Arizona Press, Tucson, pp. 1107–1148. Strobel, D.F., Zhu, X., Summers, M.E., Stevens, M.H., 1996. On the vertical thermal structure of Pluto’s atmosphere. Icarus 120, 266–289. Strobel, D.F., Cui, J., 2014. Titan’s upper atmosphere/exosphere, escape processes and rates. In: Mueller-Wodarg, I., Griffith, C., Lellouch, E., Cravens, T. (Eds.), Titan: Interior, Surface, Atmosphere and Space Environment, Cambridge Planetary Science Series. Cambridge University Press, pp. 355–375.

Zalucha, A.M., Gulbis, A.A.S., Zhu, X., Strobel, D.F., Elliot, J.L., 2011. An analysis of Pluto occultation light curves using an atmospheric radiative-conductive model. Icarus 211, 804–818. Zhu, X., et al., 2004. Radiative transfer in the middle atmosphere and planetary atmospheres. In: X., Zhu, et al. (Eds.), Observation, Theory and Modeling of Atmospheric Variability. World Scientific Pub. Co., Singapore, pp. 359–396. Zhu, X., Strobel, D.F., Erwin, J.T., 2014. The density and thermal structure of Pluto’s atmosphere and associated escape processes and rates. Icarus 228, 301–314.