© ESO, 2014

1. Introduction

Strongly irradiated by their close host star, hot Jupiters reside in extreme environments and represent a class of planets without analogue in our solar system. This type of exoplanets, the first to be discovered around main-sequence stars, remains the best available to study through observations and challenge a variety of models in the area of planetary science (see reviews on the subject by Baraffe et al. 2010; Seager & Deming 2010; Burrows & Orton 2011; Showman et al. 2011). In recent years, multiwavelength observations of transiting hot Jupiters have allowed scientists to put constraints on the physical and chemical state of their atmospheres. Among these hot Jupiters, the best characterized are probably HD 209458b and HD 189733b, which belong to some of the brightest and closest transiting systems, and for which primary transit, secondary eclipse, and phase curve measurements have been used to probe, though often with controversial interpretations, various characteristics of their atmospheres, such as the thermal structure (Deming et al. 2005, 2006; Knutson et al. 2008; Charbonneau et al. 2008), winds and day-night heat redistribution (Knutson et al. 2007, 2012; Cowan et al. 2007; Snellen et al. 2010), and mixing ratios of some of the main molecular constituents (Tinetti et al. 2007; Swain et al. 2008, 2009a, 2009b, 2010; Grillmair et al. 2008; Sing et al. 2009; Désert et al. 2009; Beaulieu et al. 2010; Gibson et al. 2011; Waldmann et al. 2012; Lee et al. 2012; Rodler et al. 2013; de Kok et al. 2013).

The ability of observations at infrared and visible wavelengths to characterize the physical and chemical state of exoplanet atmospheres has motivated the development of various types of theoretical models. On the one hand, there are those aiming at investigating the physical structure of hot-Jupiter atmospheres, either one-dimensional radiative models (Iro et al. 2005; Fortney et al. 2008; Parmentier et al. 2014) or three-dimensional general circulation models (Cho et al. 2008; Showman et al. 2009; Heng et al. 2011a,b; Dobbs-Dixon et al. 2012; Rauscher & Menou 2013; Parmentier et al. 2013). These models have shown how fascinating the climates of hot Jupiters are, with atmospheric temperatures usually in excess of 1000 K, and helped to understand some global observed trends. Some hot Jupiters are found to display a strong thermal inversion in the dayside while others do not (e.g. Fortney et al. 2008). Strong winds with velocities of a few km s^-1 develop and transport the heat from the dayside to the nightside, reducing the temperature contrast between the two hemispheres. The circulation pattern in these planets is characterized by an equatorial superrotating eastward jet. On the other hand, the chemical composition of hot Jupiters has been investigated by one-dimensional models, which currently account for thermochemical kinetics, vertical mixing, and photochemistry (Zahnle et al. 2009; Line et al. 2010, 2011; Moses et al. 2011, 2013; Kopparapu et al. 2012; Venot et al. 2012, 2014; Agúndez et al. 2014). These models have revealed the existence of three different chemical regimes in the vertical direction. A first one at the bottom of the atmosphere, where the high temperatures and pressures ensure a chemical equilibrium composition. A second one located above this, where the transport of material between deep regions and higher layers occurs faster than chemical kinetics so that abundances are quenched at the chemical equilibrium values of the quench level. And a third one located in the upper atmosphere, where the exposure to ultraviolet (UV) stellar radiation drives photochemistry. The exact boundaries between these three zones depend on the physical conditions of the atmosphere and on each species.

In addition to the retrieval of average atmospheric quantities from observations, there is a growing interest in the physical and chemical differences that may exist between different longitudes and latitudes in hot-Jupiter atmospheres, and in the possibility of probing these gradients through observations. Indeed, important temperature contrasts between different planetary sides of hot Jupiters, noticeably between day and night sides, have been predicted (Showman & Guillot 2002), observed for a dozen hot Jupiters (see Knutson et al. 2007 for the first one), qualitatively understood (Cowan & Agol 2011; Perez-Becker & Showman 2013), and confirmed by three-dimensional general circulation models (e.g. Perna et al. 2012). These temperature gradients, together with the fact that photochemistry switches on and off in the day and night sides, are at the origin of a potential chemical differentiation in the atmosphere along the horizontal dimension, especially along longitude. On the other hand, strong eastward jets with speeds of a few km s^-1 are believed to dominate the atmospheric circulation in the equatorial regions, as predicted by Showman & Guillot (2002), theorized in Showman & Polvani (2011), potentially observed by Snellen et al. (2010), and confirmed by almost all general circulation models of hot Jupiters. These strong horizontal winds are an important potential source of homogenization of the chemical composition between locations with different temperatures and UV illumination. The existence of winds and horizontal gradients in the temperature and chemical composition of hot-Jupiter atmospheres has mainly been considered from a theoretical point of view, although some of these effects can be studied through phase curve observations (Fortney et al. 2006; Cowan & Agol 2008; Majeau et al. 2012; de Wit et al. 2012), monitoring of the transit ingress and egress (Fortney et al. 2010), and Doppler shifts of spectral lines during the primary transit (Snellen et al. 2010; Miller-Ricci Kempton & Rauscher 2012; Showman et al. 2013).

The existence of horizontal chemical gradients has been addressed in the frame of a series of one-dimensional models in the vertical direction at different longitudes (e.g. Moses et al. 2011). An attempt to understand the interplay between circulation dynamics and chemistry was undertaken by Cooper & Showman (2006), who coupled a three-dimensional general circulation model of HD 209458b to a simple chemical kinetics scheme dealing with the interconversion between CO and CH₄. These authors found that, even in the presence of strong temperature gradients, the mixing ratios of CO and CH₄ are homogenized throughout the planet’s atmosphere in the 1 bar to 1 mbar pressure regime. In our team, we have recently adopted a different approach in which we coupled a robust chemical kinetics scheme to a simplified dynamical model of HD 209458b’s atmosphere (Agúndez et al. 2012). In this approach the atmosphere was assumed to rotate as a solid body, mimicking a uniform zonal wind, while vertical mixing and photochemistry were neglected. We found that the zonal wind acts as a powerful disequilibrium process that tends to homogenize the chemical composition, bringing molecular abundances at the limb and nightside regions close to chemical equilibrium values characteristic of the dayside. Here we present an improved model that simultaneously takes into account thermochemical kinetics, photochemistry, vertical mixing, and horizontal transport in the form of a uniform zonal wind. We apply our model to study the interplay between atmospheric dynamics and chemical processes, and the distribution of the main atmospheric constituents in the atmosphere of the hot Jupiters HD 209458b and HD 189733b.

2. Model

We modeled the atmospheres of HD 209458b and HD 189733b, for which we adopted the parameters derived by Southworth (2010). For the system of HD 209458 we took a stellar radius of 1.162 R_⊙, a planetary radius and mass of 1.38 R_J and 0.714 M_J (where R_J and M_J stand for Jupiter radius and mass), and an orbital distance of 0.04747 au. For the system of HD 189733 the adopted parameters are a stellar radius of 0.752 R_⊙, a planetary radius and mass of 1.151 R_J and 1.150 M_J, and a planet-to-star distance of 0.03142 au.

The atmosphere model is based on some of the outcomes of three-dimensional general circulation models (GCMs) developed for HD 209458b and HD 189733b (Showman et al. 2009; Parmentier et al. 2013, and in prep.), which indicate that circulation dynamics is dominated by a broad eastward equatorial jet. On the assumption that the eastward jet dominates the circulation pattern, it seems well justified to model the atmosphere as a vertical column that rotates along the equator, which mimicks a uniform zonal wind. The main shortcoming of this approach is that it reduces the whole circulation dynamics to a uniform zonal wind, although it has the clear advantage over more traditional one-dimensional models in the vertical direction of simultaneously taking into account the mixing and transport of material in the vertical and horizontal directions.

2.1. Pseudo two-dimensional chemical model

In one-dimensional models of planetary atmospheres, the distribution of each species in the vertical direction is governed by the coupled continuity-transport equation (1)where f_i is the mixing ratio of species i, t the time, n the total number density of particles, r the radial distance to the center of the planet, P_i and L_i the rates of production and loss, respectively, of species i, and φ_i the vertical transport flux of particles of species i (positive upward and negative downward). The first two terms on the right side of Eq. (1) account for the formation and destruction of species i by chemical and photochemical processes, while the third term accounts for the vertical transport in a spherical atmosphere. In this way, thermochemical kinetics, photochemistry, and vertical mixing can be taken into account through Eq. (1). The transport flux can be described by eddy and molecular diffusion as (2)where z is the altitude in the atmosphere with respect to some reference level (typically set at a pressure of 1 bar), T is the gas kinetic temperature, K_zz is the eddy diffusion coefficient, D_i is the coefficient of molecular diffusion of species i, H_i is the scale height of species i, H₀ is the mean scale height of the atmosphere, and α_i is the thermal diffusion factor of species i. More details on Eqs. (1) and (2) can be found, for instance, in Bauer (1973) and Yung & DeMore (1999). The coefficient of molecular diffusion D_i is estimated from the kinetic theory of gases (see Reid et al. 1988), while the factor of thermal diffusion α_i is set to −0.25 for the light species H, H₂, and He (Bauer 1973), and to 0 for the rest of species. The eddy diffusion coefficient K_zz is a rather empirical formalism to take into account advective and turbulent mixing processes in the vertical direction, and is discussed in more detail in Sect. 2.2.

To compute the abundances of the different species as a function of altitude, the atmosphere is divided into a certain number of layers and the continuous variables in Eqs. (1) and (2) are discretized as a function of altitude. After the discretization, Eq. (1) reads (3)where the superscript j refers to the jth layer, while j + 1/2 and j−1/2 refer to its upper and lower boundaries, respectively, so that layers are ordered from bottom to top. The transport fluxes of species i at the upper and lower boundaries of layer j, and , are then given by (4)where the variables evaluated at j + 1/2 and j−1/2 boundaries are approximated as the arithmetic mean of the values at layers j and j + 1 and at layers j−1 and j, respectively. We assume that there is neither gain nor loss of material in the atmosphere, and thus the transport fluxes at the bottom and top boundaries of the atmosphere are set to zero.

The rates of production and loss of each species in Eqs. (1) and (3) are given by chemical and photochemical processes. Thermochemical kinetics is taken into account with a chemical network, which consists of 104 neutral species composed of C, H, N, and O linked by 1918 chemical reactions, that has been validated in the area of combustion chemistry by numerous experiments over the 300–2500 K temperature range and the 0.01–100 bar pressure regime, and has been found suitable to model the atmospheres of hot Jupiters. Most reactions are reversed with their rate constants fulfilling detailed balance to ensure that, in the absence of disequilibrium processes such as photochemistry or mixing, thermochemical equilibrium is achieved at sufficiently long times. The reaction scheme is described in Venot et al. (2012), with some minor modifications given in Agúndez et al. (2012). As photochemical processes we consider photodissociations, whose rates depend on the incident UV flux and the relevant cross sections. The incident UV flux is calculated by solving the radiative transfer in the vertical direction for a given zenith angle, where the spherical geometry of layers is taken into account when computing the path length along each of them. Absorption and Rayleigh scattering, the latter being treated through a two-ray iterative algorithm (Isaksen et al. 1977), both contribute to the attenuation of UV light throughout the atmosphere. Absorption and photodissociation cross-sections are described in detail in Venot et al. (2012). Rayleigh-scattering cross-sections are calculated for the most abundant species from their polarizability (see e.g. Tarafdar & Vardya 1969). As UV spectrum for the host star HD 209458, we adopt the spectrum of the Sun (mean between minimum and maximum activity from Thuillier et al. 2004) below 168 nm and a Kurucz synthetic spectrum¹ at longer wavelengths. For HD 189733b, below 335 nm we adopt a UV spectrum of ϵ Eridani based on the CAB X-exoplanets archive (Sanz-Forcada et al. 2011) and observations with FUSE and HST (see details in Venot et al. 2012), and a Kurucz synthetic spectrum² above 335 nm.

In one-dimensional vertical models of planetary atmospheres, the system of differential equations given by Eq. (3), with as many equations as the number of layers times the number of species, is integrated as a function of time, starting from some initial composition, usually given by thermochemical equilibrium, until a steady state is reached. During the evolution, the physical conditions of the vertical atmosphere column remain static. In the pseudo two-dimensional approach adopted here, we consider that the vertical atmosphere column rotates around the planet’s equator, and thus the system of differential equations is integrated as a function of time with physical conditions varying with time, according to the periodic changes experienced during this travel. A vertical atmosphere column rotating around the equator mimics a uniform zonal wind, which is an idealization of the equatorial superrotating jet structure found by three-dimensional GCMs for hot-Jupiter atmospheres. This approach may be seen as a pseudo two-dimensional model in which the second dimension, which corresponds to the longitude (the first one being the altitude), is in fact treated as a time dependence in the frame of an atmosphere column rotating around the equator.

To build the pseudo two-dimensional chemical models of HD 209458b and HD 189733b the vertical atmosphere column is divided into 100–200 layers spanning the pressure range 500–10^-8 bar. The evolution of the vertical atmosphere column starts at the substellar point with an initial composition given by either thermochemical equilibrium or a one-dimensional vertical model, the latter usually resulting in shorter integration times before a periodic state is reached. The convenience of starting with the composition of the hottest substellar regions is discussed in Agúndez et al. (2012). Thermochemical equilibrium calculations were carried out using a code that minimizes the Gibbs energy based on the algorithm of Gordon & McBride (1994) and the thermochemical data described in Venot et al. (2012) for the 102 species included. A solar elemental composition (Asplund et al. 2009) was adopted for the atmospheres of both HD 209458b and HD 189733b. The planetary sphere was then discretized into a certain number of longitudes (typically 100) and the system of differential equations given by Eq. (3) was integrated as the atmosphere column moves from one longitude to the next, at a constant angular velocity. To speed up the numerical calculations, the physical variables that vary with longitude (in our case these are the vertical structures of temperature and incident UV flux) were discretized as a function of longitude, that is, they were assumed to remain constant within each discretized longitude interval. As long as there are important longitudinal temperature gradients, the atmospheric scale height also varies with longitude, so that the atmosphere expands or shrinks depending on whether it gets warmer or cooler. To incorporate this effect, which may have important consequences for transit spectra, the vertical atmosphere column was enlarged or compressed (the radius at the base of the atmosphere remaining fixed) to fulfill hydrostatic equilibrium at any longitude. The variation of the incident UV flux with longitude was taken into account through the zenith angle. At the limbs we considered a zenith angle slightly different from a right angle because of the finite apparent size of the star and because of atmospheric refraction, for which we adopted a refraction angle of half a degree as in the case of visible light at Earth.

The nonlinear system of first-order ordinary differential equations given by Eq. (3) was integrated as a function of time using a backward differentiation formula implicit method for stiff problems implemented in the Fortran solver DLSODES within the ODEPACK package³ (Hindmarsh 1983; Radhakrishnan & Hindmarsh 1993). The evolution of the vertical atmosphere column was followed during several rotation cycles until the abundances of the main atmospheric constituents achieved a periodic behavior, which for HD 209458b and HD 189733b, occurs after some tens or hundreds of rotation periods.

2.2. Atmospheric dynamics and temperature (GCMs)

The pseudo two-dimensional chemical model needs some key input data related to the zonal wind speed, thermal structure, and strength of vertical mixing. These data are calculated with the three-dimensional general circulation model SPARC/MITgcm developed by Showman et al. (2009), in which dynamics and radiative transfer are coupled. The data used here for HD 209458b are based on the simulations by Parmentier et al. (2013), while those for HD 189733b are based on calculations by Parmentier et al. (in prep.), both of which cover a pressure range from about 200 bar to 2 μbar. These GCM simulations provide a wealth of detailed information regarding the physical structure of the atmosphere, although they remain limited with respect to the chemical structure as long as the composition is assumed to be given by local chemical equilibrium.

Fig. 1 Zonal-mean zonal wind speed (positive is eastward and negative westward) as a function of latitude and pressure, as calculated with a GCM simulation of HD 209458b (Parmentier et al. 2013). Note the superrotating wind above the 1 bar pressure level in the equatorial region (±20° in latitude).

Fig. 2 Zonal-mean zonal wind speed as a function of latitude and pressure calculated with a GCM simulation of HD 189733b (Parmentier et al., in prep.). A strong superrotating equatorial jet is clearly present in the 10–10-3 bar pressure range.

2.2.1. Wind structure

Circulation dynamics in the atmospheres of hot Jupiters is dominated by a fast eastward (or superrotating) jet stream at the equator. This superrotating jet was first predicted by Showman & Guillot (2002), has been found to emerge from almost all GCM simulations of hot Jupiters (Cooper & Showman 2005; Showman et al. 2008, 2009, 2013; Dobbs-Dixon & Lin 2008; Rauscher & Menou 2010, 2012a,b; Perna et al. 2010, 2012; Heng et al. 2011a,b; Lewis et al. 2010; Kataria et al. 2013; Parmentier et al. 2013), and is also understood theoretically (Showman & Polvani 2011). A shift of the hottest point of the planet eastward from the substellar point has been directly observed in several exoplanets (Knutson et al. 2007, 2009a, 2012; Crossfield et al. 2010) and interpreted as a direct consequence of this jet. The superrotating jet in HD 209458b and HD 189733b spans over all longitudes and has a well-defined location in latitude (around ±20°) and pressure (between 1–10 bar and 10^-6–10^-3 bar), as illustrated in Figs. 1 and 2, where the zonally averaged zonal wind speed is depicted as a function of latitude and pressure for each planet. To derive a mean speed of the jet for our pseudo two-dimensional chemical model, we averaged the zonal wind speed longitudinally over the whole planet, latitudinally over ±20°, and vertically between 1 and 10^-6 bar, the latter corresponding to the top of the atmosphere in the GCM simulations. We find mean zonal wind speeds of 3.85 km s^-1 for HD 209458b and 2.43 km s^-1 for HD 189733b, in both cases in the eastward direction. These values were adopted in the pseudo two-dimensional chemical model as the speed of the zonal wind at the equator and 1 bar pressure level, thus setting the angular velocity of the rotating vertical atmosphere column (i.e., its rotation period).

At high latitudes, above 50°, the circulation is no longer dominated by the superrotating jet; the zonal-mean zonal wind is westward and the flow exhibits a complex structure, with westward and eastward winds, and a substantial day-to-night flow over the poles (Showman et al. 2009). At low latitudes, the zonal-mean zonal wind is eastward over most of the vertical structure (above the 1–10 bar pressure level), although the shape of the superrotating jet changes gradually with altitude, from a well-defined banded flow with little longitudinal variability of the jet speed in the deep atmosphere to a less banded flow with important longitudinal variations of the wind speed in the upper levels (Showman et al. 2009). It is also worth noting that according to Showman et al. (2013), the circulation regime in the atmospheres of hot Jupiters changes from a superrotating one to a high-altitude day-to-night flow when the radiative or the frictional time scales become short, as occurs at low pressures under intense insolation or strong drag forces. In this regard, we note that the GCM simulations by Parmentier et al. (2013, and in prep.) used here are based on a drag-free case, and are thus the most favorable for the presence of a strong equatorial jet. That is, we would have found a somewhat slower equatorial jet if drag forces were included in the GCM simulations (Rauscher & Menou 2012a, 2013; Showman et al. 2013).

In view of the discussion above, our pseudo two-dimensional chemical model based on a uniform zonal wind probably is a good approximation for the equatorial region (±20°) in the 1 bar to 1 mbar pressure regime, and may still provide a reasonable description of upper equatorial layers, where an eastward jet is still present although with a less uniform structure. In the polar regions our formalism may not be adequate since the circulation regime is more complex, and thus the interplay between dynamics and chemistry may lead to a very different distribution of the chemical composition from that predicted by our model. It is interesting to note that low latitudes contribute more to the projected area of the planet’s disk than polar regions, and thus planetary emission is to a large extent dominated by the equatorial regions modeled here. The same is not true for transmission spectra however, where low and high latitudes are both important.

2.2.2. Temperature structure

Among the dozen hot Jupiters for which we have good observational constraints on their atmospheric properties (Seager & Deming 2010), half of them are believed to have a strong thermal inversion at low pressure in the dayside, while the other half are thought to lack such an inversion. The presence of a stratosphere in hot Jupiters is commonly attributed to the survival in the gas phase of the strong absorbers at visible wavelengths TiO and VO (Fortney et al. 2008; Showman et al. 2009; Parmentier et al. 2013). In this theoretical framework, planets that are warm enough to have an appreciable opacity due to TiO and VO (pM class planets) host a stratosphere, while those that are cooler (pL class planets) do not develop a temperature inversion in their atmospheres (Fortney et al. 2008). This, not yet firmly established however because no unambigous detection of TiO has been obtained (Désert et al. 2008), the nature of the absorbers that cause temperature inversions in hot Jupiters is still debated. For example, photochemical products of some undetermined nature or arising from the photochemical destruction of H₂S have also been postulated as possible absorbers responsible for these stratospheres (Burrows et al. 2008; Zahnle et al. 2009). Moreover, not all planets fit into this pM/pL scheme, and other parameters such as the atmospheric elemental C/O abundance ratio (Madhusudhan 2012) or the stellar activity (Knutson et al. 2010) might control whether there are stratospheres in hot Jupiters.

Fig. 3 Temperature structure averaged latitudinally over ±20° around the equator of HD 209458b, as calculated with a GCM simulation (Parmentier et al. 2013). Note the extremely hot dayside stratosphere above the 1 mbar pressure level.

Fig. 4 Temperature structure averaged latitudinally over ±20° around the equator of HD 189733b, as calculated with a GCM simulation (Parmentier et al., in prep.).

HD 209458b and HD 189733b are good examples of these two types of hot Jupiters, the former hosting a strong thermal inversion in the dayside, while the latter does not. The temperature resulting from the GCM simulations and averaged latitudinally over an equatorial band ±20° in latitude is shown as a function of longitude and pressure in Fig. 3 for HD 209458b and in Fig. 4 for HD 189733b. This equatorial band of ±20° in latitude corresponds to the region where the equatorial jet is present in the GCM simulations (see Figs. 1 and 2). For the pseudo two-dimensional chemical model, which focuses on the equatorial region where the eastward jet develops, we adopted the temperature distribution shown in Figs. 3 and 4, assuming an isothermal atmosphere at pressures lower than 2 μbar.

In our previous study (Agúndez et al. 2012), the temperature structure of HD 209458b’s atmosphere was calculated with a one-dimensional time-dependent radiative model and resulted in an atmosphere without a strong temperature inversion. Here, the temperature structure calculated for the atmospheres of HD 209458b and HD 189733b comes from GCM simulations, which result in a strong temperature inversion for the former planet and an atmosphere without stratosphere for the latter one. This permits us to explore the chemistry of hot Jupiters with and without a stratosphere. It is also worth noting that in the case of HD 209458b, there is evidence of a dayside temperature inversion from observations of the planetary emission spectrum at infrared wavelengths (Knutson et al. 2008).

2.2.3. Vertical eddy diffusion coefficient

Fig. 5 Vertical eddy diffusion coefficient profiles for HD 209458b and HD 189733b, as calculated by following the behavior of passive tracers (solid lines; Parmentier et al. 2013, and in prep.), and as given by previous estimates based on the rms of the vertical velocity times the vertical scale height (dashed lines; Moses et al. 2011).

Another important outcome of GCM simulations is the quantification of the strength with which material is transported in the vertical direction in the atmosphere. Although this mixing is not diffusive in a rigorous sense, once averaged over the whole planet, it can be well represented by an effective eddy diffusion coefficient that varies with pressure (Parmentier et al. 2013). This variable enters directly as input into one-dimensional and pseudo two-dimensional chemical models of planetary atmospheres such as ours. The eddy diffusion coefficient is commonly estimated in the literature as the root mean square of the vertical velocity times the vertical scale height (Line et al. 2010; Moses et al. 2011). Recently, Parmentier et al. (2013) have used a more rigorous approach to estimate an effective eddy diffusion coefficient in HD 209458b by following the behavior of passive tracers in a GCM. These authors have shown that vertical mixing in hot-Jupiter atmospheres is driven by large-scale circulation patterns. There are large regions with ascending motions and large regions with descending motions, some of them contributing more to the global mixing than others. It has been also shown that a diffusion coefficient is a good representation of the vertical mixing that takes place in the three-dimensional model of the atmosphere. The resulting values for HD 209458b are 10–100 times lower than those obtained with the previous method (see Fig. 5), and are used here. The vertical profile of the eddy diffusion coefficient for HD 209458b can be approximated by the expression K_zz (cm² s^-1) = 5 × 10⁸p^-0.5, where the pressure p is expressed in bar (see Parmentier et al. 2013). In the case of HD 189733b we used the expression K_zz (cm² s^-1) = 10⁷p^-0.65, where the pressure p is again expressed in bar. This expression is based on preliminary results by Parmentier et al. (in prep.) using the method involving passive tracers, and results in values up to 1000 times lower than those obtained with the previous more crude method (see Fig. 5). In both HD 209458b and HD 189733b we considered a constant K_zz value at pressures lower than 10^-5 bar.

2.2.4. Dynamical time scales

Fig. 6 Dynamical time scales of horizontal transport () and vertical mixing () as a function of pressure in HD 209458b and HD 189733b. is computed adopting as zonal wind speed either an average over an equatorial band of latitude ±20° (solid lines) or the mean values independent of height given in Sect. 2.2.1 (dotted lines).

To assess the relative strengths of horizontal transport and vertical mixing in the atmospheres of HD 209458b and HD 189733b it is useful to argue in terms of dynamical time scales. The dynamical time scale of horizontal transport may be roughly estimated as , where R_p is the planetary radius and u the zonal wind speed, while that related to vertical mixing can be approximated as , where H is the atmospheric scale height and K_zz the eddy diffusion coefficient. If we take a zonal wind speed uniform with altitude and equal to the mean value given in Sect. 2.2.1, we find that horizontal transport occurs faster than vertical mixing over most of the vertical structure of the atmospheres of HD 209458b and HD 189733b (see dotted and dashed lines in Fig. 6). Only in the upper layers, at pressures below 10^-3–10^-4 bar, the high eddy diffusion coefficient makes vertical mixing faster than horizontal transport.

Fig. 7 Vertical cuts of the abundance distributions of some of the most abundant molecules at longitudes spanning the 0–360° range, as calculated with the pseudo two-dimensional chemical model for HD 209458b’s atmosphere.

In the deep atmosphere, however, the equatorial superrotating jet vanishes and zonal winds become slower (see Figs. 1 and 2), although horizontal transport still remains faster than or at least similar to vertical mixing (see solid and dashed lines in Fig. 6). In these deep layers, below the 1–10 bar pressure level, our assumption of a zonal wind speed uniform with altitude and with values as high as a few km s^-1 is not valid. This is clearly a limitation of the pseudo two-dimensional model, although the implications for the resulting two-dimensional distribution of atmospheric constituents are not strong because in these deep layers the temperature remains rather uniform with longitude, and molecular abundances are largely controlled by thermochemical equilibrium, which makes them quite insensitive to the strength of horizontal transport. This has been verified by running models for HD 209458b and HD 189733b with zonal wind speeds down to 1000 times slower than the nominal mean values given in Sect. 2.2.1.

In the way the pseudo two-dimensional chemical model is conceived, it clearly deals with the equatorial region of hot Jupiter atmospheres. First, the formalism adopted, in which a vertical atmosphere column rotates around the equator at a constant angular velocity, is adequate for the equatorial region (±20° in latitude), where a strong eastward jet is found to dominate the circulation according to GCM simulations (see Figs. 1 and 2). Second, the temperature structure adopted (see Figs. 3 and 4) corresponds to the average over an equatorial band of width ±20° in latitude. Third, the rotation period of the atmosphere column is calculated from the wind speed retrieved from the GCM (which is also an average over an equatorial band ±20° in latitude) and the equatorial circumference. And fourth, the longitude-dependent zenith angle adopted to compute the penetration of stellar UV photons corresponds to the equatorial latitude. The adopted formalism is therefore adequate for the equatorial region as long as circulation is dominated by an eastward jet. With these limitations in mind, we now present and discuss the chemical composition distribution resulting from the pseudo two-dimensional chemical model for the atmospheres of HD 209458b and HD 189733b.

Fig. 8 Vertical cuts of the abundance distributions of some of the most abundant molecules at longitudes spanning the 0–360° range, as calculated with the pseudo two-dimensional chemical model for HD 189733b’s atmosphere.

Fig. 9 Distribution of H2O, CO2, CH4, and HCN as a function of longitude and pressure in the equatorial band of HD 209458b’s atmosphere, as calculated with the pseudo two-dimensional chemical model.

Fig. 10 Same as Fig. 9, but for HD 189733b.

Fig. 11 Vertical distributions of the most abundant atmospheric constituents, after H2 and He, at four longitudes: substellar point (0°), evening limb (+ 90°), antistellar point (± 180°), and morning limb (−90°) in the atmosphere of HD 209458b. We show the mole fractions calculated by the pseudo two-dimensional chemical model (solid lines), by a model that neglects vertical mixing (horizontal transport case; dashed-dotted lines), by a one-dimensional vertical model that neglects horizontal transport (vertical mixing case; dashed lines), and by local thermochemical equilibrium (dotted lines). Photochemistry is taken into account in all cases but the last.

3. Distribution of atmospheric constituents

3.1. Overview

A first glance at the calculated distribution of the chemical composition with altitude and longitude in the atmospheres of HD 209458b and HD 189733b can be obtained by examining the ranges over which the vertical abundance profiles vary with longitude. This information is shown in Figs. 7 and 8 for some of the most abundant species, after H₂ and He. We can see that some molecules such as CO, H₂O, and N₂ show little abundance variation with longitude, while some others such as CH₄, CO₂, NH₃, and HCN experience important changes in their abundances as longitude varies. Abundance variations are usually restricted to the upper regions of the atmosphere (above the 10^-1–10^-3 bar pressure level, depending on the molecule) but not to the lower atmosphere, where molecules maintain rather uniform abundances with longitude. On the one hand, longitudinal gradients in the temperature and incident stellar UV flux drive the abundance variations with longitude, while on the other, the zonal wind tends to homogenize the chemical composition in the longitudinal direction, resulting in the complex abundance distributions shown in Figs. 7 and 8. These results agree with the predictions of Cooper & Showman (2006) concerning the CO distribution in HD 209458b’s atmosphere. These authors coupled a GCM to a simple chemical kinetics scheme dealing with the interconversion between CO and CH₄ and found that CO shows a rather homogeneous distribution with longitude and latitude in spite of the strong variations predicted by chemical equilibrium. We also find a rather homogeneous distribution of CO with longitude, although the same is not true for other molecules that display important longitudinal abundance gradients.

Fig. 12 Same as Fig. 11, but for HD 189733b.

In Figs. 9 and 10 we show the atmospheric distribution of selected molecules that may influence planetary spectra as a function of longitude and pressure. Water vapor illustrates the case of a molecule with a rather uniform distribution throughout the atmosphere of both planets, except for a slight enhancement at high pressures (>10 bar) and a small depletion, which in HD 209458b occurs at about 10^-5 bar eastward of the substellar point and is induced by the warm stratosphere, and in HD 189733b takes place in the upper dayside layers (above the 10^-7 bar pressure level) through photochemical destruction. Carbon monoxide also has a quite uniform distribution and is not shown in Figs. 9 and 10. Carbon dioxide is perhaps the most abundant molecule showing important longitudinal abundance variations, with a marked day-to-night contrast. In HD 209458b this molecule is enhanced in the cooler nightside, where it is thermodynamically favored. In HD 189733b the nightside enhancement is only barely apparent in the 10^-5–10^-1 bar pressure range, while in upper layers the situation is reversed and CO₂ becomes depleted in the nightside regions because of a complex interplay between chemistry and dynamics. In the atmospheres of both planets, CO₂ maintains a mixing ratio between a few 10^-8 and a few 10^-5. The hydrides CH₄ and NH₃ show important abundance variations in the vertical direction, their abundance decrease when moving toward upper low-pressure layers, and also some longitudinal variability, which is only important at low abundance levels, however. In HD 209458b, methane is largely suppressed above the 1 mbar pressure level because of the stratosphere. In HD 189733b it is present at a more important level, except in the very upper layers where its depletion in the warmer dayside regions is propagated by the jet to the east, contaminating the nightside regions to a large extent. Hydrogen cyanide also shows important abundance variations with both longitude and altitude. This molecule is greatly enhanced by the action of photochemistry, and thus becomes quite abundant in the upper dayside regions of HD 189733b and to a lower extent in the upper dayside layers of HD 209458b, where photochemistry is largely supressed by the presence of the stratosphere (Moses et al. 2011; Venot et al. 2012). The distribution of HCN in the upper atmosphere shows that the eastward jet results in a contamination of nightside regions with HCN formed in the dayside.

As long as there is an important departure from chemical equilibrium in the atmospheric composition of both planets, the assumption of local chemical equilibrium in the GCM simulations may be an issue and one potential source of inconsistency between the GCM and the chemical model. Much of the thermal budget of these atmospheres, however, is controlled by water vapor, whose abundance is rather uniform and close to chemical equilibrium. This fact may justify to some extent the assumption of local chemical equilibrium in GCMs. We note however that other atmospheric constituents such as CO and CO₂ can also play an important role in the thermal balance of hot-Jupiter atmospheres, especially for elemental compositions far from solar, in which case the hypothesis of chemical equilibrium usually adopted in GCMs may not be adequate. Obviously, a more accurate and self-consistent approach would be to couple a robust chemical kinetics network to a GCM, although this is a very challenging computational task.

3.2. Comparison with limiting cases

To obtain insight into the predicted distribution of molecules in the atmospheres of HD 209458b and HD 189733b, a useful and pedagogical exercise is to compare the abundance distributions calculated by the pseudo two-dimensional chemical model with those predicted in various limiting cases. A first one in which vertical mixing is neglected and therefore the only disequilibrium processes are horizontal advection and photochemistry (horizontal transport case), a second one consisting of a one-dimensional vertical model including vertical mixing and photochemistry, which neglects horizontal transport (vertical mixing case), and a third one which is given by local thermochemical equilibrium. Figures 11 and 12 show the vertical distributions of some of the most abundant species, after H₂ and He, in the atmospheres of HD 209458b and HD 189733b, respectively, at four longitudes (substellar and antistellar points, and evening and morning limbs⁴), as calculated by the pseudo two-dimensional model and the three aforementioned limiting cases. We may summarize the effects of horizontal transport (modeled as a uniform zonal wind) and vertical mixing (modeled as an eddy diffusion process) by saying that horizontal transport tends to homogenize abundances in the horizontal direction, bringing them close to chemical equilibrium values of the hottest dayside regions, while vertical mixing tends to homogenize abundances in the vertical direction, bringing them close to chemical equilibrium values of hot bottom regions.

The effect of horizontal transport is perfectly illustrated in the case of methane. In both HD 209458b and HD 189733b, the abundance profile of CH₄ given by the horizontal transport case (blue dashed-dotted lines in Figs. 11 and 12) almost perfectly resembles the chemical equilibrium profile at the substellar point (blue dotted lines in upper left panel of Figs. 11 and 12), and remains almost invariant with longitude in spite of the important abundance enhancement predicted by chemical equilibrium in the cooler nightside and morning limb regions. The existence of a strong stratosphere in HD 209458b introduces important differences with respect to HD 189733b. The hot temperatures in the upper dayside layers of HD 209458b result in short chemical time scales and therefore allows chemical kinetics to mitigate to some extent the horizontal quenching induced by the zonal wind. This is clearly seen for CO₂ (magenta dashed-dotted lines in Figs. 11 and 12), whose abundance varies longitudinally within 2–3 orders of magnitude in HD 209458b, while in HD 189733b it shows an almost uniform distribution with longitude. The abundance distributions obtained in the horizontal transport case are qualitatively similar to those presented by Agúndez et al. (2012), although there are some quantitative differences due to the lack of photochemistry and a temperature inversion for HD 209458b in that previous study. Photochemistry plays in fact an important role in the horizontal transport case (dashed-dotted lines in Figs. 11 and 12), as it causes molecular abundances to vary with longitude in the upper layers due to the switch on/off of photochemistry in the day and night sides, as the wind surrounds the planet. Note also that the lack of vertical mixing in this case causes the photochemically active region to shift down to the level where, in the absence of vertical transport, chemical kinetics is able to counterbalance photodissociations, that is, to synthesize during the night the molecules that have been photodissociated during the day. Another interesting consequence of photochemistry in the horizontal transport case is that molecules such as HCN (gray dashed-dotted lines in Figs. 11 and 12), which are formed by photochemistry in the upper dayside regions, remain present in the upper nightside regions as a consequence of the continuous horizontal transport, and can in fact increase their abundances through the molecular synthesis ocurring during the night.

In the extreme case where vertical transport completely dominates over any kind of horizontal transport, the homogenization is produced in the vertical, and not longitudinal, direction. The value at which a given molecular abundance is quenched vertically corresponds to the chemical equilibrium abundance at the altitude where the rates of chemical reactions and vertical transport become similar, the so-called quench region. This quench region may be located at a different altitude for each species, although in hot Jupiters such as HD 209458b and HD 189733b it is usually located in the 10–10^-2 bar pressure range (Moses et al. 2011; Venot et al. 2012; also this study). Assuming the strength of vertical mixing does not vary with longitude (as done in this study), the vertical mixing case would yield uniform abundances with longitude if temperatures do not vary much with longitude in the range of altitudes where abundances are usually quenched vertically. In this case, the quench region for a given species would be the same at all longitudes, and so would the vertically quenched abundance. According to the GCM simulations of HD 209458b and HD 189733b, the temperature varies significantly with longitude above the 1 bar pressure level, and thus the exact values at which the abundances of the different species are quenched vertically vary with longitude. The temperature constrast between day and nightside regions is therefore one of the main causes of abundance variations with longitude, as illustrated by CH₄ in both planets (blue dashed lines in Figs. 11 and 12). Another factor that drives longitudinal abundance gradients in the vertical mixing case is photochemistry, which switches on and off in the day and nightsides, respectively. Without horizontal transport that connects the day and nightsides, abundances become rather flat in the vertical direction in the nightside, where photochemistry is suppressed, and display more complicated vertical abundance profiles in the dayside, where photochemistry causes molecules such as NH₃ to be depleted while some others such as HCN are enhanced. Note that because we used eddy diffusion coefficients significantly below those adopted in previous studies (e.g. Moses et al. 2011; Venot et al. 2012), the vertical quench of abundances in the dayside is not as apparent because it is strongly counterbalanced by photochemistry.

In the pseudo two-dimensional model, in which both horizontal transport and vertical mixing are simultaneously taken into account, the distribution of atmospheric constituents (solid lines in Figs. 11 and 12) results from the combined effect of various processes that tend to drive the chemical composition to a variety of distributions. On the one hand, chemical kinetics proceeds to drive the composition close to local chemical equilibrium. On the other hand, horizontal transport tends to homogenize abundances longitudinally, while vertical mixing does the same in the vertical direction. Finally, stellar UV photons tend to photodissociate molecules in the upper dayside layers, and new molecules are formed through chemical reactions involving the radicals produced in the photodissociations. Among these processes, horizontal transport and vertical mixing compete in homogenizing the chemical composition in the longitudinal and vertical directions, respectively. In the atmospheres of HD 209458b and HD 189733b horizontal transport occurs faster than vertical mixing below the ~1 mbar pressure level (see Sect. 2.2.4), and therefore molecular abundances are strongly homogenized in the longitudinal direction in this region. In upper layers the competition of mixing and photochemical processes results in a more complex distribution of atmospheric constituents.

Molecular abundances show a wide variety of behaviors when both horizontal transport and vertical mixing are considered simultaneously. The abundances of molecules such as CH₄, NH₃, and HCN tend to follow those given by the vertical mixing case at the substellar region, but at other longitudes the situation is quite different depending on the molecule (blue, green, and gray solid lines in Figs. 11 and 12). At the antistellar point, for example, the abundance profiles of CH₄ and NH₃ are closer to those predicted by the pure horizontal transport case than by the vertical mixing one, but HCN does follow a behavior completely different from each of these two limiting cases. The abundances of CO, H₂O, and N₂ show little variation with longitude or altitude and are therefore not affected by whether horizontal transport or vertical mixing dominates. Nevertheless, the coupling of horizontal transport and vertical mixing results in some curious behaviors, such as that of water vapor at the substellar point of HD 209458b (red lines in Fig. 11). The two limiting cases of pure horizontal transport and pure vertical mixing predict a decline in its abundance in the upper layers because of photodissociation and because of a low chemical equilibrium abundance at these low pressures. However, horizontal and vertical dynamics working simultaneously bring water from more humid regions so that there is no decline in its abundance up to the top of the atmosphere (at 10^-8 bar in our model). In summary, taking into account both horizontal transport and vertical mixing produces complex abundance distributions that in many cases cannot be predicted a priori.

Fig. 13 Vertical cuts of the abundance distributions of some of the most abundant molecules in HD 209458b’s atmosphere at longitudes spanning the 0–360° range, as calculated (from top to bottom) with the pseudo two-dimensional model, in the horizontal transport and vertical mixing cases, and under local chemical equilibrium.

Fig. 14 Same as Fig. 13, but for HD 189733b.

We may have a different view of the situation by looking at the ranges over which the vertical abundance profiles vary with longitude in the various limiting cases (see Figs. 13 and 14). Our attention first focuses on the fact that local chemical equilibrium predicts strong variations of the chemical composition with longitude in the atmospheres of both HD 209458b and HD 189733b. This is especially true for CH₄ and CO in the latter planet, where methane becomes more abundant than carbon monoxide in the cooler nightside regions. Disequilibrium processes, however, in particular horizontal transport and vertical mixing, reduce to a large extent the longitudinal variability of molecular abundances. As already stated, although perhaps more clearly seen in Figs. 13 and 14, horizontal transport tends to homogenize abundances with longitude. The effect of a purely horizontal transport is perfectly illustrated in HD 189733b’s atmosphere, where, except for the photochemically active region in the upper layers, the distribution of molecules is remarkably homogeneous with longitude (see horizontal transport panel in Fig. 14). In the atmosphere of HD 209458b, on the other hand, a pure horizontal transport allows for some longitudinal variability in the abundances of CO₂ and NH₃ above the 10^-3 bar pressure level (see horizontal transport panel in Fig. 13), mainly because of the activation of chemical kinetics in the dayside stratosphere and its ability to counterbalance the homogenization driven by horizontal transport. In the vertical mixing case (i.e., no horizontal transport), abundances are more uniform in the vertical direction but show important longitudinal variations, with a marked day/night asymmetry characterized by rather flat vertical abundance profiles in the nightside and abundances varying with altitude in the dayside because of the influence of photochemistry (see e.g. CH₄, NH₃, and HCN in vertical mixing panels of Figs. 13 and 14). When horizontal transport and vertical mixing are considered simultaneously (top panels of Figs. 13 and 14), the distribution of molecules in the lower atmosphere of both HD 209458b and HD 189733b, below the 10^-3 bar pressure level, remains remarkably homogeneous with longitude and close to that given by the pure vertical mixing case at the substellar regions. That is, the chemical composition of the hottest dayside regions propagates to the remaining longitudes, which indicates that the zonal wind transports material faster than vertical mixing processes do. In the upper atmosphere the abundance profiles become more complicated because of the combined effect of the photochemistry that takes place in the dayside and the mixing of material ocurring in both the vertical and horizontal directions.

3.3. Comparison with previous one-dimensional models

Fig. 15 Effect of eddy coefficient profile and 1D/2D character of the model on the vertical abundance profiles of some of the most abundant molecules in HD 209458b. We show abundances at the substellar point and at the evening and morning limbs, as given by the pseudo 2D model using the nominal eddy coefficient profile (solid lines), by the pseudo 2D model using the Moses et al. 2011’s eddy profile (dashed lines), and by a one-dimensional vertical model using the Moses et al. 2011’s eddy profile (dotted lines).

Fig. 16 Same as Fig. 15, but for HD 189733b.

It is interesting to compare the results obtained with the pseudo two-dimensional model with previous results from one-dimensional vertical models (Moses et al. 2011; Venot et al. 2012). There are two main differences between our model and these previous ones. The first is related to the eddy diffusion coefficients adopted, which are noticeably lower in this study because they are calculated by following the behavior of passive tracers in GCM simulations, while those used previously were also estimated from GCM simulations but as the root mean square of the vertical velocity times the vertical scale height. The second is related to the very nature of the model, which in our case is a pseudo two-dimensional model that simultaneously takes into account horizontal transport and vertical mixing, while in these previous studies horizontal transport is neglected. To isolate the differences caused by each of these factors we compare in Figs. 15 and 16 the vertical abundance distributions of some of the most abundant molecules at the substellar point and at the two limbs, as calculated with our pseudo two-dimensional model using the nominal vertical profile of the eddy diffusion coefficient (see Sect. 2.2.3), as given by the same pseudo two-dimensional model but using the high eddy diffusion coefficient profiles derived by Moses et al. (2011), which are about 10–100 times higher than ours for HD 209458b and about 10–1000 times higher than ours for HD 189733b, and as computed with a one-dimensional vertical model using the high eddy diffusivity values of Moses et al. (2011).

The main effect of increasing the strength of vertical mixing in the frame of a pseudo two-dimensional model is that the quench region shifts down to lower altitudes. For molecules such as CH₄, NH₃, and HCN, this implies that their vertically quenched abundances increase (compare solid and dashed lines in Figs. 15 and 16). If horizontal transport is completely supressed, that is, moving from a pseudo two-dimensional model to a one-dimensional vertical model (from dashed to dotted lines in Figs. 15 and 16), the horizontal homogenization of abundances is completely lost and thus the abundances of species such as CH₄ experience more important variations with longitude. Another interesting consequence of suppressing horizontal transport concerns water vapor, carbon monoxide, and molecular nitrogen, whose abundances decrease in the upper layers of the dayside regions (see upper panel in Figs. 15 and 16). In HD 209458b the depletion of these molecules is caused by the hot stratosphere, where neutral O and C atoms are favored over molecules, while in HD 189733b it is caused by photodissociation by UV photons. This loss of H₂O, CO, and N₂ molecules in the upper dayside layers is shifted to higher altitudes when horizontal transport, which brings molecules from other longitudes, is taken into account.

The vertical abundance profiles calculated with the one-dimensional vertical model at the substellar point and at the two limbs (dotted lines in Figs. 15 and 16) can be compared with the one-dimensional results of Moses et al. (2011) using averaged thermal profiles for the dayside and terminator regions (see also Venot et al. 2012). There is a good overall agreement between our substellar point results and their dayside results, on the one hand and on the other, between our results at the two limbs and their results at the terminator region, except for CH₄, NH₃, and HCN, for which we find vertically quenched abundances lower by about one order of magnitude. Although there are some differences in the adopted elemental abundances, stellar UV spectra, and zenith angles, the main source of the discrepancies is attributed to the different temperature profiles adopted. On the one hand there are slight differences between the GCM results of Showman et al. (2009), adopted by Moses et al. (2011) and Venot et al. (2012), and those of Parmentier et al. (2013, and in prep.), which are adopted here. On the other, and more importantly, the temperature profiles have a different nature. They are averages over the dayside and terminator regions in their case, while in ours they correspond to specific longitudes. The dayside average temperature profile of Moses et al. (2011) is cooler than our substellar temperature profile by about 100 K around the 1 bar pressure level in both planets, which results in vertically quenched abundances of CH₄ and NH₃ higher than ours by about one order of magnitude (part of the abundance differences are also due to the different chemical network adopted; see Fig. 7 of Venot et al. 2012). This serves to illustrate how relatively small changes of temperature in the 0.1−10 bar pressure regime – the quench region for most molecules – may induce important variations in the vertically quenched abundance of certain molecules. This also raises the question of whether it is convenient to use a temperature profile averaged over the dayside in one-dimensional chemical models that aim at obtaining a vertical distribution of molecules representative of the dayside. Although it may be a reasonable choice if one is limited by the one-dimensional character of the model, averaging the temperature over the whole dayside masks the temperatures of the hottest regions, near the substellar point, which are in fact the most important as they control much of the chemical composition at other longitudes if horizontal transport becomes important.

In summary, the main implications of using a pseudo two-dimensional approach and of the downward revision of the eddy values in the atmospheres of HD 209458b and HD 189733b are that, on the one hand, the longitudinal variability of the chemical composition is greatly reduced compared with the expectations of pure chemical equilibrium or one-dimensional vertical models and, on the other hand, the mixing ratios of CH₄, NH₃, and HCN are significantly reduced compared with results of previous one-dimensional models (by one order of magnitude or more with respect to the results of Moses et al. 2011), down to levels at which their influence on the planetary spectra are probably minor.

4. Calculated vs. observed molecular abundances

We now proceed to a discussion in which we compare the molecular abundances calculated with the pseudo two-dimensional chemical model and those derived from observations. Our main aim here is to evaluate whether or not the calculated composition, which is based on plausible physical and chemical grounds, is compatible with the mixing ratios derived by retrieval methods used to interpret the observations. The molecules H₂O, CO, CO₂, and CH₄ have all been claimed to be detected in the atmospheres of HD 209458b and HD 189733b either in the terminator region of the planet from primary transit observations, in the dayside from secondary eclipse observations, or in both regions using the two methods. Although we are not in a position to cast doubt on any of these detections, given the controversial results often found by different authors in the interpretation of spectra of exoplanets it is advisable to be cautious when using the derived mixing ratios to argue in any direction. Having this in mind, hereafter we use the term detection instead of claim of detection.

Water vapor and carbon monoxide are calculated with nearly their maximum possible abundances in both planets and show a rather homogeneous distribution as a function of both altitude and longitude (see Figs. 7 and 8). Adopting a solar elemental composition, as done here, the calculated mixing ratios of both H₂O and CO are around 5 × 10^-4. Water vapor being the species that provides most of the atmospheric opacity at infrared wavelengths, it was the first molecule to be detected in the atmosphere of an extrasolar planet, concretely in the transmission spectrum of HD 189733b (Tinetti et al. 2007), and H₂O mixing ratios derived from observations for both HD 189733b and HD 209458b are usually in the range of the calculated value of 5 × 10^-4 (Tinetti et al. 2007; Grillmair et al. 2008; Swain et al. 2008, 2009a,b; Madhusudhan & Seager 2009; Beaulieu et al. 2010; Lee et al. 2012; Line et al. 2013; Deming et al. 2013). Carbon monoxide, although less evident than water vapor, has also been detected in both planets and the mixing ratios derived are in the range of the values inferred for H₂O and expected from the chemical model (Swain et al. 2009a; Désert et al. 2009; Madhusudhan & Seager 2009; Lee et al. 2012; Line et al. 2013).

The calculated mixing ratio of carbon dioxide in the two hot Jupiters is in the range 10^-7–10^-6 depending on the pressure level, with a more important longitudinal variation in the atmosphere of HD 209458b than in that of HD 189733b (see Figs. 7 and 8). This molecule has been also detected through secondary-eclipse observations in the dayside of HD 189733b, with mixing ratios spanning a wide range from 10^-7 up to more than 10^-3 (Swain et al. 2009a; Madhusudhan & Seager 2009; Lee et al. 2012; Line et al. 2013), and in the dayside of HD 209458b, with a mixing ratio in the range 10^-6–10^-5 (Swain et al. 2009b). Taking into account the uncertainties associated with the values retrieved from observations, the agreement with the calculated abundance is reasonably good for CO₂.

The most important discrepancies between calculated and observed abundances are probably found for methane. This molecule is predicted to be very abundant in the cooler nightside regions of both planets, especially in HD 189733b, according to chemical equilibrium (see lower panels in Figs. 13 and 14), but reaches quite low abundances everywhere in the atmosphere according to the pseudo two-dimensional non-equilibrium model (see upper panels in Figs. 13 and 14). In both hot Jupiters, the calculated mixing ratio of CH₄ is in fact significantly lower than the predictions of previous one-dimensional models (Moses et al. 2011; Venot et al. 2012), a finding that strengthens the conflict with observations. We find that the mixing ratio of CH₄ above the 1 bar pressure level is below 10^-7 in HD 209458b and below 10^-6 in HD 189733b, whatever the side of the planet.

In HD 209458b, secondary-eclipse observations have been interpreted as evidence of methane being present in the dayside with a mixing ratio between 2 × 10^-5 and 2 × 10^-4 (Swain et al. 2009b), or within the less constraining range 4 × 10^-8–3 × 10^-2 (Madhusudhan & Seager 2009). In fact, the abundance of CH₄ retrieved in these studies is similar or even higher than that retrieved for H₂O, which is clearly not the case according to our predictions. It seems difficult to reconcile the low abundance of CH₄ calculated by the pseudo two-dimensional chemical model with the high methane content inferred from observations, which points to some fundamental problem in either of the two sides. As concerns the chemical model, an enhancement of the vertical transport to the levels adopted by Moses et al. (2011) or the supression of horizontal transport would increase the abundance of CH₄ only slightly (see Fig. 15). Photochemistry, which might potentially enhance the abundance of CH₄, is largely supressed by the stratosphere in the dayside atmosphere of HD 209458b. An elemental composition of the planetary atmosphere far from the solar one with, for example, an elemental C/O abundance ratio higher than 1, or some unidentified disequilibrium process, which might be related to, for instance, clouds or hazes, might lead to a high methane content in the warm atmospheric layers of HD 209458b’s dayside. Some problems on the observational side cannot be ruled out, taking into account the difficulties associated to the acquisition of photometric fluxes of exoplanets and the possibility of incomplete spectroscopic line lists of some molecules relevant to the interpretation of spectra of exoplanets (see e.g. the recently published line list for hot methane by Hargreaves et al. 2012).

In HD 189733b, contradictory results exist on the detection of methane in both the terminator and dayside regions. Swain et al. (2008) reported the detection of CH₄ through primary-transit observations, with a derived mixing ratio of about 5 × 10^-5, although Sing et al. (2009) did not find evidence of its presence in the transmission spectrum. These contradictory results obtained using NICMOS data could point to non-negligible systematics in the data (e.g., Gibson et al. 2012). Controversial results also exist on the detection of CH₄ in the dayside emission spectrum of HD 189733b (Swain et al. 2009a, 2010; Madhusudhan & Seager 2009; Waldmann et al. 2012; Lee et al. 2012; Line et al. 2013; Birkby et al. 2013). Until observations can draw more reliable conclusions it is difficult to decide whether or not observations and models are in conflict regarding the abundance of CH₄ in HD 189733b.

5. Variations in the planetary spectra

Fig. 17 Pressure level probed by transmission and emission spectra as a function of wavelength for HD 209458b and HD 189733b. Dashed lines correspond to a model with a mean vertical profile averaged over the terminator and show the pressure level at which the tangential optical depth equals 2/3, which is in fact a transmission spectrum expressed in terms of atmospheric pressure instead of planetary radius. Solid lines correspond to a model with a mean vertical profile averaged over the dayside and indicate the pressure level at which the optical depth in the vertical outward direction equals 2/3, which is an approximate location of the region from where most of the planetary emission arises, the regions below being opaque and those above being translucent.

Fig. 18 Calculated emission spectra for HD 209458b (left panel) and HD 189733b (right panel) at 4 different phases in which the planet faces the observer the day and night sides, and the two sides in between, centered on the evening and morning limbs. Spectra have been smoothed to a resolving power R = 300. The adopted vertical profiles of temperature and mixing ratios are projected area-weighted averages over the corresponding emitting hemisphere, with the distribution of molecules given by our nominal pseudo two-dimensional chemical model. The planetary flux is shown relative to that of the star, for which the Kurucz synthetic spectrum and stellar radius given in Sect. 2 are adopted. Gray dashed lines correspond to planetary blackbody temperatures of 2000, 1500, and 1000 K (from top to bottom in HD 209458b’s panel) and of 1000 K (in HD 189733b’s panel). The inset in HD 209458b’s panel compares the emission spectrum around 4.3 μm as calculated for the dayside and as computed using the mean dayside temperature profile and the mean chemical composition of the nightside. The two spectra are nearly identical except for a slight difference at 4.3 μm due to CO2.

The calculated distribution of molecules in the atmospheres of HD 209458b and HD 189733b may be probed by observations. Instead of comparing synthetic spectra and available observations of these two planets, we are here mainly interested in evaluating whether the longitudinal variability of the chemical composition may be probed by observations. For example, the monitoring of the emission spectrum of the planet at different phases during an orbital period would probe the composition in the different sides of the planet. In addition, the observation of the transmission spectrum at the ingress and egress during primary transit conditions would allow one to probe possible chemical differentiation between the morning and evening limbs of the planet’s terminator. Planetary emission and transmission spectra were computed using the line-by-line radiative transfer code described in Appendix A. Since the code is currently limited because it is one-dimensional in the vertical direction, we adopted mean vertical profiles by averaging the temperature structure in longitude and latitude given by the GCM simulations of Parmentier et al. (2013, and in prep.) and the longitudinal distribution of abundances obtained with the pseudo two-dimensional chemical model. In the case of emission spectra, we adopted a weighted average profile of temperature and of mixing ratios over the hemisphere facing the observer (weighted by the projected area on the planetary disk to better represent the situation encountered by an observer), where mixing ratios were assumed to be uniform with latitude. In transmission spectra, vertical profiles are simply obtained by averaging over the whole terminator, or over the morning or evening limb. After adopting an average pressure-temperature profile, the planetary radius (see values in Sect. 2) is assigned to the 1 bar pressure level and the altitude of each layer in the atmosphere is computed according to hydrostatic equilibrium. We note that similarly to the case of one-dimensional chemical models, the use of average vertical profiles in calculating planetary spectra is an approximation that masks the longitudinal and latitudinal structure of temperature and chemical composition and may result in non-negligible inaccuracies in the appearance of the spectra, which we plan to investigate in the future.

Is is useful to begin our discussion on planetary spectra with a pedagogical plot that shows the pressure level probed by transmission and emission spectra for HD 209458b and HD 189733b (see Fig. 17). We first note that transmission and emission spectra probe different pressure levels, with transmission spectra being sensitive to upper atmospheric layers than emission spectra. At infrared wavelengths (1–30 μm), and for the thermal and chemical composition profiles adopted by us for these two hot Jupiters, emission spectra probe pressures between 10 and 10^-2 bar, while transmission spectra probes the 1–10^-3 bar pressure regime. A second aspect worth noting is that there are strong variations with wavelength in both types of spectra, which implies that observations at different wavelengths are sensitive to the physical and chemical conditions of different pressure levels. It is always useful to keep these ideas in mind when analyzing the vertical distribution of molecules calculated with a chemical model, because only a very specific region of the atmosphere becomes relevant to planetary spectra.

5.1. Variation of emission spectra with phase

A modulation of the planetary emission with the orbital phase has been observed for HD 189733b by monitoring the photometric flux in the 8 μm band of Spitzer IRAC during a good part of the orbit of the planet (Knutson et al. 2007). This has served to evidence the important temperature contrast between the different sides of the planet, noticeably between day and night, and indirectly the presence of strong winds that can redistribute the energy from the day to the night side, due to an observed shift between the hot spot and the substellar point. Various theoretical studies have also been interested in predicting the variation of the planetary flux with the orbital phase in HD 209458b and HD 189733b using the temperature structure calculated with GCM simulations (Fortney et al. 2006; Showman et al. 2008, 2009; Rauscher et al. 2008; Burrows et al. 2010; Rauscher & Menou 2013). Most previous studies have focused on the link between light curves and variations of temperature between the different planetary sides, and on the comparison between predicted and observed photometric fluxes. Here we are instead interested in discussing the influence of the temperature but also that of the chemical composition (assumed to be given by chemical equilibrium in previous studies) on the variation of the planetary emission with phase.

Fig. 19 Calculated transmission spectra for the evening and morning limbs of HD 209458b (left panel) and HD 189733b (right panel), where the vertical structure is obtained by averaging the temperature over each limb and adopting the abundance profiles at each limb from the nominal pseudo two-dimensional chemical model. Spectra have been smoothed to a resolving power R = 300. The transit depth is simply calculated as (Rp(λ)/R∗)2, where Rp(λ) is the calculated radius of the planet as a function of wavelength and R∗ is the stellar radius (see values in Sect. 2). The absolute scale of the transmission spectrum is set by our choice of assigning the value of the planetary radius given in Sect. 2 to the 1 bar pressure level. Since no attempt has been made to reproduce the absolute scale indicated by primary transit observations, calculated transit depths are somewhat higher than given by observations. The insets in both panels compare the transmission spectrum around 4.3 μm as calculated for the evening limb and as computed using the mean temperature profile of the evening limb and the chemical composition corresponding to the morning limb. The most important differences between both spectra occur around 4.3 and 15 μm, due to CO2.

We show in Fig. 18 how the calculated emission spectra of HD 209458b and HD 189733b vary with the phase of the planet. Important variations with phase are apparent in HD 209458b, whose strong dayside stratosphere causes the dayside emission spectrum to be significantly brighter and to have a noticeably different spectral shape than at other phases. In HD 189733b, the modulation of the flux and the variation of the spectral shape with phase are also important although less pronounced. Emission spectra are controlled on the one hand, by the vertical temperature structure, and on the other, by the abundances of the main atmospheric constituents providing opacity. The sensitivity of emission spectra to the thermal structure is illustrated in HD 209458b, whose dayside (facing a temperature inversion to the observer) shows some spectral features that appear in emission and not in absorption, as occurs for the other planetary sides of HD 209458b and for HD 189733b. These differences in the spectra can be used to infer whether there is a stratosphere in the atmosphere of a hot Jupiter from observations of its dayside emission spectrum (Knutson et al. 2008, 2009b; Burrows et al. 2008; Machalek et al. 2008, 2010; Todorov et al. 2010; Madhusudhan & Seager 2010). It is also interesting to note how similar the emission spectra of night and morning sides are in the two planets, as are the day and evening sides in the case of HD 189733b. This is a consequence of the eastward transport of energy by the superrotating jet, which shifts the hottest and coldest regions to the east of the substellar and antistellar points, respectively. In the calculated emission spectra of both HD 209458b and HD 189733b, most of the atmospheric opacity along the 1–30 μm wavelength range is provided by water vapor, with carbon monoxide contributing at 2.3 and 4.6 μm, CO₂ at 4.3 and 15 μm, and collision-induced absorption by H₂-H₂ in certain wavelength ranges below 4 μm. No other species leaves appreciable signatures in the calculated emission spectra of HD 209458b, although in that of HD 189733b CH₄ contributes around 3.3 and 7.7 μm, NH₃ around 10.6 μm, and HCN at 14 μm.

An interesting question that arises from the change in the emission spectrum with phase is whether it is entirely caused by the variation of temperature in the different sides of the planet or whether the longitudinal variation of the chemical composition contributes to an important extent. To illustrate this point we compare in the inset of HD 209458b’s panel in Fig. 18 the dayside emission spectrum with a synthetic spectrum calculated using the mean vertical temperature structure of the dayside and the mean chemical composition of the nightside. The two spectra are nearly identical except for a slight difference around 4.3 μm, a spectral region where atmospheric opacity is to a large extent dominated by CO₂. Similar models in which the temperature structure and the chemical composition are adopted from different planetary sides indicate that variations of HD 209458b’s emission spectrum with phase are almost entirely caused by changes of temperature, with the only effect that can be purely adscribed to variations in the chemical composition being restricted to the tiny variation (less than 0.02% in the planet-to-star flux ratio) at 4.3 μm, which is caused by the longitudinal variation of about one order of magnitude in the abundance of CO₂ (see Fig. 7). The reasons of the small impact of the chemical composition on the variation of emission spectra with phase are related to the important longitudinal homogenization of the abundances driven by the zonal wind in HD 209458b (see Fig. 7). In fact, most of the atmospheric opacity affecting the emission spectrum comes from H₂O, CO, and CO₂, in order of decreasing importance, and the two former molecules show remarkably uniform abundances with longitude, while only the abundance of the latter molecule experiences some longitudinal variation, leading to a slight variation of the planetary flux with phase around 4.3 μm. In HD 189733b, the homogenization of the chemical composition with longitude is even more marked than in HD 209458b because of the lack of a stratosphere and the rather low eddy coefficient values (compare Figs. 7 and 8). Because the abundances of H₂O, CO, CO₂, CH₄, NH₃, and HCN (the main molecules providing opacity, in order of decreasing importance) vary little between the different sides of HD 189733b at the pressures probed by emission spectra (>10^-2 bar), the impact of the chemical composition on the change of the emission spectrum with phase becomes almost negligible, even around 4.3 μm because of the reduced longitudinal variation of the abundance of CO₂.

5.2. Transmission spectra of evening and morning limbs

Fig. 20 Contributions of the different sources of opacity to the transmission spectra of HD 209458b (left panel) and HD 189733b (right panel) smoothed to a resolving power R = 300. The H2-He continuum, whose contribution is similar in shape to the H2-H2 continuum but lower because of the lower abundance of He with respect to H2, is not shown. The spectra have been calculated using the temperature from GCM simulations and the chemical composition from the nominal pseudo two-dimensional chemical model averaged over the whole terminator. Each line represents the transmission spectrum that results from a model in which only the opacity provided by each source is taken into account.

Variations in the composition of the atmosphere between the different sides of the planet may also be probed by transmission spectroscopy. Indeed, it is a priori possible to probe differences in the thermal and chemical structure of the two limbs if observations are able to obtain the transmission spectrum during the first half of the primary transit ingress, which would probe the leading or morning limb, and during the second half of the transit egress, which would probe the trailing or evening limb. Although a non-zero impact parameter during the transit would complicate the situation somewhat and such observations are very challenging today, they may be feasible in the near future. The subject has been addressed theoretically for hot Jupiters such as HD 189733b and HD 209458b in some studies in which the differences between the transmission spectra of leading and trailing limbs are evaluated under different assumptions for the chemical composition of each of the two limbs, either chemical equilibrium or some disequilibrium estimation (Fortney et al. 2010; Burrows et al. 2010). Here we revisit the subject in the light of the molecular abundances calculated in this study with the pseudo two-dimensional chemical model.

To illustrate the possibility that transmission spectroscopy might be able to distinguish between the two limbs of HD 209458b and HD 189733b we show in Fig. 19 the transmission spectrum calculated by adopting the chemical composition and mean temperature of the evening and morning limbs of these two exoplanets. Since we are mainly interested in comparing the spectra at the two limbs and not in comparing with observations, we set the absolute scale of transmission spectra by simply assigning the value of the planetary radius given in Sect. 2 to the 1 bar pressure level and made no attempt to reproduce the absolute scale of the photometric transit depths derived from observations. In HD 209458b and HD 189733b, the transmission spectrum of the evening limb shows a higher degree of absorption and also the variations of the transit depth with wavelength have a larger amplitude than at the morning limb, whose transmission spectrum is flatter. These differences are mainly due to the different temperature profile of the two limbs. Because the atmosphere at the morning limb is cooler and thus has a smaller scale height than at the evening limb, it becomes more compact, resulting in smaller apparent radii at all wavelengths and a flatter transmission spectrum. In addition to this dependence of the transmission spectrum on temperature, which causes it to shift up or down and to have a more elongated or flattened overall shape, the spectral structure is controlled by the relative abundances of the main species that provide opacity in the atmosphere. Figure 20 shows the relative contributions of the different sources of opacity taken into account in calculating the transmission spectra. Similarly to the emission spectra, in the calculated transmission spectra of HD 209458b and HD 189733b most of the atmospheric opacity at infrared wavelengths is provided by H₂O, with CO being important at 2.3 and 4.6 μm, CO₂ at 4.3 and 15 μm, the H₂-H₂ continuum at certain wavelengths below 3 μm, and, in the case of HD 189733b, CH₄ having some contribution around 3.3 and 7.7 μm, NH₃ around 10.6 μm, and HCN at 14 μm.

Similarly to emission spectra, we evaluated to which extent differences in the chemical composition of evening and morning limbs contribute to the change of the transmission spectrum from one limb to the other. To this purpose, we computed transmission spectra in which we switched the temperature and chemical profiles between the two different limbs. As an example we compare in the insets of left and right panels in Fig. 19 the transmission spectrum of the evening limb with a synthetic spectrum calculated using the temperature structure of the evening limb and the chemical composition of the morning limb. In HD 209458b the two spectra are very similar, except for a different degree of absorption around 4.3 and 15 μm, which is due to the difference of nearly one order of magnitude in the abundance of CO₂ between the two limbs (see Fig. 11). Similarly to emission spectra, the longitudinal homogenization driven by horizontal transport is at the origin of the weak impact of other molecules on the variation of transmission spectra between both limbs. Because the abundances of CO and H₂O are very similar in both limbs, and other molecules such as CH₄, NH₃, and HCN contribute little to the atmospheric opacity at infrared wavelengths because of their rather low abundances, the only chemical effect contributing to the change of the transmission spectrum from one limb to the other of HD 209458b is restricted to carbon dioxide. In HD 189733b, the even stronger longitudinal homogenization of abundances (compare Figs. 7 and 8) diminishes the extent of chemical effects, which are now restricted to a very weak change of the absorption around 4.3 and 15 μm (see inset in HD 189733b’s panel of Fig. 19), again due to a slight increase in the abundance of CO₂ when moving from the evening limb to the morning one (see Fig. 12).

6. Summary

We have developed a pseudo two-dimensional model of a planetary atmosphere that takes into account thermochemical kinetics, photochemistry, vertical mixing, and horizontal transport, and allows one to calculate the distribution with altitude and longitude of the main atmospheric constituents. Horizontal transport was modeled through a uniform zonal wind and thus the model is best suited for studying the atmosphere of planets whose circulation dynamics is dominated by an equatorial superrotating jet, as is expected to be the case of hot Jupiters. We therefore applied the model to study the atmospheres of the well-known exoplanets HD 209458b and HD 189733b. We used the temperature structure from GCM simulations and parameterized the turbulent mixing in the vertical direction using an eddy coefficient profile, which was calculated by following the behavior of passive tracers in GCM simulations, a method that results in substantially lower eddy values, by a factor of 10–100 in HD 209458b and of 10–1000 in HD 189733b, than previous estimates based on cruder methods.

Molecular abundances homogenized with longitude to values typical of the hottest dayside regions. – We found that the distribution of molecules in the atmospheres of HD 209458b and HD 189733b is quite complex because of the interplay of the various (photo)chemical and dynamical processes at work, which form, destroy, and transport molecules throughout the atmosphere. Much of the distribution of the atmospheric constituents is driven by the strong zonal wind, which reaches speeds of a few km s^-1, and the limited extent of vertical transport, with relatively low eddy diffusion coefficients below 10⁹ cm² s^-1 around the 1 bar pressure level, resulting in an important homogenization of molecular abundances with longitude, in particular in the atmosphere of HD 189733b, which lacks a stratosphere and has quite low eddy diffusion coefficients. Moreover, molecular abundances are quenched horizontally to values typical of the hottest dayside regions, and therefore the composition of the cooler nightside regions is highly contaminated by that of warmer dayside regions. In hot Jupiters with a temperature inversion, such as HD 209458b, the longitudinal homogenization of molecular abundances is not as marked as in planets lacking a stratosphere, such as HD 189733b. In general, the cooler the planet, the stronger the homogenization of the chemical composition with longitude. Furthermore, in cooler planets such as hot Neptunes orbiting M dwarfs (e.g., GJ 436b) the temperature contrast between day and nightsides decreases because the cooling rate scales with the cube of temperature (e.g., Lewis et al. 2010), and therefore the composition is expected to be even more homogeneous with longitude than in warmer planets such as HD 209458b and HD 189733b. However, unlike hot Jupiters, hot Neptunes may have an atmospheric metallicity much higher than solar (Line et al. 2011; Moses et al. 2013; Agúndez et al. 2014; Venot et al. 2014), which makes it interesting to investigate the extent of the spatial variation of molecular abundances in their atmospheres.

Low methane content. – A major consequence of our pseudo two-dimensional chemical model is that methane reaches quite low abundances in the atmospheres of HD 209458b and HD 189733b, lower than the values calculated by previous one-dimensional models. The main reason for the low CH₄ abundance is that most of the atmosphere is contaminated by the hottest dayside regions, where the chemical equilibrium abundance of CH₄ is the lowest. The calculated mixing ratio of CH₄ in the dayside of HD 209458b is significantly below the values inferred from observations, which points to some fundamental problem in either the chemical model or the observational side. If the strength of vertical transport is substantially higher than in our nominal model, the calculated abundance of some molecules such as CH₄ and NH₃ would experience significant enhancement, especially in HD 189733b, although a conflict with observations would still exist regarding CH₄ in the dayside of HD 209458b.

Variability of planetary spectra driven by thermal, rather than chemical, gradients. – An important consequence of the strong longitudinal homogenization of molecular abundances in the atmospheres of HD 209458b and HD 189733b is that the variability of the chemical composition has little effect on the way the emission spectrum is modified with phase and on the changes of the transmission spectrum from the transit ingress to the egress. Temperature variations and not chemical gradients are therefore at the origin of these types of variations in the planetary spectra. Only the longitudinal variation of the abundance of CO₂, of nearly one order of magnitude, in the atmosphere of HD 209458b, is predicted to induce variations in the planetary spectra around 4.3 and 15 μm. We note, however, that an inhomogenous distribution of clouds and/or hazes (none of them included in our model) may induce important variations in the emission spectra with phase and in the transmission spectra from one limb to the other. These variations are best characterized at short wavelengths. Indeed, there is evidence of the presence of hazes in the atmosphere of HD 189733b (Lecavelier des Etangs et al. 2008; Sing et al. 2009), and an inhomogeneous distribution of clouds has recently been inferred for the hot Jupiter Kepler 7b (Demory et al. 2013).

The main drawback of our pseudo two-dimensional chemical model is the oversimplification of atmospheric dynamics, which is probably adequate for equatorial regions, but not at high latitudes. Ideally, GCM simulations coupled to a robust chemical network would provide an even more realistic view of the distribution of molecules in the atmospheres of HD 209458b and HD 189733b, but such calculations are very challenging from a computational point of view. Telescope facilities planned for the near or more distant future, such as the James Webb Space Telescope, Spica, and EChO, will be able to test some of the predictions of our pseudo two-dimensional model, in particular the low abundance of methane in the two planets and the important longitudinal homogenization of the chemical composition.

Acknowledgments

We thank our anonymous referee for insightful comments which helped to improve this article. We acknowledge Adam P. Showman and Jonathan J. Fortney for the use of the SPARC/MITgcm code, Vincent Hue for useful discussions on photochemical models, Sergio Blanco-Cuaresma and Christophe Cossou for their help with Python and Fortran, and Vincent Eymet and Philip von Paris for kindly helping to validate the line–by–line radiative transfer code. M.A. and F.S. acknowledge support from the European Research Council (ERC Grant 209622: E₃ARTHs). O.V. acknowledges support from the KU Leuven IDO project IDO/10/2013 and from the FWO Postdoctoral Fellowship Program. Computer time for this study was provided by the computing facilities MCIA (Mésocentre de Calcul Intensif Aquitain) of the Université de Bordeaux and of the Université de Pau et des Pays de l’Adour.

References

Appendix A: Calculation of planetary spectra

To investigate the influence of the physical and chemical structure of the atmospheres of HD 209458b and HD 189733b on their transmission and emission spectra, we developed a line-by-line radiative transfer code that is independent of the pseudo two-dimensional chemical code. Currently, the code is limited because it is one-dimensional in the sense that the atmosphere is divided into various layers in the vertical direction (typically 60 spanning the 10–10^-6 bar pressure range) and each layer is assumed to be homogeneous with longitude and latitude. Therefore, each layer is characterized by a given pressure, temperature, and chemical composition, and longitudinal and latitudinal gradients are neglected. For transmission spectra, the physical and chemical profile in the vertical direction at either the east or west limb, or a mean of the profiles at both limbs, can be used. For emission spectra, the limitations caused by the one-dimensional character of the code can be partially alleviated by adopting thermal and chemical vertical profiles averaged in some manner over the hemisphere facing the observer.

It is common in infrared spectroscopy to use the wavenumber with units of cm^-1, instead of frequency or wavelength, and we therefore adopt this choice hereafter as well. At this stage, there are various sources of opacity included in the code. On the one hand, we consider collision induced absorption (CIA) by H₂-H₂, for which available absorption coefficients cover the wavelength range 10–25 000 cm^-1 and temperatures between 60 and 7000 K (Borysow et al. 2001; Borysow 2002), and by H₂-He, in which case absorption coefficients in the wavelength range 10–25 000 cm^-1 and for temperatures in the range 100–7000 K are available (Borysow et al. 1989, 1997; Borysow & Frommhold 1989). CIA absorption coefficients scale with the square of pressure and thus become the dominant source of opacity at high pressures, usually above 1 bar. On the other hand, we consider spectroscopic transitions (mostly ro-vibrational transitions lying at infrared wavelengths) of H₂O, CO, and CO₂, whose data are taken from HITEMP (Rothman et al. 2010), and of CH₄, NH₃, and HCN, for which data from HITRAN (Rothman et al. 2009) are adopted.

The spectral region of interest is divided into a certain number of spectral bins, whose widths are determined by the spectral resolution imposed. In each layer of the atmosphere, the contribution of a spectroscopic transition j (centered on a wavenumber and which belongs to a species i) to the absorption coefficient in a spectral bin l (centered on a wavenumber and having a width ), which we may label as , can be expressed as (A.1)where S_j(T) is the line intensity of the spectroscopic transition j, which depends on the temperature T and is usually given with units of cm^-1/(molecule cm^-2) in the HITRAN and HITEMP databases, n_i is the number density of species i in the atmospheric layer where the absorption coefficient is to be evaluated, and φ_j is the line profile function of transition j, which has units of inverse of wavenumber (i.e., cm) and must be normalized such that the integral of from to yields unity. The absorption coefficient has units of cm^-1. The integral in Eq. (A.1) extends between the lower and upper wavenumber edges of the spectral bin l. The line profile function is taken as a Voigt profile, which results from the convolution of a Gaussian and a Lorentzian profile, and thus accounts for the Doppler and pressure broadening of spectral lines in each layer of the atmosphere due to thermal motions and collisions, respectively. The Voigt profile function is calculated numerically with a routine based on an implementation of Humlicek’s algorithm by Kuntz (1997).

Ideally, extremely high spectral resolution would be desirable to properly resolve the narrowest line profiles, although in practice this is too expensive in terms of computing time. For the calculations carried out here we adopted a spectral resolution of 0.03 cm^-1, which is on the order of the line widths in the layers that provide most atmospheric opacity and has been found to be high enough to yield relative errors below 1% in the computed spectra. Spectroscopic transitions lying farther away than 50 cm^-1 of a given spectral bin l have not been taken into account when computing the absorption coefficient . In addition, a cutoff in the line intensity S_j(296 K) of 10^-40 cm^-1/(molecule cm^-2) was adopted to neglect weak lines and speed up the calculations. This is perhaps the most delicate aspect because weak lines are numerous, especially in the HITEMP line lists, and at high temperatures the sum of all them results in non-negligible opacity enhancements in certain spectral regions. Calculations carried out with different line intensity cutoffs in selected spectral regions indicate that the relative error in the calculated spectra is at most 10% in the hottest case studied here (the dayside emission spectrum of HD 209458b), and lower than 5% in the remaining computed spectra. Currently, the code does not take into account light scattering, which becomes important at wavelengths shorter than ~1 μm, and therefore calculated spectra are reliable at infrared wavelengths, but not in the visible region of the electromagnetic spectrum. The code does not take into account the Doppler shift of spectral lines due to atmospheric winds either, an effect that may be observable in high-resolution spectra of hot Jupiters, where winds are strong (Snellen et al. 2010; Miller-Ricci Kempton & Rauscher 2012; Showman et al. 2013).

At each atmospheric layer, the absorption coefficient in each wavenumber bin l is calculated as the sum of the contributions from all the spectroscopic transitions of the various absorbing species included, together with the contributions of the CIA couples, that is, (A.2)where the sum in m extends to the H₂-H₂ and H₂-He couples.

After calculating the absorption coefficient in each wavenumber bin and at each atmospheric layer, computing the transmission and emission planetary spectra becomes straightforward, provided scattering is not considered. To calculate the transmission spectrum, the optical depth along a tangential line of sight intersecting the planet’s atmosphere is computed as a function of the impact parameter b in each spectral interval l as (A.3)where Δℓ_h(b) is the path length (in cm) intersected by layer h along the tangential line of sight at impact parameter b, and the sum in h extends to all atmospheric layers N(b) intersected by the tangential line of sight. The apparent radius of the planet in each wavenumber interval is then retrieved as the impact parameter for which the optical depth along the tangential line of sight becomes 2/3. This latter value is rather arbitrary and is mainly chosen for similarity with the definition of a stellar photosphere, although it is very close to the value of 0.56 inferred by Lecavelier des Etangs et al. (2008) and is not critical to derive the apparent radius of the planet. To obtain the emission spectrum, the emergent specific intensity along the observer’s line of sight is computed as a function of impact parameter in each wavenumber bin. We thus need to solve the equation of radiative transfer along the various paths pointing toward the observer that pass through the planetary atmosphere at different impact parameters. As long as the different atmospheric layers are homogeneous, the equation of radiative transfer can be solved sequentially from the back to the front for each of the atmospheric layers intersected by the path. For each intersected layer, the equation of radiative transfer reads (see e.g. Rybicki & Lightman 2004) (A.4)where the subscript indicates wavenumber dependence, and are the incoming and outgoing specific intensities that enter and emerge, respectively, from the current layer along a given path, is the optical depth along the path within the current layer, and is Planck’s function. The final emission spectrum is then computed by averaging the wavenumber-dependent specific intensity, calculated as a function of the impact parameter, over the projected area of the emitting hemisphere.

The line-by-line radiative transfer code was checked against the suite of radiative transfer tools “kspectrum”⁵, which has been widely used to model the atmosphere of solar system planets such as Venus (Eymet et al. 2009).

All Figures

	Fig. 1 Zonal-mean zonal wind speed (positive is eastward and negative westward) as a function of latitude and pressure, as calculated with a GCM simulation of HD 209458b (Parmentier et al. 2013). Note the superrotating wind above the 1 bar pressure level in the equatorial region (±20° in latitude).
In the text

	Fig. 2 Zonal-mean zonal wind speed as a function of latitude and pressure calculated with a GCM simulation of HD 189733b (Parmentier et al., in prep.). A strong superrotating equatorial jet is clearly present in the 10–10-3 bar pressure range.
In the text

	Fig. 3 Temperature structure averaged latitudinally over ±20° around the equator of HD 209458b, as calculated with a GCM simulation (Parmentier et al. 2013). Note the extremely hot dayside stratosphere above the 1 mbar pressure level.
In the text

	Fig. 4 Temperature structure averaged latitudinally over ±20° around the equator of HD 189733b, as calculated with a GCM simulation (Parmentier et al., in prep.).
In the text

	Fig. 5 Vertical eddy diffusion coefficient profiles for HD 209458b and HD 189733b, as calculated by following the behavior of passive tracers (solid lines; Parmentier et al. 2013, and in prep.), and as given by previous estimates based on the rms of the vertical velocity times the vertical scale height (dashed lines; Moses et al. 2011).
In the text

	Fig. 6 Dynamical time scales of horizontal transport () and vertical mixing () as a function of pressure in HD 209458b and HD 189733b. is computed adopting as zonal wind speed either an average over an equatorial band of latitude ±20° (solid lines) or the mean values independent of height given in Sect. 2.2.1 (dotted lines).
In the text

	Fig. 7 Vertical cuts of the abundance distributions of some of the most abundant molecules at longitudes spanning the 0–360° range, as calculated with the pseudo two-dimensional chemical model for HD 209458b’s atmosphere.
In the text

	Fig. 8 Vertical cuts of the abundance distributions of some of the most abundant molecules at longitudes spanning the 0–360° range, as calculated with the pseudo two-dimensional chemical model for HD 189733b’s atmosphere.
In the text

	Fig. 9 Distribution of H2O, CO2, CH4, and HCN as a function of longitude and pressure in the equatorial band of HD 209458b’s atmosphere, as calculated with the pseudo two-dimensional chemical model.
In the text

	Fig. 10 Same as Fig. 9, but for HD 189733b.
In the text

Fig. 11 Vertical distributions of the most abundant atmospheric constituents, after H2 and He, at four longitudes: substellar point (0°), evening limb (+ 90°), antistellar point (± 180°), and morning limb (−90°) in the atmosphere of HD 209458b. We show the mole fractions calculated by the pseudo two-dimensional chemical model (solid lines), by a model that neglects vertical mixing (horizontal transport case; dashed-dotted lines), by a one-dimensional vertical model that neglects horizontal transport (vertical mixing case; dashed lines), and by local thermochemical equilibrium (dotted lines). Photochemistry is taken into account in all cases but the last.

In the text

	Fig. 12 Same as Fig. 11, but for HD 189733b.
In the text

	Fig. 13 Vertical cuts of the abundance distributions of some of the most abundant molecules in HD 209458b’s atmosphere at longitudes spanning the 0–360° range, as calculated (from top to bottom) with the pseudo two-dimensional model, in the horizontal transport and vertical mixing cases, and under local chemical equilibrium.
In the text

	Fig. 14 Same as Fig. 13, but for HD 189733b.
In the text

Fig. 15 Effect of eddy coefficient profile and 1D/2D character of the model on the vertical abundance profiles of some of the most abundant molecules in HD 209458b. We show abundances at the substellar point and at the evening and morning limbs, as given by the pseudo 2D model using the nominal eddy coefficient profile (solid lines), by the pseudo 2D model using the Moses et al. 2011’s eddy profile (dashed lines), and by a one-dimensional vertical model using the Moses et al. 2011’s eddy profile (dotted lines).

In the text

	Fig. 16 Same as Fig. 15, but for HD 189733b.
In the text

Fig. 17 Pressure level probed by transmission and emission spectra as a function of wavelength for HD 209458b and HD 189733b. Dashed lines correspond to a model with a mean vertical profile averaged over the terminator and show the pressure level at which the tangential optical depth equals 2/3, which is in fact a transmission spectrum expressed in terms of atmospheric pressure instead of planetary radius. Solid lines correspond to a model with a mean vertical profile averaged over the dayside and indicate the pressure level at which the optical depth in the vertical outward direction equals 2/3, which is an approximate location of the region from where most of the planetary emission arises, the regions below being opaque and those above being translucent.

In the text

Fig. 18 Calculated emission spectra for HD 209458b (left panel) and HD 189733b (right panel) at 4 different phases in which the planet faces the observer the day and night sides, and the two sides in between, centered on the evening and morning limbs. Spectra have been smoothed to a resolving power R = 300. The adopted vertical profiles of temperature and mixing ratios are projected area-weighted averages over the corresponding emitting hemisphere, with the distribution of molecules given by our nominal pseudo two-dimensional chemical model. The planetary flux is shown relative to that of the star, for which the Kurucz synthetic spectrum and stellar radius given in Sect. 2 are adopted. Gray dashed lines correspond to planetary blackbody temperatures of 2000, 1500, and 1000 K (from top to bottom in HD 209458b’s panel) and of 1000 K (in HD 189733b’s panel). The inset in HD 209458b’s panel compares the emission spectrum around 4.3 μm as calculated for the dayside and as computed using the mean dayside temperature profile and the mean chemical composition of the nightside. The two spectra are nearly identical except for a slight difference at 4.3 μm due to CO2.

In the text

Fig. 19 Calculated transmission spectra for the evening and morning limbs of HD 209458b (left panel) and HD 189733b (right panel), where the vertical structure is obtained by averaging the temperature over each limb and adopting the abundance profiles at each limb from the nominal pseudo two-dimensional chemical model. Spectra have been smoothed to a resolving power R = 300. The transit depth is simply calculated as (Rp(λ)/R∗)2, where Rp(λ) is the calculated radius of the planet as a function of wavelength and R∗ is the stellar radius (see values in Sect. 2). The absolute scale of the transmission spectrum is set by our choice of assigning the value of the planetary radius given in Sect. 2 to the 1 bar pressure level. Since no attempt has been made to reproduce the absolute scale indicated by primary transit observations, calculated transit depths are somewhat higher than given by observations. The insets in both panels compare the transmission spectrum around 4.3 μm as calculated for the evening limb and as computed using the mean temperature profile of the evening limb and the chemical composition corresponding to the morning limb. The most important differences between both spectra occur around 4.3 and 15 μm, due to CO2.

In the text

Fig. 20 Contributions of the different sources of opacity to the transmission spectra of HD 209458b (left panel) and HD 189733b (right panel) smoothed to a resolving power R = 300. The H2-He continuum, whose contribution is similar in shape to the H2-H2 continuum but lower because of the lower abundance of He with respect to H2, is not shown. The spectra have been calculated using the temperature from GCM simulations and the chemical composition from the nominal pseudo two-dimensional chemical model averaged over the whole terminator. Each line represents the transmission spectrum that results from a model in which only the opacity provided by each source is taken into account.

In the text