P Lavvas

A Arfaux

Impact of photochemical hazes and gases on exoplanet atmospheric thermal structure

Keywords:

We investigate the impact of photochemical hazes and disequilibrium gases on the thermal structure of hot-Jupiters, using a detailed 1-D radiative-convective model. We find that the inclusion of photochemical hazes results in major heating of the upper and cooling of the lower atmosphere. Sulphur containing species, such as SH, S 2 and S 3 provide significant opacity in the middle atmosphere and lead to local heating near 1 mbar, while OH, CH, NH, and CN radicals produced by the photochemistry affect the thermal structure near 1 µbar. Furthermore we show that the modifications on the thermal structure from photochemical gases and hazes can have important ramifications for the interpretation of transit observations. Specifically, our study for the hazy HD 189733 b shows that the hotter upper atmosphere resulting from the inclusion of photochemical haze opacity imposes an expansion of the atmosphere, thus a steeper transit signature in the UV-Visible part of the spectrum. In addition, the temperature changes in the photosphere also affect the secondary eclipse spectrum. For HD 209458 b we find that a small haze opacity could be present in this atmosphere, at pressures below 1 mbar, which could be a result of both photochemical hazes and condensates. Our results motivate the inclusion of radiative feedback from photochemical hazes in general circulation models for a proper evaluation of atmospheric dynamics.

INTRODUCTION

Observations of hot-Jupiters at UV-visible wavelengths reveal a slope in their transit depth. While for some of these planets the observed slope is consistent with the anticipated molecular scattering, for others the required extinction along the line of sight is much higher, suggesting the presence of aerosols (Sing et al. 2016). Similarly, observations of Neptune-type exoplanets also suggest a variable degree of haziness, with its presence better manifested in transit depth measurements at the 1.4 µm H2O band (Crossfield & Kreidberg 2017). However, current data does not allow for a clear distinction between an aerosol layer or a high mean molecular weight atmosphere -a problem that may be ameliorated with future observations with broader wavelength coverage and/or higher precision. Stellar variability can potentially affect the interpretation of short wavelength observations (e.g. [START_REF] Rackham | [END_REF]). On the other hand, the persistence of the aerosol-indicator features ⋆ E-mail: panayotis.lavvas@univ-reims.fr (PL) over multiple star/planet systems suggests that aerosols are an important component of exoplanet atmospheres.

The exact nature of these aerosols is still under debate, with two types of heterogeneous components -photochemical hazes and condensates (clouds) -being the two main anticipated families (Morley et al. 2013(Morley et al. , 2015;;Kempton et al. 2017;Gao et al. 2020), reflecting our understanding of gas giant atmospheres established through observations of the solar system [START_REF] West | Jupiter, Jovian Clouds and Haze[END_REF][START_REF] West | [END_REF]. However, the composition of exoplanet aerosols could be drastically different from those in our neighbourhood. In exoplanets condensates include silicates and other exotic compounds that are present in the gas phase of these hot environments (Burrows & Sharp 1999;Wakeford & Sing 2015;Mbarek & Kempton 2016;[START_REF] Woitke | [END_REF], while photochemical hazes are expected to have a carbon (He et al. 2018;Fleury et al. 2019) or sulfur (Zahnle et al. 2009;Gao et al. 2017) based structure with possible inclusions of other atoms. Although the relative contributions of condensates and photochemical hazes is not yet fully clarified, the impact of such particles on the thermal structure could be significant.

Previous studies on the impact of aerosols on the atmospheric thermal structure of exoplanets have mainly focused on the role of condensates. Studies on hot-Jupiters considered condensates of MgSiO3 and Fe and found that they could have an impact on the thermal structure and the emitted thermal radiation, but its magnitude will vary among different planets and according to the assumptions made regarding the properties of the cloud particles (Fortney et al. 2005;Lines et al. 2018). Typically, the presence of optically thick clouds results in cooling of the atmosphere below the cloud deck and heating above, as the altitude range of photon absorption is modified accordingly. Nevertheless, the impact of cloud particles on hot-Jupiter's thermal structure seems to be moderate and most studies of cloud formation so far have focused on the reverse effect of how horizontal temperature variations could affect the cloud formation (Parmentier et al. 2013;Helling et al. 2019). For sub-Neptunes, studies on GJ 1214b demonstrate that the KCl and ZnS clouds anticipated in this atmosphere (Morley et al. 2013) generate a similar effect: temperatures, depending on the cloud particle radius considered, change by 50 K for 1 µm particles up to ∼300 K for 0.1 µm particles (Charnay et al. 2015). Investigations of the microphysical properties of cloud particles suggest radii of a few µm for silicate clouds in hot-Jupiters (Powell et al. 2018) and KCl/ZnS clouds in sub-Neptunes [START_REF] Gao | [END_REF]Ohno & Okuzumi 2018) at the pressure levels where the cloud optical depth is one, implying that a moderate effect is anticipated on the thermal structure from cloud particles in both planet types. On the other hand, more recent studies find significant radiative feedbacks by clouds on hot-Jupiters (Roman & Rauscher 2019;Roman et al. 2020).

Photochemical hazes are typically produced at lower pressures than those of cloud formation. Microphysics models for such particles suggest particle radii of a few nm in the upper atmosphere, based on the constraints imposed by the observed transit depths and the anticipated formation yields from the atmospheric photochemistry (Lavvas & Koskinen 2017;Gao et al. 2017;Kawashima & Ikoma 2018;Kawashima et al. 2019;Lavvas et al. 2019;Adams et al. 2019;Gao et al. 2020). Depending on their abundance and optical properties photochemical hazes can result in heating or cooling of the atmosphere (McKay et al. 1991;Arney et al. 2016;Zhang et al. 2017). Previous estimates on the role of photochemical hazes on GJ1214b based on particles with 0.1 µm radii suggest a significant radiative feedback (Morley et al. 2015). We explore further this interaction in our study using a radiative-convective model and detailed microphysics modelling of the photochemical haze formation. We focus only on the role of photochemical hazes, thus do not consider the contribution of cloud opacity in our simulations.

Most of the previous studies on exoplanet thermal structure have assumed thermochemical equilibrium for the abundances of gaseous species. However, detailed models of disequilibrium chemistry driven by atmospheric mixing and photochemistry (Zahnle et al. 2009;Moses et al. 2011;Lavvas et al. 2014;[START_REF] Vardavas | Radiation and Climate[END_REF] demonstrate that the divergence between the two conditions can be significant for the abundances of the species affecting the thermal structure. An evaluation of the impact of disequilibrium chem-istry on the thermal structure suggested temperature differences up to ∼100 K (Drummond et al. 2016). Thus, we further explore the implications of chemical disequilibrium along with the role of photochemical hazes.

As photochemical hazes are mainly affecting the upper atmosphere's opacity (p < 1 mbar), divergence from local thermodynamic equilibrium (LTE) conditions is anticipated. Here, we do not consider non-LTE effects for the gas phase, which although may be important, are beyond our current goal of evaluation of photochemical haze impact and require a separate dedicated study. However, we do consider the temperature disequilibrium between the particles and the gas envelope that is important in the upper atmosphere (Lavvas & Koskinen 2017).

Finally, we perform our simulations in 1D so that we can explore in detail the role of photochemical hazes and gases. Nevertheless, dynamics should affect the horizontal distribution of photochemical hazes (Komacek et al. 2019;Steinrueck et al. 2020). Our goal with this study is to evaluate the impact of these components on the thermal structure of the upper atmosphere and motivate the inclusion of their radiative effects in general circulation models of exoplanet atmospheres.

In Section 2 we present a description of the radiativeconvective model used in the study (2.1), the opacities and gas composition considered (2.2), the photochemical haze properties used (2.3) and how the dis-equilibrium between particles and gas temperature is evaluated (2.4). In Section 3 we present results for two reference case exoplanets, those of HD 209458 b (3.1) and HD 189733 b (3.2). As explained further below, for these calculations we consider as input disequilibrium chemical compositions and photochemical haze size distributions that were derived in previous studies, assuming a temperature profile from GCM simulations. In Section 4 we further explore how our simulated temperature profiles impact the interpretation of transit and secondary eclipse observations for these exoplanets, as well as, present a preliminary evaluation of how the revised temperature profiles may affect the disequilibrium chemistry and photochemical haze properties. Section 5 presents our final conclusions.

METHODS

This section provides an overview of the model used in the study, of the gaseous and particulate opacities considered, and how we treat the possible disequilibrium between the particles and the gas background.

Radiative-Convective model

We use a 1D radiative-convective model to simulate the exoplanet thermal structure from the deep atmosphere (10 3 bar) to the lower thermosphere (∼10 -10 bar). RC models provide rapid calculations for the atmospheric thermal structure and have been applied to the Earth's atmosphere (Vardavas & Carver 1984), atmospheres in the solar system (McKay et al. 1989;Guillot et al. 1994;Lavvas et al. 2008), as well as to exoplanet atmospheres (Burrows et al. 1997;Fortney et al. 2003;Iro et al. 2005;Amundsen et al. 2014;Drummond et al. 2016). The model used here was described in detail in its previous application to Titan's atmosphere (Lavvas et al. 2008), while further details on the radiative transfer methods used can be found in classical textbooks (e.g. Chandrasekhar 1960;[START_REF] Goody | Atmospheric Radiation[END_REF][START_REF] Vardavas | Radiation and Climate[END_REF]. As the topic is extensively covered in the existing literature, here we provide only a short description of the necessary components.

RC models separate the atmosphere into layers and solve for the vertical temperature profile considering two domains; one where radiative equilibrium dominates the atmospheric energy balance, and one deeper in the atmosphere, where convective mixing provides a more stable energy exchange relative to the radiative equilibrium. In the radiative balance regime, equilibrium requires that the net flux in each atmospheric layer is zero. Thus, the atmospheric temperature is adjusted until the net thermal emission balances the net stellar radiation within each layer. In the convective regime the temperature gradient follows the atmospheric lapse rate, Γ, and the temperature magnitude adjusts to the value that provides the same energy flux in the transition region, as the pure radiative equilibrium solution would. The location of the transition from a radiative to a convective regime depends on the planet/star case studied and is not certain to always occur within the pressure region under investigation here. Thus, we evaluate its presence based on the calculated temperature profile and adjust for it if the simulated temperature profile becomes steeper than the lapse rate. Previous studies suggested that detached convective layers may exist in the high pressure regimes of exoplanet gas-giants (Guillot et al. 1994;Burrows et al. 1997). However, such transitions occur at higher pressures than those considered here or when the effective temperature of the atmosphere is much smaller than those addressed in our study. Hence, we consider a single convective region in our work.

As implied in the above description, the typical approach for the evaluation of stellar and planetary (thermal) fluxes is to perform independent calculations that provide the flux at each atmospheric level for these two components. For the stellar energy deposition we consider wavelengths ranging between 0.5 Å to 10 µm, in an expanding grid of 660 wavelength bands. Similarly for the thermal fluxes we consider a grid of 600 wavenumber bands of 50 cm -1 width extending from 0 to 30000 cm -1 . Although the radiative transfer equation can be solved analytically in the case of no scattering (Chandrasekhar 1960), the contribution of molecular scattering at short wavelengths and haze scattering over the whole spectrum necessitates numerical solvers. Here we use the delta-Eddington approximation and solve for the three moments of the radiative transfer equation for both stellar and thermal radiation, thus allowing for the calculation of the net flux at each atmospheric level [START_REF] Vardavas | Radiation and Climate[END_REF].

We perform our nominal calculations considering a disk average geometry, thus we divide the incoming stellar radiation by a factor of 4 to take into account the diurnal and horizontal averaging. Although the planets under consideration are expected to be tidally locked to their parent star, we find that disk average insolation conditions are in better agreement with 3D circulation results regarding the thermal structure, compared to cases where the incoming radiation is averaged only over the day-hemisphere (division by factor of 2). For the stellar fluxes simulation we furthermore consider that any radiation reaching the lower boundary is locally absorbed. Practically near the lower boundary both direct and diffuse stellar components attain very small flux values due to the local high opacity. For the thermal radiation fluxes we consider no incoming flux at the top of the atmosphere. At the lower boundary we consider a net outgoing flux characterised by a black body temperature Tint that represents energy released from the processes occurring in the interior of the gas giants. Such an addition was found necessary for explaining the observed mass-radius ratios of multiple exoplanets (Guillot & Showman 2002), with typical values of Tint ranging between 100 and 300 K (Iro et al. 2005;Fortney et al. 2008). Subsequent studies suggest that Tint may attain higher values for hot-Jupiters (Thorngren et al. 2019), however in the simulations presented below we consider Tint = 100 K to facilitate comparison with previous thermal structure evaluations.

Gaseous opacities

Disequilibrium chemical abundances are calculated with a kinetic model extending from 10 3 bar to 10 -10 bar that solves the continuity equation for each species taking into account atmospheric mixing (eddy), molecular diffusion, and (photo)chemical reactions (Lavvas et al. 2014). The eddy profile used is based on general circulation model results [START_REF] Showman | [END_REF]. For the radiation transfer, our gaseous opacity list includes numerous atomic and molecular absorbers the contributions of which varies significantly across the spectrum. At high energies (λ<0.3µm), our scheme is basically an inheritance from the photochemistry calculations for exoplanet atmospheres, thus it includes absorption cross sections from multiple photochemically active species. At longer wavelengths opacities are calculated based on the correlated-k method (Lacis & Oinas 1991) utilising molecular line lists from the Exomol, Theorets, and MoLLIST groups. These include opacities for H2O (Polyansky et al. 2018), CO (Li et al. 2015), CH4 (Rey et al. 2017), NH3 (Coles et al. 2019), H2S (Azzam et al. 2016), CO2 Rothman et al. (2010), HCN (Barber et al. 2013), C2H2 (Chubb et al. 2020), TiO (an updated in 2012 version of Plez 1998), VO (McKemmish et al. 2016), He + 3 (Mizus et al. 2017) and HeH + (Amaral et al. 2019). In the current simulations we do not treat condensation and cloud formation thus do not include opacities for silicate species. As our simulations benefit from chemical composition abundances calculated for disequilibrium chemistry, we additionally included opacities for abundant radicals present in the upper atmosphere from the photolysis of major molecules. These include OH (Yousefi et al. 2018), CH (Masseron et al. 2014), NH (Fernando et al. 2018), SH [START_REF] Gorman | [END_REF] and CN (Brooke et al. 2014). The cross sections of these radicals present strong absorption bands at visible wavelengths (Fig. 1), thus we investigate their possible contribution to the heating of the upper atmosphere. Na and K opacities are calculated as in Lavvas et al. (2014), while bound-free and free-free opacities by H -and H - 2 are based on John (1988), Bell & Berrington (1987), and (Bell 1980). Furthermore we include Rayleigh scattering by the most abundant species, Thompson scattering by free electrons and collision induced absorption by H2-H2 and H2-He pairs.

Haze impact

Haze particles can efficiently absorb and scatter electromagnetic radiation therefore affect the radiation field in the atmosphere. The degree of this interaction depends on the size distribution of the particles and their optical properties (refractive index). For the particle size distribution we use results from our previous photochemical haze microphysics studies for the atmosphere of HD 189733 b (Lavvas & Koskinen 2017), which provided particle size distributions able to satisfy the observed transit depth constraints. In a nutshell, these studies evaluated the particle size distribution assuming a photochemical haze mass flux in the upper atmosphere the magnitude of which was based on photochemistry calculations. The model solves for the particle size distribution (explicitly described in size bins) by solving the continuity equation for each particle size taking into account particle settling, atmospheric mixing (eddy) and particle coagulation and ablation. For the considered photochemical haze optical properties (see below), the particle mass flux is estimated based on the photolysis of the major precursors CH4, HCN and C2H2 assuming a yield for the production of photochemical hazes that is adjusted to fit the observed transit spectra. Further details on the transit spectra simulation and the evaluation of the photochemical haze yield based on observations and laboratory measurements can be found in the relevant studies (Lavvas & Koskinen 2017;Lavvas et al. 2019).

The refractive index of the photochemically produced haze particles in exoplanet atmospheres is currently unknown as the available observational constraints are not sufficient to allow for the retrieval of this parameter, while no laboratory estimations on this parameter have yet been published. Titan's photochemical hazes for which both observational and experimental constraints are available (see review by [START_REF] West | Titan[END_REF], are often used as analog for exoplanet photochemical hazes. However, Titan's haze pyrolise at temperatures below 600 • C (Israel et al. 2005;Morisson et al. 2016), which is far lower than the anticipated temperatures in hot-Jupiters. Thus, this type of material cannot form or survive at such conditions. A different chemical structure is anticipated for exoplanet photochemical hazes, which would impose further changes in the resulting refractive index of the material. This conclusion leads to the suggestion that a soot-type composition is a more likely representative for hot exoplanets as these materials have a much higher resistance to elevated temperatures (Lavvas & Koskinen 2017). Major differences exist among these two compositions, with the soot composition providing a significantly stronger absorption compared to Titan's analogs, across the whole spectrum (Fig. 2). The properties of exoplanet photochemical hazes may sufficiently change among different planets, reflecting changes in the atmospheric chemical composition and the energy input to the atmosphere (He et al. 2020). As we cannot with certainty constraint this parameter yet, we consider these two (Titan and soot) composition types as limiting cases that potentially envelope the range of photochemical haze compositions in exoplanet atmospheres.

Disequilibrium aspects

In the interaction of haze with radiation, particles are heated therefore their temperature (Tp) increases above that of their surrounding environment (Tg). The rate at which the excess particle energy is transferred to the atmosphere depends on the efficiency of collisions of gas species on the particle surface, as well as, the atmospheric opacity at the wavelengths where the excess (thermal) energy of the particles is radiated. The energy balance at steady state for a single spherical particle of radius r, takes the form:

πr 2 QνF ⊙ ν dν-4πr 2 Qν[πBν (Tp)-πBν(Tg)]dν-Ėcol = 0 (1)
where the first term is the rate at which the particle absorbs stellar radiation, the second is the net rate of energy loss through thermal emission, and the third is the rate of energy exchange through collisions with gas molecules. In the above, Qν is the absorption efficiency of the spherical particle, which through Kirchhoff's law, equals its emission effi- ciency, F ⊙ ν is the attenuated stellar flux at the location under investigation, and Bν is the black body emission. The collisional term, Ėcol , depends on the atmospheric density. At locations where the gas mean free path, λg, is much larger than the particle radius (free molecular regime, Kn=λg/r>> 1) and the particle velocity is small compared to the thermal velocity of the gas, the exchange rate is given by [START_REF] Gombosi | Gaskinetic Theory[END_REF]:

Ėcol = αT πr 2 γ + 1 γ -1 1 2 pVT 1 - Tp Tg ( 2
)
where αT is the thermal accommodation coefficient (assumed unity in our calculations), γ the heat capacity ratio, p the gas pressure and VT the gas thermal velocity. At the opposite limit of continuous flow (Kn << 1), the energy exchange rate is described through the thermal conduction, κ, of the gas:

Ėcol = 4πrκ(Tp -Tg) (3)
There is no analytical solution able to describe the transition between these two regimes. Instead, various approximations have been suggested. For example, [START_REF] Sherman | A Survey of Experimental Results and Methods for the Transition Regime of Rarefied Gas Dynamics[END_REF] suggested a simple harmonic mean formula Ė = Ėc Ėf /( Ėc+ Ėf ) to describe the transition region, while similar approaches where further suggested by other researchers based on experimental measurements (see review by Filippov & Rosner 1999). Here, we use an interpolation formula given by [START_REF] Pruppacher | Microphysics of clouds and precipitation[END_REF]:

Ė = 4πr κ δ (Tp -Tg) with δ = r r + λg + VT κ αT rρcp (4)
which converges to the above two limiting cases at the appropriate conditions of Kn values. For small temperature differences between the particle and the gas the interpolation formulae are in agreement with detailed stochastic simulations (Filippov & Rosner 1999).

The above equations allow for the evaluation of the particle equilibrium temperature at different atmospheric levels. As the gas collision term becomes less efficient with decreasing pressure, particles warm up and their excess energy is balanced by their thermal emission. Previous evaluations for hot-Jupiters have demonstrated how the particle temperature starts to diverge from the gas temperature at pressures below ∼10 mbar and with increasing magnitude as particle radius increases. Thus, in the middle and upper atmosphere the particle heating efficiency decreases significantly below unity (Lavvas & Koskinen 2017). However, using this particle heating efficiency as a multiplicative factor of the haze extinction could be misleading for the evaluation of the atmospheric temperature. Even if the collisional energy transfer to the gas phase is small, stellar photons at different energies are absorbed by the haze, thus are not available to heat deeper levels of the atmosphere. In addition, the thermally emitted photons from the particles could be absorbed by the gas (and haze) envelope at different atmospheric locations instead of escaping the atmosphere, therefore could further lead to heating.

To properly address the impact of disequilibrium between the gas and particle temperatures on the thermal structure we included the excess radiative energy of the hot particle in the radiative transfer scheme and simulated the thermal energy fluxes deposited at different atmospheric levels only from this component. The source function at each atmospheric level is calculated from:

Sν = 1 χν p [Bν (Tp) -Bν (Tg)]σpnp (5)
with χν the local extinction coefficient and σp and np the emission cross section and number density of size bin p particles. The calculated thermal fluxes from the particles are subsequently added to the solar and planetary components for the evaluation of the net flux.

APPLICATION TO REFERENCE CASES

Here we describe the application of the model to the reference exoplanet atmospheres of HD 209458 b and HD 189733 b. As described above we use as input for the gas and photochemical haze composition results from our previous evaluations. However, as explained below, these studies considered as input thermal structure profiles derived from independent evaluations that assumed chemical equilibrium and clear atmosphere conditions (i.e. no hazes). Therefore, the results presented below serve to demonstrate the implications of disequilibrium chemistry and photochemical haze opacity on the thermal structure and to motivate the need for a self-consistent description in the future.

HD 209458 b

We start with the case of an apparently clear atmosphere, that of HD 209458 b. This planet has been a paradigm of exoplanetary investigations and has been the focus of multiple studies of thermal structure due to availability of observations (Guillot & Showman 2002;Iro et al. 2005;Fortney et al. 2006;[START_REF] Showman | [END_REF]). Therefore, it provides an opportunity to validate our model against observations and compare with previous studies, before investigating the role of photochemical hazes. HD 209458 is a G0V star with similar spectral characteristics to the sun. Thus, we consider a similar spectrum to the sun scaled to the stellar radius of HD 209458 for the stellar insolation. For the atmospheric chemical composition, we consider our previous evaluation of disequilibrium composition (Lavvas et al. 2014) that assumed solar elemental abundances. These computations where performed using the day-average temperature profile evaluated from a general circulation model [START_REF] Showman | [END_REF]). For clarity these inputs are presented as dashed lines in Figs. 3 &7.

A demonstration of the different species contributing to the atmospheric optical depth at different pressure levels reveals the major opacity sources (Fig. 4). We focus here mainly on the wavelength range relevant to the deposition of stellar energy. At pressures below 1µbar and λ<0.1µm the atmospheric opacity is driven mainly by H2O and atomic/molecular hydrogen (not shown), while at longer wavelengths we find that H -provides the strongest absorption contribution, which however has a small optical depth (τ <10 -6 ) and is smaller than the Thompson scattering opacity by electrons. Near 1 µbar, the densities of primary radicals (dominated by OH, CH and CN) become significant enough to provide a clear signature in the opacity at visible wavelengths (however they remain too small to be observable), while in the near IR the CO contribution is evident. Deeper in the atmosphere, at 10 µbar, the absorption bands of H2O in the near IR start to have a clear signature in the total opacity, while the Na and K signatures in UV (photoionization) and visible also clearly appear.

At higher pressures, another radical that starts to contribute near 0.1 mbar is SH that is formed from the photolysis of H2S. The contribution of SH was considered potentially important according to the preliminary evaluations by Zahnle et al. (2009). However, their SH cross section was based on reported line widths from the solar chromosphere that where subsequently revised and the novel constraints to the SH line profiles results in a significantly weaker cross section (by factors >10 at visible wavelengths [START_REF] Gorman | [END_REF]). Nevertheless, even with the revised cross section SH does have a dominant contribution to the opacity near 0.3 µm at pressures below ∼0.1 mbar. Similarly, S2 is another major opacity source below 1 mbar, centred in the boundary between UV and visible radiation. Previous studies suggested the theoretical cross section evaluated by Heijden & Mullen (2001). Yet these calculations correspond to high pressure plasma conditions that are not representative of the relevant atmospheric conditions under investigation here. Instead, we use a cross section evaluated from an experimentally derived linelist and applied to higher temperatures appropriate to exoplanet conditions (R.V. Yelle personal communication).

Moving to deeper pressure levels of HD 209458 b's atmosphere the main opacity sources remain the same but their individual contributions vary. Near 1 bar H2O dominates over the near-IR, Na and K in the visible and SH, S2 near 0.3 µ. Contributions by H2S are also important, while at p > 1 bar opacity by collision induced absorption starts to be significant.

The above opacity distribution imposes a direct signature to the stellar fluxes penetrating the atmosphere (Fig. 5). At λ < 900 Å (the H ionization limit) photons are absorbed at very low pressures (p < 1µbar), while at longer wavelengths the H2O, SH, and S2 define the penetration depths in the FUV. Photons at these energies reach progressively deeper pressure levels down to ∼1 mbar at 0.2 µm and ∼1 bar at 0.3 µm. At visible wavelengths Na and K dominate the penetration of radiation, while at near-IR H2O dominates the extinction of the stellar fluxes with stellar photons completely absorbed by ∼10 bar.

For the planetary emission, the corresponding thermal opacity is dominated by H2O at most IR wavelengths, while as discussed above, Na and K opacities dominate at visible. The planetary emission spectrum is characteristic of a black body of ∼1500 K with clear signatures of H2O bands and Na and K absorption at the wings and emission at the cores of the atomic lines (Fig. 6). Our results demonstrate that the majority of the outgoing thermal radiation originates from pressure levels between 1 bar and 1 mbar (Fig. 6). The situation only changes at the cores of the alkali lines which are optically thick at low pressure levels and their emission is significantly affected by the upper atmosphere temperature. The scattered stellar light flux leaving the atmosphere (Fig. 6) becomes comparable to the thermally emitted radiation near the K line and exceeds the thermal component at shorter wavelengths.

The balance between the stellar and planetary fluxes results in the simulated thermal profile (Fig. 7). We first compare our calculations with previous 1D studies considering thermochemical equilibrium composition with condensation effects included in the abundances of gaseous components. Specifically we compare with profiles generated under the assumption of complete redistribution (4π) of the incoming radiation and an internal temperature boundary condition of Tint=100K (Iro et al. 2005;Fortney et al. 2006). Our simulated temperature profile is similar to the previous evaluations in the lower atmosphere demonstrating that differences in the chemical composition between equilibrium and dis-equilibrium conditions do not have a major impact on the thermal structure below ∼ 0.1 bar for this planet.

At lower pressures however, the contribution of SH and S2 visible opacities based on the disequilibrium chemistry induce a significant atmospheric heating between 10 mbar and 0.1 mbar, which was not present in previous evaluations. Without these opacities, the temperature near 1 mbar Fortney et al. 2006). The purple line presents our simulated temperature profile without the contribution of photochemical radicals in the simulations (see text). All these simulations assume a 4π redistribution of incoming stellar energy and T int =100 K. The green line presents our nominal thermal structure assuming energy deposition on the day side only (2π redistribution). The black dashed line presents the thermospheric temperature profile from escape models (Koskinen et al. 2013), and the grey dotted lines the equilibrium saturation curves for MgSiO 3 and Mg 2 SiO 4 (Visscher et al. 2010). For all simulations, convection does not occur in the simulated pressure range. The inset compares our nominal 2π (black) and 4π (green) temperature profiles with the profiles retrieved by Line et al. (2016) based on secondary eclipse observations. The red and blue profiles with associated uncertainties are based on different parameterizations for the temperature structure used in their analysis. Right: Temperature profiles including TiO/VO at equilibrium density (see text). The thick lines present our simulations for 4π (blue) and 2π (red) stellar energy redistribution and the dashed line corresponds to the day-average temperature profile from 3D simulations [START_REF] Showman | [END_REF], which is the temperature profile used for the evaluation of the disequilibrium chemistry.

drops by ∼ 100K relative to the nominal profile, thus approaching further the Iro et al. (2005) profile. At pressures below 1µbar our profile demonstrates a sharp temperature increase with temperatures reaching values above 2500 K, which brings our simulated profile close to evaluations of the thermospheric temperature structure (Koskinen et al. 2013). This rise reflects the loss of the main atmospheric coolant, H2O, at these pressures. Our simulation does not include all the involved processes at these low pressures (e.g. cooling by advection and radiative recombination) for our results to be directly compared to the escape simulations. We note however that the radicals produced by the photolysis of the main species do affect to a small degree the thermal structure in the lower thermosphere. Our calculations suggest that removing their opacity contribution results in a maximum ∼120 K drop in the atmospheric temperature near 10µbar relative to the nominal profile. Thus the contribution of these species, as well as, of other radiative active molecules may affect the evaluation of the lower thermosphere conditions in escape simulations.

We also compare our calculations with the results from 3D simulations. For this comparison we consider the day-average profile derived from the GCM results of [START_REF] Showman | [END_REF]. The latter authors considered a thermochemical equilibrium composition but furthermore included the abundances of TiO and VO without any condensation effects. The contribution of TiO and VO is particularly important in the upper atmosphere as there the equilibrium abundances demonstrate a maximum, which in combination with the strong absorption bands at visible wavelengths for these molecules facilitates a strong local energy deposition. This effect is demonstrated by the dotted line in Fig. 4. As anticipated, considering the TiO-VO equilibrium densities based on the GCM temperature profile along with the disequilibrium chemical composition for the other species, our simulated temperature demonstrates a strong increase in the upper atmosphere with a corresponding cooling of the lower atmosphere, relative to the nominal profile (Fig. 7). The resulting thermal structure is rather similar to the GCM results. Our simulated profile assuming complete redistribution is slightly cooler than the day-average results from the GCM simulation. Assuming that the energy is deposited only on the day side (2π) we obtain a profile that is substantially warmer than the GCM results. Thus, we conclude that assuming a complete redistribution for the 1D simulations provides more representative results for the average temperature profile. Nevertheless, as we discuss below, the secondary eclipse observations are better reproduced assuming energy deposition on the day side.

Note that the thermal structure below 10 bar would require an extensively long simulation time for convergence, which is unrealistic for a GCM. Thus, the temperature profile there is forced to converge to an 1D simulation (Fortney et al. 2006), which however is warmer from the 1D results presented in the nominal case. This suggests that the GCM profile at p > 10 bar probably corresponds to a higher Tint value, thus cannot be directly compared to our results. Moreover, the GCM simulations do not extend to pressures smaller than 2µbar and the presented profile was interpolated with a hydrodynamic escape simulation temperature profile (see Moses et al. 2011, for details). The local minimum in our simulated profile near 1µbar is strongly affected by CO emission.

Our results demonstrate that the differences between Figure 9. Average particle radius and corresponding particle density in the atmosphere of HD 189733 b for two temperature profiles corresponding to the cases with and without the impact of soot type photochemical hazes. The particle mass flux is set at 10 -11 gcm -2 s -1 .

equilibrium and disequilibrium composition can have major implications for the thermal structure of the upper atmosphere, thus affect the local chemical composition. We explore these implications further in section 4.

HD 189733 b

HD 189733 b is considered a paradigm of a hazy hot-Jupiter. For the evaluation of its thermal structure we consider the disequilibrium chemical composition and photochemical haze size distribution derived previously (Lavvas & Koskinen 2017). These chemical composition results were based on a general circulation temperature profile and eddy mixing coefficient [START_REF] Showman | [END_REF]) assuming thermochemical equilibrium (with condensation, i.e. no TiO and VO) and solar elemental abundances. The photochemical haze particle size distribution was derived based on the observed transit depth and a detailed microphysics model. The input chemical composition and the temperature profile used for its evaluation are shown with dashed lines in Figs. 8 & 12, respectively. The particle average properties are shown in Fig. 9. HD 189733 is a K star and we consider for its spectrum a composite of ǫ-Eridani and a PHOENIX model at the stellar parameters (Lavvas & Koskinen 2017).

As the qualitative results for this atmosphere are similar to that of HD 209458 b, we focus here on the impact of disequilibrium chemistry and photochemical hazes.

The gaseous opacity distribution of HD 189733 b is similar to that of HD 209458 b presented earlier. However, for HD 189733 b more molecules contribute substantially to the total opacity such as NH3 (near 0.2µm) and S3 (near 0.4µm), while HCN and CH4 contribute in the H2O bands troughs in the near-IR (Fig. 10). This difference is a result of the survival of H2O to higher altitudes in this atmosphere, allowing other species such as CH4 and NH3 to survive (through UV screening and through reduction of atomic hydrogen abundance), and contribute to the atmospheric opacity. Moreover, the photolysis of the latter species leads to a higher population of primary radicals such as OH, CH, NH, SH and CN in the upper atmosphere of HD 189733 b compared to that of HD 209458 b, with a corresponding signature on the opacity spectrum.

Nevertheless, the major difference of HD 189733 b compared to the opacity of HD 209458 b is the photochemical haze contribution. The latter is the dominant opacity source in the upper atmosphere for pressure levels above 0.1 mbar, while it significantly contributes to the visible opacity between the atomic/molecular lines down to ∼1 bar (Fig. 10). This evaluation is based on the soot refractive index. However, assuming the Titan-type refractive index does not modify qualitatively the picture; although the total opacity reduces, haze remains the dominant absorber in the upper atmosphere.

Considering only the gaseous composition, our simulated profile for complete redistribution (4π) and Tint = 100 K is similar to previous evaluations of thermochemical equilibrium (black lines in Fig. 12). Specifically, our calculations are in good agreement in the lower atmosphere below 10 mbar, with equilibrium composition 1D simulations (Fortney et al. 2006), as well as with the day-average profile from the [START_REF] Showman | [END_REF] general circulation model. However, as observed in the case of HD 209458 b, the disequilibrium contributions of SH, S2 (and also S3 for HD 189733 b) provide a significant opacity that drives a local temperature increase centred at 0.5 mbar, which was not present in previous studies. At the peak, the temperature difference with the case where SH, S2 and S3 opacities are not included (grey line in Fig. 12) is ∼220 K. Moreover, it is important to highlight that the difference in the abundances of these species between equilibrium and disequilibrium evaluation is much larger than in the case of HD 209458 b, and results in a more clear signature in the penetrating stellar fluxes in the atmosphere of HD 189733 b near 0.3 µm (Fig. 11). Therefore, in this cooler planet consideration of disequilibrium processes has a pronounced impact on the thermal structure.

At pressures above 0.1µbar the lack of efficient atmospheric cooling due to the destruction of H2O results in a rapid temperature increase as identified in hydrodynamic escape models (Koskinen et al. 2013).

Inclusion of the haze (soot) opacity in the thermal structure simulation results in major heating of the upper atmosphere at pressures lower than 0.1 bar and with the maximum temperature increase reaching ∼400 K at 1µbar (solid red line in Fig. 12). As more photons are absorbed in the upper atmosphere, heating of the lower atmosphere is reduced resulting in a temperature drop of ∼130 K in the isothermal part of the atmosphere near 10 bar. Assuming the Titan-type refractive index, the corresponding temperature changes are still significant (solid blue line in Fig. 12). In this case the presence of haze particles leads to a similar temperature increase of ∼400 K near 1 µbar, but the vertical extend of the temperature rise is smaller and the profile merges with the clear atmosphere profile below ∼1 mbar. This is a direct consequence of the lower absorptivity of Titan-type material relative to soot that allows more of the visible/near-IR radiation to penetrate to deeper atmospheric layers. As a result the corresponding cooling of the Figure 12. Temperature profiles for HD 189733 b for 4π (left) and 2π (right) energy redistribution. Solid lines present our simulated profiles of disequilibrium chemistry without hazes (black), with soot opacity haze (red), and with Titan-type opacity haze (blue). The grey line presents the case of no haze and without the opacity from sulfur species (see text). Previous studies from 1D (dash-dotted, Fortney et al. 2006) and 3D (dashed, [START_REF] Showman | [END_REF]) simulations are also presented, with the latter corresponding to the temperature profile used for the evaluation of the dis-equilibrium chemical composition. The dash-dotted lines present the average particle temperature demonstrating the disequilibrium between particles and gas, while the dotted lines present the resulting atmospheric temperature if the disequilibrium between particles and gas is not taken into account. For all simulations convection is limited to the first layer near the lower boundary (10 3 bar).

lower atmosphere is limited to ∼45 K relative to the clear atmosphere case. Nevertheless, within the range of optical properties considered here, we find that photochemical hazes will have a significant impact on the atmospheric thermal structure.

We must highlight at this point that we considered the same particle size distribution for both refractive index temperature evaluations. This size distribution derived from the microphysics model was constrained by fitting the observed transit spectrum with a particle mass flux of 5×10 -12 g cm -2 s -1 and considering the soot refractive index and a temperature profile that did not include the photochemical haze impact (Lavvas & Koskinen 2017). As discussed below, assuming the Titan-type refractive index with the same size distribution results in a transit spectrum that falls short of the observed, therefore a higher particle mass flux would need to be considered for this type of particle composition. The required higher opacity will increase further the role of such particles in heating the upper atmosphere, as identified above for the soot composition.

The strong interaction of stellar radiation with the haze particles leads to their rapid heating (Fig. 12). For the soot composition assumption, the average particle temperature diverges from the gas temperature at p < 0.1 mbar and remains ∼100 K higher than the atmospheric temperature (with haze heating included) at pressure levels above 1µbar. For the Titan-type composition particles become hotter than the gas above ∼10µbar and reach a maximum temperature difference of ∼300K at lower pressures. Although soot composition particles absorb more efficiently, their resulting temperature is lower from those assuming Titan-type composition. This feature relates to the higher atmospheric density of the upper atmosphere for the soot case, which results in a more efficient energy transfer from particles to gas, relative to the Titan-type case. Note that at all cases the temperature difference between the particles and the Figure 13. Stellar heating (dashed line) and thermal cooling rates at the atmospheric temperature (solid line) and the excess emission from soot type hot-particles (red line). The thin grey line overlapping the dashed line is the sum of thermal rates, thus demonstrates the energy balance in the calculations.

gas remains small enough (Tp/Tg<1.3) for the description of energy exchange rate to be accurate.

The strong difference between gas and particle temperatures demonstrates that the inclusion of the excess particle thermal emission is important for the correct evaluation of the energy balance in the upper atmosphere. For the soot composition, particle emission peaks between 0.1 µbar and 1 µbar but its global impact to the atmospheric energy balance is important above 0.1 mbar (Fig. 13). For Titan-type particles the effect is reduced but the picture is qualitatively the same. If particles were assumed to be in equilibrium with the gas the resulting atmospheric temperature near 1 µbar would be higher by ∼100-200 K depending on the particle refractive index assumed (see dotted lines in Fig. 12). The radiation emitted by the hot particles reflects the refractive index assumed. The high absorptivity of soots over the whole spectrum results in a black-body emission spectrum, while for Titan-type particles the strong absorption features in the near IR leave a strong signature to the emitted spectrum (Fig. 14). However, at the top of the atmosphere the particle emission is negligible compared to the gas emission from deeper layers, thus it would not be observable. As with the case of HD 209458 b, the emitted radiation dominantly originates from levels between 1 bar and 0.1 mbar. At the top of the atmosphere thermally emitted radiation dominates at wavelengths longer than the K line and reflected stellar radiation prevails at shorter wavelengths.

IMPLICATIONS

Our results demonstrate that the impact of photochemical hazes and disequilibrium chemistry on the atmospheric temperature can be significant. These effects could have further ramifications for the interpretation of transit and secondary eclipse observations. Thus, we explore these implications here.

HD 209458 b

Our previous evaluation of the transit signature of HD 209458 b suggested a rather clear atmosphere, at least for the part of the atmosphere probed by the observations (Lavvas et al. 2014). However that conclusion was based on the GCM temperature profile that included the thermochemical equilibrium abundances of TiO and VO, i.e. corresponding to a significantly warmer atmosphere than the one we calculate here. This was actually a critical factor in suggesting that the higher temperature of HD 209458 b relative to that of HD 189733 b was responsible for the lack of photochemical hazes in the former atmosphere: the photochemical composition simulations demonstrated that the main gas species considered as precursors of the photochemical hazes were significantly reduced at high temperatures. The cooler temperature profile under the disequilibrium chemical composition questions this conclusion, thus necessitates a re-evaluation of the chemical composition and photochemistry mass fluxes.

The observational constraints for the atmosphere of HD 209458 b also demand a re-evaluation of the presence of hazes. A new analysis of HST/STIS observations (Sing et al. 2016) revealed the presence of the missing potassium (Lavvas et al. 2014), but also resulted in a significantly stronger transit signature in the UV-Visible part of the spectrum relative to the previous analyses (Sing et al. 2008;Deming et al. 2013). Therefore, the updated constraints may be consistent with the presence of hazes in this atmosphere, after all (Pinhas et al. 2019).

To address these concerns we re-evaluated the disequilibrium chemical composition considering the new temperature profile evaluated in this study, all other parameters (e.g. stellar flux, eddy mixing profile) remaining the same as in our previous evaluation (Lavvas et al. 2014). We consider the complete redistribution profile for this evaluation as the disk average conditions it corresponds to, are likely more representative of the atmospheric limbs probed in transit. Clearly the cooler temperature affects the chemical composition of the upper atmosphere with most prominent modification the transition from H2 dominance to H dominance moving from 0.1 mbar in the hot temperature profile to ∼1 µbar in the cooler temperature profile (Fig. 3). This modification has further ramifications for other species that are affected by the abundance of atomic hydrogen such as H2O, N2 and CO for which their mixing ratios increase in the upper atmosphere. Using the revised thermal structure and chemical composition, the simulated transit spectrum (black line in Fig. 15) demonstrates a stronger signature in the UV-Visible slope compared to the spectrum based on the equilibrium thermal structure (red line). Most importantly, the revised observations demonstrate an even higher transit depth in the same spectral region, suggesting that a haze contribution could be possible.

The main photochemical haze precursors such as CH4, HCN and C2H2 demonstrate higher abundances at the disequilibrium temperature. The most significant increase appears in the profile of HCN that attains a mixing ratio of ∼10 -6 in the upper atmosphere, while the profiles of CH4 and C2H2 also increase but their mixing ratios remain in the order of 10 -9 . The corresponding mass fluxes from the photolysis of these species reflects the changes in their abun- 2016). The red solid line presents the simulated spectrum assuming the equilibrium composition temperature profile (i.e. utilises the chemistry and temperature profiles presented by dashed lines in Figs. 3 and7), the black solid line assumes the disequilibrium composition temperature profile (i.e utilises the chemistry profiles presented by solid lines in Fig. 3 and black line temperature profile in Fig. 7), and the green solid line assumes the disequilibrium + photochemical haze profile discussed in the text. For the eclipse spectrum observations are from Deming et al. (2005) 2016) (orange crosses, see also inset). Dashed lines present theoretical planet/star flux ratios for planetary black body emissions at the presented temperatures. The solid line is our simulated spectrum for the disequilibrium chemical composition assuming energy deposition on the day side, while the dash-dotted line spectrum assumes complete redistribution.

dances (Table 1). The mass flux generated by the photolysis of HCN at pressure levels above 10µbar is ∼60× larger from our previous evaluation based on the equilibrium temperature profile (Lavvas et al. 2014), providing a mass flux of 4.8×10 -13 gcm -2 s -1 . Similarly, photolysis of CH4 and C2H2 generate mass fluxes that are ∼10 6 × and ∼10 3 × larger than previous estimates, respectively. Previously, we estimated the mass flux for soot particles from the photolysis of these three species considering a moderate yield of 1 per cent. Such an estimate for the new conditions investigated results in a mass flux of 5.2×10 -15 gcm -2 s -1 , which is ∼70× larger than the previous evaluation. Our previous conclusion regarding the role of high temperature as a modulator of the photochemical haze production remains valid: the revised mass flux in the atmosphere of HD 209458 b is still significantly smaller than that in the cooler atmosphere of HD 189733 b. Nevertheless, the updated mass flux magnitude may suggest that photochemical hazes could be present in the atmosphere of HD 209458 b.

We therefore considered how a photochemical haze mass flux with the above magnitude would affect the transit spectrum. Our preliminary evaluation of the particle size distribution suggests that the effect would be negligible. However, for this evaluation we assumed the eddy mixing profile derived from the atmospheric circulation simulations that considered the equilibrium chemical composition with TiO/VO [START_REF] Showman | [END_REF], thus is not representative of the conditions investigated here. Assuming a weaker eddy mixing, as also suggested by subsequent evaluations with GCM models (Parmentier et al. 2013;Zhang & Showman 2018), allows for a more efficient particle growth (Lavvas & Koskinen 2017). A 100× weaker eddy profile results in the green transit spectrum in Fig. 15 that agrees better to the observed UV-Visible transit depth. However, a higher formation yield for photochemical hazes is also possible (Lavvas & Koskinen 2017). Moreover, photolysis of CO and CO2 may further con- tribute to the mass flux of photochemical hazes, as recent laboratory experiments suggest, albeit of a smaller yield relative to the pathways based on methane (He et al. 2018).

On the other hand, the revised observations indicate the presence of other features in the transit spectrum. Particularly features at the Na wings and at shorter wavelengths are reminiscent of gaseous opacities that would suggest there are missing components in the gaseous abundances considered. In the past we evaluated the contributions of heavier atoms such as Al, Ca, Mg, Fe and Si in the transit spectrum (see fig. 10 in Lavvas et al. 2014), assuming that these are not lost due to condensation, with the results suggesting that they could contribute significantly to the opacity. However, such contributions are expected to be far more significant in ultra hot exoplanets (Lothringer et al. 2020). On the contrary, the low temperature conditions we calculate here suggest that heavy metals should remain condensed in the atmosphere of HD 209458 b and should not contribute to the atmospheric opacity. This conclusion is further supported by a recent reanalysis of the NUV spectrum of HD 209458 b that did not detect the presence of Mg/Mg + in the upper atmosphere (Cubillos et al. 2020), although the detection of Fe beyond the planetary roche lobe is a concern. Moreover, the observed transit shape of H2O at 1.4 µm appears also affected by the haze, which is not reproduced by the few nm radius particles resulting from our microphysics simulations. This characteristic further indicates that condensate clouds could also be affecting the transit depth at near IR wavelengths. Based on our transit simulation the pressures probed at wavelengths λ<2µm are larger than 1 mbar, which makes the possibility of photochemical haze interaction with cloud formation very likely (see condensation curves for silicates in Fig. 7). All these uncertainties demonstrate the complexity of the problem, as well as, the need for further investigation of the interaction mechanisms among different types of particles and their feedback to the gas background.

We need to keep in mind that there are further complications arising from the 3D aspect of the atmosphere. These complications are better demonstrated by the comparison with the secondary eclipse observations. Our simulations assuming complete redistribution fall short of the observations, but assuming energy deposition only on the day side provides a spectrum that is in better agreement (Fig. 15). The corresponding temperature profile (green line in Fig. 7) is ∼200-300 K hotter than the nominal case of complete redistribution, and consistent with the Line et al. (2016) temperature retrievals based on secondary eclipse observations (see inset in Fig. 7). However, even with the 2πredistribution the simulated eclipse spectrum does demonstrate a lower emission near 10 µm relative to the observations. Our calculations suggest that photons in this spectral region mainly originate from pressure levels between 1 mbar and 1 bar, i.e. where also our transit simulations suggest that possible hazes should reside. The presence of photochemical hazes could locally heat the atmosphere and increase the thermally emitted radiation, bringing the simulated emission spectrum in better agreement with the observations, as demonstrated below for the case of HD 189733 b. The same effect could also result from the presence of clouds. Thus, the secondary eclipse observations may further support the presence of aerosols in the atmosphere of HD 209458 b.

HD 189733 b

For HD 189733 b, we note that the transit at short wavelengths (λ 1µm) probes the upper atmosphere at pressures between ∼1 mbar and ∼1 µbar (Lavvas & Koskinen 2017). As demonstrated above, the inclusion of photochemical haze opacity affects this part of the atmosphere and the involved changes to the atmospheric scale height could modify the simulated transit spectra. Therefore, we evaluated the difference in the transit spectrum of HD 189733 b assuming the temperature profiles with and without haze (Fig. 16). The results reveal that the hotter upper atmosphere with the inclusion of soot haze is significantly more extended yielding a stronger signature to the transit depth by increasing the UV-Visible slope (red line). The impact of hazes on the ther-mal structure, may therefore reduce the requirement for a strong eddy mixing to keep the particles aloft at high altitudes (Ohno & Kawashima 2020) and potentially allow for a convergence between the eddy mixing profile necessary to explain observations and those estimated by GCM models (Parmentier et al. 2013;Zhang & Showman 2018).

On the other hand, considering the Titan-type composition demonstrates the requirement for a re-evaluation of the appropriate mass flux and size distribution; the resulting transit depth assuming the same particle size distribution as for soot particles falls short of the observations (blue line) due to the lower opacity demonstrated by the Titan-type composition. A preliminary evaluation suggests that a factor of ∼10 higher opacity would be required to match the transit spectrum at visible and UV wavelengths assuming a Titantype refractive index. Such an increase in opacity would also bring the temperature closer to the profile derived for the soot composition particles. Thus, the important conclusion from this exercise is that the correct interpretation of transit spectra for hazy exoplanets requires proper consideration of the haze impact on the thermal structure.

Our temperature profiles presented above furthermore demonstrate that inclusion of photochemical hazes also affects the photosphere therefore could modify the secondary eclipse signature. Our simulated eclipse spectrum for HD 189733 b is consistent with the observations for all simulated cases assuming energy depositions on the day side (Fig. 16). Nevertheless, the temperature changes imposed by the inclusion of haze opacity reduces the difference between model and observation (particularly for λ>10µm), with the impact being more pronounced for the soot composition haze that has a stronger influence to the part of the atmosphere affecting the thermally emitted radiation (see Fig. 12). Although this trend is clear, a definite conclusion regarding the role of photochemical hazes on the eclipse spectrum requires further investigation; photochemical hazes could act as nucleation sites for the formation of clouds in the photosphere, which will modify both the size distribution and optical properties of the resulting particles. Therefore a complete evaluation of the impact of hazes on the atmospheric structure and the interpretation of eclipse observations requires a coupled photochemical haze-cloud simulation. With the advent of more accurate observations in the forthcoming years through JWST, ARIEL, and other observatories, it will be possible to set more accurate constraints on the contribution of clouds and hazes.

The temperature increase in the upper atmosphere of HD 189733 b imposed by the photochemical haze opacity will have an impact on both the atmospheric chemistry and the formation and evolution of the haze particles, which we have not yet taken into account. Therefore, our work here presents only part of the implications from the inclusion of such particles. Considering only gas phase chemistry Drummond et al. (2016) presented self-consistent simulations of disequilibrium chemical composition and temperature and found that the coupling reduced the departure from equilibrium chemistry composition. However, the haze opacity in our simulations is orders of magnitude higher than the gaseous opacity at visible wavelengths (Fig. 10). This dramatic difference in opacity sources demonstrates why the potential changes in the chemical composition from the self- 2016). Lines present simulations assuming the thermochemical equilibrium temperature profile (black) and the disequilibrium temperature profiles evaluated here under the soot (red) and the Titan-type (blue) composition assumptions (for all cases the chemical composition is described by the dashed profiles in Fig. 8). For the secondary eclipse spectrum observations are from Charbonneau et al. (2008) (red crosses), Knutson et al. (2007Knutson et al. ( , 2009) (green asterisks), Grillmair et al. (2008) (blue diamonds), Swain et al. (2009) (cyan triangles). Dashed lines are theoretical planet/star flux ratios for planetary black body emissions at the presented temperatures. Solid lines present the simulated spectrum for the cases presented earlier (black: no haze, red: soot haze, blue: Titan-type haze) evaluated assuming energy deposition on the day side.

consistent description are not likely to modify our conclusion on the impact of hazes on the atmospheric thermal structure.

As a demonstration of the anticipated ramifications, our disequilibrium chemical composition simulations using the disequilibrium temperature profile with soot haze opacity reveal non-negligible changes in the profiles of the main gaseous species (Fig. 8). These changes reflect modifications of the chemical kinetic rates and diffusion properties of the species through their temperature dependence. Note that temperature also affects the thermochemical equilibrium that controls the lower boundary conditions for the main species. Thus, results in further modifications of the main abundances depending on the elemental composition and Tint value assumed [START_REF] Fortney | Beyond Equilibrium Temperature: How the Atmosphere/Interior Connection Affects the Onset of Methane, Ammonia, and Clouds in Warm Transiting Giant Planets[END_REF]. The disequilibrium temperature profile results in higher abundances for H2O, CH4 and NH3 in the photosphere region, while the H2 to H transition starts at lower pressures near 0.1 mbar. CH4 and C2H2 have smaller abundances in the upper atmosphere, compared to the profiles under the equilibrium temperature profile, while HCN demonstrates a strong increase in the whole atmosphere. The mass fluxes from the photolysis of the main species are indeed different from the previous evaluation based on the equilibrium composition temperature profile (Table 1) suggesting a higher production of photochemical hazes.

However, even without the changes in the background gas composition, the higher temperature in the upper atmosphere will independently affect the microphysical evolution of the hazes. As a demonstration of the impact of temperature on the evolution of the haze particles we compared their microphysical evolution under the clear atmosphere and soot-composition haze scenarios (Fig. 9). The results show that the higher temperature drives a more efficient collision rate of the particles that yields a faster particle growth in the upper atmosphere relative to the case of the clear atmosphere temperature profile. Note that the temperature increase is not sufficient to sublimate the particles, as seen in the lower atmosphere (Lavvas & Koskinen 2017), therefore does not pose a problem for the presence of the particles at the simulated conditions.

Finally, we can note from the simulated transit spectra (Fig. 16) that a refractive index similar to the soot, but with a lower absorptivity towards the near IR would be more consistent with the observations. However, a re-evaluation of the microphysics taking into account the impact on the thermal structure is required before any definite conclusions. The revised chemical composition and haze distribution will further affect the thermal structure necessitating a new temperature evaluation as well as a new evaluation of the optimum haze mass fluxes to match the transit spectrum. Therefore it is clear that a self-consistent description of thermal structure, atmospheric chemistry and haze microphysics is necessary for a proper evaluation of the feedback among each domain, and the interpretation of current and future observations.

CONCLUSIONS

We describe a detailed radiative-convective model for the investigation of the atmospheric thermal structure of hot exoplanets and investigate the impact of photochemical hazes. We find that the inclusion of such particles result in major heating of the upper and cooling of the lower atmosphere. The magnitude of this effect depends on the refractive index assumed for the optical properties of the particles, but our simulations demonstrate that the effect is important for the highly absorbing soot composition, as well as, for the relatively more scattering Titan-type composition.

We also considered the role of disequilibrium chemical composition relative to the usual assumption of thermochemical equilibrium (with or without rainout) in previous thermal structure evaluations. We find that sulphur containing species, such as SH, S2 and S3 provide significant opacity in the middle atmosphere and lead to local heating near 1 mbar. Similarly OH, CH, NH, and CN radicals produced by the photochemistry affect the thermal structure near 1µbar.

Furthermore we show that the modifications on the thermal structure from photochemical gases and hazes can have important ramifications for the interpretation of transit observations. Our study for the hazy HD 189733 b shows that the hotter upper atmosphere resulting from the inclusion of photochemical haze opacity imposes an expansion of the atmosphere, thus a steeper transit signature in the UV-Visible part of the spectrum. In addition, the temperature changes in the photosphere also affect the secondary eclipse spectrum. For HD 209458 b we find that the cooler atmospheric thermal structure under disequilibrium chemical composition provides a transit spectrum that is no longer consistent with the latest observational constraints. This conclusions suggests that some haze opacity could be present in this atmosphere, at pressures below 1 mbar, which could be a result of both photochemical hazes and condensates. However further investigation is required with a selfconsistent coupled description.

The above conclusions for the role of hazes are particularly important for general circulation models that investigate the horizontal distribution of particles by the winds and are able to evaluate the whole phase curve variation of the planetary emission. Currently most studies that involve particles follow a tracer approach that reveals how fixed size particles are distributed in the atmosphere but do not include the feedback of particle absorption to the thermal structure and the winds field. Our study motivates the inclusion of this feedback in future GCM simulations.

Figure 1 .

 1 Figure 1. Comparison of absorption cross sections of radicals OH, CH, NH, and SH at 1000 K. All cross sections have been smoothed for clarity.

Figure 2 .

 2 Figure 2. Comparison of refractive indices of soot type (red) and Titan-type (blue) compositions. Soots demonstrate a far more absorbing behaviour than the Titan-type composition, manifested over the whole spectrum. On the contrary the Titan-type material absorbs mainly in the short visible and UV wavelengths and at narrow bands in the near IR.

Figure 3 .

 3 Figure 3. Dis-equilibrium chemical composition for the atmosphere of HD 209458 b assuming the Showman et al. (2009) dayaverage GCM temperature profile (dashed lines) and the disequilibrium composition temperature profile calculated in this work for complete (4π) redistribution (solid lines).

Figure 4 .

 4 Figure4. Gaseous opacity contributions across wavelength for stellar energy deposition calculations and at different pressures in the atmosphere of HD 209458 b. Different colours present the absorption optical depth contribution by different gases as labeled in each panel. The solid black line corresponds to the total absorption optical depth based on the dis-equilibrium composition, while the dotted line is the total opacity when the contributions of TiO and VO (evaluated at thermodynamic equilibrium) are also taken into account. The dash-dotted lines present the total scattering opacity (molecular and electrons).

Figure 5 .Figure 6 .

 56 Figure 5. Stellar flux penetration in the atmosphere of HD 209458 b at different wavelengths.

Figure 7 .

 7 Figure 7. Temperature profiles for HD 209458 b. Left: Temperature profile under nominal conditions of disequilibrium chemistry (thick solid black line) compared to other 1D simulations assuming thermochemical equilibrium composition (dash-dotted,Iro et al. 2005), (dashtriple-dottedFortney et al. 2006). The purple line presents our simulated temperature profile without the contribution of photochemical radicals in the simulations (see text). All these simulations assume a 4π redistribution of incoming stellar energy and T int =100 K. The green line presents our nominal thermal structure assuming energy deposition on the day side only (2π redistribution). The black dashed line presents the thermospheric temperature profile from escape models(Koskinen et al. 2013), and the grey dotted lines the equilibrium saturation curves for MgSiO 3 and Mg 2 SiO 4(Visscher et al. 2010). For all simulations, convection does not occur in the simulated pressure range. The inset compares our nominal 2π (black) and 4π (green) temperature profiles with the profiles retrieved byLine et al. (2016) based on secondary eclipse observations. The red and blue profiles with associated uncertainties are based on different parameterizations for the temperature structure used in their analysis. Right: Temperature profiles including TiO/VO at equilibrium density (see text). The thick lines present our simulations for 4π (blue) and 2π (red) stellar energy redistribution and the dashed line corresponds to the day-average temperature profile from 3D simulations[START_REF] Showman | [END_REF], which is the temperature profile used for the evaluation of the disequilibrium chemistry.

Figure 8 .

 8 Figure 8. Dis-equilibrium chemical composition for the atmosphere of HD 189733 b assuming the Showman et al. (2009) dayaverage GCM temperature profile (dashed lines) and the disequilibrium composition temperature profile with soot haze calculated in this work for complete (4π) redistribution (solid lines).

Figure 10 .

 10 Figure10. Gas and photochemical haze opacity contributions across wavelength for stellar energy deposition calculations at different pressures in the atmosphere of HD 189733 b. Different colours present the absorption optical depth contribution by different gases based on the dis-equilibrium chemical composition. The dashed black line presents the absorption optical depth by the haze particles and the solid black line the total absorption optical depth (gas and haze). The dash-dotted lines present the total scattering optical depth.

Figure 11 .

 11 Figure 11. Stellar flux penetration in the atmosphere of HD 189733 b at different wavelengths.

Figure 14 .

 14 Figure14. Planetary thermal emission from gases (black) and particles (blue) and stellar reflected radiation (red) at the top of the simulated atmosphere, assuming 4π-redistribution. Dashed lines present black body emission at 1000 K, 1250 K and 1500 K. The top panel corresponds to soot composition particles and the bottom to Titan-type composition.

Figure 15 .

 15 Figure 15. Observed (symbols) and simulated (lines) transit (left) and secondary eclipse (right) spectra of HD 209458 b. Transit observations are fromSing et al. (2016). The red solid line presents the simulated spectrum assuming the equilibrium composition temperature profile (i.e. utilises the chemistry and temperature profiles presented by dashed lines in Figs.3 and 7), the black solid line assumes the disequilibrium composition temperature profile (i.e utilises the chemistry profiles presented by solid lines in Fig.3and black line temperature profile in Fig.7), and the green solid line assumes the disequilibrium + photochemical haze profile discussed in the text. For the eclipse spectrum observations are from Deming et al. (2005) (red cross), Swain et al. (2008) (cyan triangles), Swain et al. (2009) (navy asterisks), Crossfield et al. (2012) (green asterisk), Diamond-Lowe et al. (2014) (blue diamonds), Zellem et al. (2014) (purple square), and Line et al. (2016) (orange crosses, see also inset). Dashed lines present theoretical planet/star flux ratios for planetary black body emissions at the presented temperatures. The solid line is our simulated spectrum for the disequilibrium chemical composition assuming energy deposition on the day side, while the dash-dotted line spectrum assumes complete redistribution.

  Figure 15. Observed (symbols) and simulated (lines) transit (left) and secondary eclipse (right) spectra of HD 209458 b. Transit observations are fromSing et al. (2016). The red solid line presents the simulated spectrum assuming the equilibrium composition temperature profile (i.e. utilises the chemistry and temperature profiles presented by dashed lines in Figs.3 and 7), the black solid line assumes the disequilibrium composition temperature profile (i.e utilises the chemistry profiles presented by solid lines in Fig.3and black line temperature profile in Fig.7), and the green solid line assumes the disequilibrium + photochemical haze profile discussed in the text. For the eclipse spectrum observations are from Deming et al. (2005) (red cross), Swain et al. (2008) (cyan triangles), Swain et al. (2009) (navy asterisks), Crossfield et al. (2012) (green asterisk), Diamond-Lowe et al. (2014) (blue diamonds), Zellem et al. (2014) (purple square), and Line et al. (2016) (orange crosses, see also inset). Dashed lines present theoretical planet/star flux ratios for planetary black body emissions at the presented temperatures. The solid line is our simulated spectrum for the disequilibrium chemical composition assuming energy deposition on the day side, while the dash-dotted line spectrum assumes complete redistribution.

Figure 16 .

 16 Figure16. Observed (symbols) and simulated (lines) transit (left) and secondary eclipse (right) spectra of HD 189733 b. For the transit spectrum observations are fromSing et al. (2016). Lines present simulations assuming the thermochemical equilibrium temperature profile (black) and the disequilibrium temperature profiles evaluated here under the soot (red) and the Titan-type (blue) composition assumptions (for all cases the chemical composition is described by the dashed profiles in Fig.8). For the secondary eclipse spectrum observations are fromCharbonneau et al. (2008) (red crosses),Knutson et al. (2007Knutson et al. ( , 2009) (green asterisks),Grillmair et al. (2008) (blue diamonds),Swain et al. (2009) (cyan triangles). Dashed lines are theoretical planet/star flux ratios for planetary black body emissions at the presented temperatures. Solid lines present the simulated spectrum for the cases presented earlier (black: no haze, red: soot haze, blue: Titan-type haze) evaluated assuming energy deposition on the day side.

Table 1 .

 1 Photolysis mass fluxes (in gcm -2 s -1 ) of major species in the upper atmosphere (p<10µbar) of HD 209458 b and HD 189733 b under the two temperature profiles of equilibrium (E) and disequilibrium (DE) composition of each case. The soot mass flux is estimated from the sum of CH 4 , C 2 H 2 and HCN mass fluxes assuming a moderate 1% efficiency(Lavvas & Koskinen 2017).

		HD 209458 b	HD 189733 b
	Species	E p-T	DE p-T	E p-T	DE p-T
	CH 4	1.1(-22)	6.2(-16) 7.1(-16)	2.8(-16)
	HCN	7.8(-15)	4.8(-13) 2.9(-11)	2.1(-10)
	C 2 H 2	2.8(-17)	4.1(-14) 1.8(-11)	8.6(-12)
	CO	4.8(-11)	5.9(-12) 4.8(-12)	5.4(-13)
	N 2	4.5(-12)	5.1(-13) 7.0(-13)	2.7(-12)
	NH 3	3.5(-16)	2.9(-13) 1.6(-15)	3.0(-15)
	H 2 S	1.8(-11)	1.5(-11) 5.9(-13)	6.3(-13)
	S 3	1.8(-19)	2.0(-12) 1.4(-13)	3.2(-14)
	Soot (1%)	7.8(-17)	5.2(-15) 4.7(-13)	2.2(-12)

© 2020 The Authors

MNRAS 000, 1-18 (2020)

ACKNOWLEDGEMENTS

We acknowledge support by the Programme National de Planétologie through the project TISSAGE. We thank R.V. Yelle for providing the S2 linelist and T.T. Koskinen for comments on the manuscript. We also thank our reviewer for constructive comments and suggestions that resulted in an improved manuscript.

DATA AVAILABILITY

The data underlying this article will be shared on reasonable request to the corresponding author.