A nongrey analytical model for irradiated atmospheres^{⋆}
I. Derivation
^{1}
Laboratoire Lagrange, UMR7293, Université de Nice SophiaAntipolis, CNRS,
Observatoire de la Côte d’Azur,
06300
Nice,
France
email:
vivien.parmentier@oca.eu
^{2}
Department of Astronomy and Astrophysics, University of
California, Santa
Cruz, CA
95064,
USA
Received:
22
July
2013
Accepted:
25
November
2013
Context. Semigrey atmospheric models (with one opacity for the visible and one opacity for the infrared) are useful for understanding the global structure of irradiated atmospheres, their dynamics, and the interior structure and evolution of planets, brown dwarfs, and stars. When compared to direct numerical radiative transfer calculations for irradiated exoplanets, however, these models systematically overestimate the temperatures at low optical depths, independently of the opacity parameters.
Aims. We investigate why semigrey models fail at low optical depths and provide a more accurate approximation to the atmospheric structure by accounting for the variable opacity in the infrared.
Methods. Using the Eddington approximation, we derive an analytical model to account for lines and/or bands in the infrared. Four parameters (instead of two for the semigrey models) are used: a visible opacity (κ_{v}), two infrared opacities, (κ_{1} and κ_{2}), and β (the fraction of the energy in the beam with opacities κ_{1}). We consider that the atmosphere receives an incident irradiation in the visible with an effective temperature T_{irr} and at an angle μ_{∗}, and that it is heated from below with an effective temperature T_{int}.
Results. Our nongrey, irradiated line model is found to provide a range of temperatures that is consistent with that obtained by numerical calculations. We find that if the stellar flux is absorbed at optical depth larger than τ_{lim} = (κ_{R}/κ_{1}κ_{2})(κ_{R}κ_{P}/3)^{1/2}, it is mainly transported by the channel of lowest opacity whereas if it is absorbed at τ ≳ τ_{lim} it is mainly transported by the channel of highest opacity, independently of the spectral width of those channels. For low values of β (expected when lines are dominant), we find that the nongrey effects significantly cool the upper atmosphere. However, for β ≳ 1/2 (appropriate in the presence of bands with a wavelengthdependence smaller than or comparable to the width of the Planck function), we find that the temperature structure is affected down to infrared optical depths unity and deeper as a result of the socalled blanketing effect.
Conclusions. The expressions that we derive can be used to provide a proper functional form for algorithms that invert the atmospheric properties from spectral information. Because a full atmospheric structure can be calculated directly, these expressions should be useful for simulations of the dynamics of these atmospheres and of the thermal evolution of the planets. Finally, they should be used to test full radiative transfer models and to improve their convergence.
Key words: radiative transfer / planets and satellites: atmospheres / stars: atmospheres / planetary systems
A FORTRAN implementation of the analytical model is available at the CDS via anonymous ftp to cdsarc.ustrasbg.fr (130.79.128.5) or via http://cdsarc.ustrasbg.fr/vizbin/qcat?J/A+A/562/A133
© ESO, 2014
1. Introduction
The discovery of numerous starplanet systems and the possibility of characterizing the planets’ atmospheric properties has led to a great many publications using radiative transfer calculations, often taken “offtheshelf” from numerical models. Given the infinite amount of possible compositions for mostly unknown exoplanetary atmospheres, it is highly valuable to be able to perform very fast calculations and also to understand what determines the thermal structure of an irradiated atmosphere.
Analytical radiative transfer solutions for atmospheres have been calculated with a variety of assumptions and in different contexts (e.g. Eddington 1916; Chandrasekhar 1935, 1960; King 1956; Matsui & Abe 1986; Weaver & Ramanathan 1995; Pujol & North 2003; Chevallier et al. 2007; Shaviv et al. 2011). However, the discovery of superEarths, giant exoplanets, brown dwarfs, and lowmass stars close to a source of intense radiation has prompted the need for solutions that account for both an outside and an inside radiation field, and properly link low and high optical depths levels. Hubeny et al. (2003), Rutily et al. (2008), Hansen (2008), Guillot (2010), Robinson & Catling (2012), and Heng et al. (2012) provide these solutions in the framework of a semigrey model, with one opacity for the incoming irradiation (generally mostly at visible wavelengths), and one opacity for the thermal radiation field (generally mostly at infrared wavelengths). These approximations have been used in hydrodynamical models of planetary atmospheres (e.g. Heng et al. 2011; Rauscher & Menou 2013), planetary evolution models (e.g. MillerRicci & Fortney 2010; Guillot & Havel 2011; Budaj et al. 2012), planet synthesis models (Mordasini et al. 2012a,b), retrieval methods (Line et al. 2012), and a variety of other applications.
Fig. 1 Optical depth vs. atmospheric temperature in units of the effective temperature. A numerical solution obtained from Fortney et al. (2008) (thick black line) is compared to the semigrey analytical solutions of Guillot (2010) for values of the greenhouse factor ranging from 0.01 to 100 (black to red lines). Low values of γ_{v} are redder. We used , T_{irr} = 1250K corresponding to the dayside average profile of a planet at 0.05 au from a sunlike star, and accounting for the albedo obtained in the numerical model. The internal temperature is T_{int} = 125 K and gravity 25ms^{2}. The effective temperature of the studied slice of atmosphere^{1} is defined by . For the numerical solution, the relation between pressure and optical depth was calculated using Rosseland mean opacities; TiO and VO opacities were not included. 

Open with DEXTER 
As shown in Fig. 1 for an atmosphere irradiated from above with a flux and heated from below with a flux , while semigrey models provide solutions that are wellbehaved when compared to full numerical solutions at optical depths larger than about unity, the temperatures at lowoptical depths appear to be systematically hotter than in the numerical solutions. Most importantly, this occurs regardless of the choice of the two parameters of the problem, i.e. the thermal (infrared) opacity κ_{th} and the ratio of the visible to infrared opacity γ_{v} ≡ κ_{v}/κ_{th}. For hot Jupiters, as in the example in Fig. 1, the real temperature profiles at low optical depths can be several hundreds of Kelvins cooler than predicted by the semigrey solutions.
The levels probed both by transit spectroscopy and by the observations of secondary eclipses of exoplanets often correspond to lowoptical depth levels (e.g. Burrows et al. 2007; Fortney et al. 2008; Showman et al. 2009), i.e. where semigrey models seem to systematically overestimate the temperatures. Furthermore, because the problem persists regardless of the main parameters, this implies that the functional form of the semigrey solutions is probably not appropriate for inversion models. Nongrey effects are known to facilitate the cooling of the upper atmosphere (see Pierrehumbert 2010, for a qualitative explanation) so they must be included. This is the purpose of the present paper.
We first describe previous analytical methods used to solve the radiative transfer problem analytically. In Sect. 3, we then derive an analytical nongrey line model and compare it to previous models in Sect. 4. In Sect. 5 we study the role of nongrey effects in shaping the atmospheric thermal structure. Eventually we apply our model to the structure of irradiated giant planets in Sect. 6. We note here that while we focus the discussion on exoplanets, we believe that this model is applicable to a much wider variety of problems, as long as an atmosphere is irradiated both from above and below. Our method can also be used to solve the radiative transfer equations in other geometries, such as the thermal structure of protoplanetary disk. We provide our conclusions in Sect. 7.
2. Assumptions and previous analytical models
2.1. Setting
2.1.1. The equation of radiative transfer
Following Guillot (2010), we will consider the problem of a planeparallel atmosphere in local thermodynamic equilibrium that receives from above a collimated flux at an angle θ_{∗} = cos^{1}(μ_{∗}) from the vertical, and from below an isotropic flux . The total energy budget of the modelled atmosphere is then set by , which defines the effective temperature in this paper^{1}.The irradiation and intrinsic fluxes are generally characterized by very different wavelengths. Although this is not required in the solution that we propose, it is convenient to think of them as being emitted preferentially in the visible and in the infrared, respectively. Scattering processes can influence both the thermal and visible radiation. As shown by Heng et al. (2012), the solution including symmetrical scattering for the incoming radiation with the Eddington approximation is equivalent to the one obtained without scattering when the irradiation flux is reduced by a factor of (1 − A) and the visible opacity is reduced by a factor 1/, where A is the Bond albedo of the planet and ξ is the ratio of the absorption to the extinction opacity (see also Meador & Weaver 1980, for a review of the different two stream methods including scattering). In this paper, the irradiation flux and the visible opacity are treated as parameters, thus we implicitly take into account symmetrical scattering of the incoming radiation (e.g. Rayleigh scattering). Scattering of the thermal radiation, however, is neglected.
In order to solve the radiative transfer problem for a plane parallel atmosphere in local thermodynamic equilibrium, one has to solve the following equation for all frequencies ν and all directions μ (Chandrasekhar 1960), (1)where I_{μν} is the specific intensity at the wavelength ν propagating with an angle θ = cos^{1}(μ) with the vertical, κ_{ν} is the opacity at a given wavelength, B_{ν} is the Planck function, and dm = ρdz is the mass increment along the path of the radiation. As usual, T, ρ, and z are the atmospheric temperature, density and height, respectively. The main difficulty in solving Eq. (1)lies in its triple dependence on μ, ν, and T and its additional dependence on m. An analytical solution requires simplifications of the opacities and of the dependence of the radiation intensity on angle.
2.1.2. Opacities and optical depth
The need for simplification implies that mean opacities must be used. The most common one is the Rosseland mean, defined as (2)When at all wavelengths the mean free path of photons is small compared to the scale height of the atmosphere, the radiative gradient obeys its welldefined diffusion limit and (unless convection sets in) the temperature gradient become that obtained from a grey atmosphere in which the opacity is set to the Rosseland mean (Mihalas & Mihalas 1984, p. 350). We hence define the optical depth τ on the basis of the Rosseland mean opacity, such that, along the vertical direction (3)Assuming hydrostatic equilibrium, the relation between pressure and optical depth can be found by integrating Eq. (3): (4)The optical depth thus becomes the natural variable to account for the dependence on depth in the radiative transfer problem. For any strictly positive Rosseland mean opacities, Eq. (4)is a bijection relating pressure and optical depth. Thus, solution of the radiative transfer equations in terms of optical depth can be converted to a solution in term of pressure for any functional form of the Rosseland mean opacities.
The second mean opacity that is traditionally used for radiative transfer is the socalled Planck mean: (5)We use the ratio of the Planck and Rosseland mean opacities to quantify the nongreyness of the atmosphere: (6)While the value of the Rosseland mean opacity is dominated by the lowest values of the opacity function κ_{ν}, the Planck mean opacity is dominated by its highest values. Thus, it can be shown that γ_{P} = 1 for a grey atmosphere and γ_{P} > 1 for a nongrey atmosphere (King 1956).
In irradiated atmospheres, a collimated flux coming from the star is absorbed at different atmospheric levels. We name κ_{v} the opacity relevant to the absorption of the stellar flux. As will be shown in Sect. 3.8, the absorption of the visible flux appears linearly in the radiative transfer equations. Thus, a solution can be found using multiple visible opacity bands κ_{v1}, κ_{v2}, etc.
We further define the ratio of the visible opacity to the mean (Rosseland) thermal opacity: (7)In order to solve the radiative transfer problem analytically, we suppose that γ_{v} is constant with optical depth. Once γ_{v} is chosen, we can solve the equations for the visible radiation independently from the final thermal structure of the atmosphere. Of course, purely grey models are such that γ_{v} = 1.
2.1.3. The picketfence model
Fig. 2 Simplified thermal opacities for the picketfence model. β = δν/Δν is the equivalent bandwidth (see text). 

Open with DEXTER 
It is important to note at this point that two sets of opacities with different wavelength dependences may have the same Rosseland and Planck means. We must constrain the problem further, and to this intent, we now consider the simplest possible line model, known as the picketfence model (Mihalas 1978), where the thermal opacities can take two different values κ_{1} and κ_{2} (see Fig. 2) such that (8)We define an equivalent bandwidth by: (9)The characteristic width of the Planck function can be defined as . When choosing Δν ≪ Δν_{P}, the Planck function can be considered constant over Δν and we get β = δν/Δν. The Planck and Rosseland mean opacities then become (see Eqs. (2)and (5)) We also define the following ratios: Following Chandrasekhar (1935), we also define a limit optical depth (15)The role of τ_{lim} in shaping the final temperature profile is discussed in Sect. 5 and its variations with R and β is pictured in Fig. 8.
2.2. The method of discrete ordinates for the nonirradiated problem
2.2.1. The grey case
An approximate method to solve Eq. (1) including the angular dependency has been developed by Chandrasekhar (1960) in the case of a nonirradiated atmosphere (T_{irr} = 0). The idea is to replace the integrals over angle in Eq. (1) by a Gaussian sum over μ. It can then be solved to an arbitrary precision by increasing the number of terms in the sum. The boundary condition at the top of the atmosphere is simply given as I_{μ < 0}(0) = 0. The expansion to the fourth term yields the temperature profile (16)with Q = 0.706920, L_{1} = −0.083921, L_{2} = −0.036187, L_{3} = −0.009461, k_{1} = 1.103188, k_{2} = 1.591778, and k_{3} = 4.45808 (Chandrasekhar 1960, Table VIII)^{2}. One of the important results from this formalism is that the skin temperature of the planet (defined as the temperature at zero optical depth) is independent of the order of expansion and therefore corresponds to the exact value (17)This expression is exact only in the limit of a grey, nonirradiated atmosphere.
2.2.2. The nongrey case
Chandrasekhar (1960) also developed a perturbation method in order to include nongrey thermal opacities. This method was improved by Krook (1963). However, these perturbation methods either work for small departures from the grey opacities or involve a fastidious iterative procedure (e.g. Unno & Yamashita 1960; Avrett & Krook 1963) and are no longer fully analytical. However, considering that the variations in the opacities are small compared to the variations of the Planck function, analytical solutions can be found for an arbitrarily large departure from the grey opacities. Noting the similar role of μ and κ_{ν} in Eq. (1), King (1956), following Münch (1946), used the method of discrete ordinates in order to turn the integrals over frequency into Gaussian sums. For the picketfence model defined in Sect. 2.1.3 and the second approximation for the angular dependency, King’s method leads to the following temperature profile (18)As in the grey case, the method of discrete ordinates leads to an exact relation for the skin temperature, whatever the dependency of κ_{ν} on frequency (but no dependence on pressure or temperature): (19)In the grey limit, γ_{P} = 1 and we recover Eq. (17). Otherwise, γ_{P} > 1, implying that for a nonirradiated atmosphere, nongrey effects always tend to lower the atmospheric skin temperature.
2.3. Moment equation method
2.3.1. Equations for the momentum of the radiation intensity
A simpler way to solve the radiative transfer equation has been carried out by Eddington (1916). The idea is to solve the equation using the different momentum of the intensity defined as (20)Then, integrating over μ Eq. (1)and μ times Eq. (1)one gets the momentum equations Assuming the atmosphere to be in radiative equilibrium, we can write (23)For a grey atmosphere (κ_{ν} = κ_{R} ∀ν), Eqs. (21)−(23) can be integrated over frequency, leading to an equation for J, H, K, and B, the frequencyintegrated versions of J_{ν}, H_{ν}, K_{ν}, and B_{ν}. The radiative equilibrium equation becomes (24)The frequency integrated versions of Eqs. (21)−(23) are a set of three equations with four unknowns. The system is not closed because by integrating Eq. (1) over all angles we have lost the information on the angular dependency of the irradiation. A closure relationship that contains this angular dependency is therefore needed. A common closure relationship, known as the Eddington approximation is (25)This relationship is exact in two very different cases: when the radiation field is isotropic (I_{μ} independent of μ), and in the twostream approximation ( and ). Although this seems to be a very restrictive approximation, it is relevant for the deep layers of the atmosphere because of the quasiisotropy of the radiation field there. It is also good for the top of the atmosphere, where the flux comes mainly from the τ ≈ 1 layer. Indeed, the exact solution gives a ratio J/K that differs by no more than 20% from the 1/3 ratio over the whole atmosphere and leads to a temperature profile that is correct to 4% in the grey case (see the plain blue line in Fig. 4).
2.3.2. Top boundary condition
Although in the method of discrete ordinates the boundary condition at the top of the atmosphere is intuitive, in the momentum equations method it is less obvious and different choices have been made by different authors. Usually, the expression for J is known and some integration constants need to be found. Two equations are needed, one for J(0) and one for H(0). Four possibilities are widely used in the literature from which one has to choose two:

1.
The radiative equilibrium equation that relates the emergent flux at the top of the atmosphere to the internal flux from the planet and the incident flux from the star.

2.
An adhoc relation between H(0) and J(0) at τ = 0: H(0) = f_{H}J(0), where f_{H} is often called the second Eddington coefficient.

3.
A calculation of H(0) from the second moment equation (Eq. (22)) and the Eddington approximation.

4.
A calculation of H(0) from the integration of the source function through the entire atmosphere, known as the Milne equation (Mihalas & Mihalas 1984, p. 347): .
For grey and semigrey models, the first condition is natural and so it was used by Hansen (2008) and Guillot (2010). For the other part of the top boundary condition, Guillot (2010) chose to use the second and Hansen (2008) the fourth condition (see Appendix A of Guillot 2010 for a comparison of the two expressions).
In the case of a nongrey model, the first condition cannot be implemented (at least directly) because it is a constraint on the total thermal flux, but it provides no information on how the thermal flux is split between the opacity bands that are considered. Chandrasekhar (1935) therefore uses conditions 2 and 3 in each of the opacity bands for his nongrey, nonirradiated model. He also notes that using condition 4 instead of condition 3 should yield better results, but it leads to more complex expressions. In this work, because an accurate treatment of the flux is needed for the nongrey irradiated model, we will use conditions 2 and 4 in each of the opacity bands. All these models are discussed in the next sections and summarized in Table 1.
2.3.3. Nonirradiated, grey case
In this section we consider the case T_{irr} = 0. Under the grey approximation, using the conditions 1 and 4, the temperature profile is given by (Mihalas & Mihalas 1984, p. 357) (26)which leads to the same solution as assuming condition 2 with . The skin temperature is then (27)which differs from the exact solution (Eq. (17)) by a factor of /2. Assuming / is thus tempting, as it leads to the correct skin temperature. Unfortunately, it leads to a temperature profile which is less accurate around τ ≈ 1.
2.3.4. Nonirradiated nongrey picketfence model
Chandrasekhar (1935) provides solutions to the moment equations for the picketfence model presented in Sect. 2.1.3. He assumes that the relation with , valid in the grey case under the Eddington approximation (see Eq. (27)), holds for the two thermal channels separately. Using this condition together with condition 3, he obtains a temperature profile (28)and an equation for the skin temperature (29)As expected, this equation reduces to Eq. (27) in the limit γ_{P} = 1. For high values of γ_{P}, this relation differs by a factor of 4/3 from the exact one derived with the method of discrete ordinates [19]. As happens for the grey case, using / would lead to the exact solution for the skin temperature, but at the expense of the accuracy of the profile at deeper levels. Again, we note that, in the nonirradiated case, the temperature at the top of the atmosphere is determined by a single parameter, γ_{P}, representing the nongreyness of the atmosphere. A comparison of the different expressions for the skin temperature is provided in Sect. 4.
2.3.5. Irradiated semigrey model
In the case of irradiated atmospheres, the presence of an incoming collimated flux at the top of the atmosphere breaks the angular symmetry of the equations. The radiative transfer problem thus can no longer be solved analytically (at least not in a simple way) through the discrete ordinates technique. The momentum method is thus required.
To solve the problem, the radiation field is split into two parts: the incoming, collimated radiation field on one hand, the thermal radiation field on the other. The radiative equilibrium equation (Eq. (23)) links the two streams, as can be seen in Sect. 3.1 (see also Hansen 2008; Guillot 2010; Robinson & Catling 2012). As mentioned in Sect. 2.1.2, when the incident radiation is at a much shorter wavelength than the thermal emission of the atmosphere, the two streams correspond to different characteristic wavelengths and may often be labelled as visible and infrared. This is not a requirement, however: the solutions apply if the radiation field corresponds to other wavelengths or if they overlap.
As discussed previously (Sect. 2.3.2), the boundary condition at the top of the model can be chosen in several ways. When using condition 2, Guillot (2010) lets the value of be either 1/2 or . Those values are based on the nonirradiated case: is the value that arises from the calculation of the angle dependence between H(0) and J(0) in the isotropic case, but provides a skin temperature that agrees with the exact value. The two solutions differ by ≈3% at most (see Fig. 5 hereafter), and choosing one over another is not crucial. In any case, for an easier comparison, we provide here the solution of Guillot (2010) for (30)where μ_{∗} is the cosine of the angle of the incident radiation. The skin temperature is (31)For γ_{v} → 0, the incident radiation is absorbed in the deep layers of the atmosphere and the skin temperature converges to the skin temperature of a grey model with an effective temperature . The semigrey model depends only on the parameter γ_{v}.
As discussed in the introduction, the semigrey model predicts minimum temperatures that are generally higher than the numerical solutions for irradiated exoplanets, independent of the choice of γ_{v} (see Fig. 1). In fact, similarly to the skin temperature, the minimum temperature of a semigrey atmosphere, shown in Fig. 3, depends only on the values of T_{eff, μ∗} and γ_{v}. It is the lowest and equal to T_{eff, μ∗}/2^{1/4} both in the γ_{v} → 0 and γ_{v} → ∞ limits. This lower bound for the semigrey temperature profile is hotter than the one obtained by numerical calculations taking into account the full set of opacities. The discrepancy is much larger than the variations resulting from the approximation of the momentum method. Clearly, nongrey effects must be invoked to explain the low temperatures obtained by numerical models at low optical depths.
Fig. 3 Minimum temperature of the semigrey model in terms of the effective temperature as a function of γ_{v}/μ_{∗}. 

Open with DEXTER 
Summary of the different models compared in this paper.
3. An analytical irradiated nongrey picketfence model
3.1. Equations
We now derive the equations for an irradiated atmosphere in local thermodynamic equilibrium with infrared line opacities as described in Sect. 2.1.2. Thus, our model contains three different opacities: κ_{1} and κ_{2} for the thermal radiation and κ_{v} for the incoming radiation of the star. As explained before, the difference between the thermal and the visible channel depends on the angular dependency of the radiation and not on the frequency. Although the method of discrete ordinates is shown to lead to more exact results, it is difficult to adapt to the irradiated case. Therefore, following Chandrasekhar (1935) and Guillot (2010), we solve the radiative transfer equations using the momentum equations. Integrating Eqs. (21) and (22) over each thermal band we obtain
where the subscript indicates the integrated quantities over the given thermal band. Thus, for a quantity X_{ν} we have:
The Planck function is considered constant over each bin of frequency Δν and therefore B_{1} = βB and B_{2} = (1 − β)B.
We now assume that the Eddington approximation is valid in the two bands separately: J_{(1,2)} = 3K_{(1,2)}. Equations (32)and (33)can be combined into and the radiative equilibrium equation becomes (38)where the quantities with subscript v are the momentum of the incident stellar radiation. Assuming that the incoming stellar radiation arrives as a collimated flux and hit the top of the atmosphere with an angle θ_{∗} = cos^{1}μ_{∗} they can are given by (39)with the total incident intensity.
The absorption of the stellar irradiation can be treated separately from the thermal radiation and J_{v} is given by Eq. (13) of Guillot (2010), (40)where we have simplified the notation by introducing the parameter .
Equations (36) to (38) are a set of three coupled equations with three unknowns J_{1}, J_{2}, and B. In order to decouple these equations we define two new variables: (41)Conversely, we can come back to the original variables: (42)Using the combination of Eqs. (36) + (37) and (38)we get (43)The combination of Eqs. (36)+(37) yields (44)Noting that , Eq. (38) becomes (45)Equations (43)−(45) are now a set of two uncoupled differential equations and a linear equation.
3.2. Boundary conditions
To solve the differential equations we need to specify the boundary conditions. When τ → + ∞ we want to fulfill the diffusion approximation J_{ν} = B_{ν} (Mihalas & Mihalas 1984, p. 350). In our case this translates to J_{1} = βB and J_{2} = (1 − β)B. Furthermore, at these levels, the gradient of B should also obey the diffusion approximation (Mihalas & Mihalas 1984) (46)where is the thermal flux coming from the interior of the planet. Using the system of Eqs. (41) and noting that , we can derive a condition on J_{γ} and J_{γ3}: For τ → 0 we specify the geometry of the intensity by setting: (49)Furthermore, we calculate the flux at the top of the atmosphere in each band using Eq. (79.21) from Mihalas & Mihalas (1984). From the assumption of local thermodynamic equilibrium, the source function in the two bands is S_{1}(τ_{1}) = βB(τ/γ_{1}) and S_{2}(τ_{2}) = (1 − β)B(τ/γ_{2}). The upper boundary condition on the flux of the two bands thus becomes (50)
3.3. Solution
The solution of a secondorder differential equation with constant coefficient is the sum of the solutions of the homogeneous equation and a particular solution of the complete equation. Thus, solutions of Eq. (43) must be of the form: (51)Applying the boundary condition Eq. (47), we get C_{2} = 3H. For τ = 0 we obtain (52)Using Eq. (45)to eliminate B and replacing J_{γ} by its solution, Eq. (44) becomes: (53)Again, solutions of this differential equation must be the sum of the solutions of the homogeneous equation and one solution of the complete equation. The homogeneous solution must have the form (54)where C_{3} and C_{4} are constants of integration to be determined using the boundary conditions. We look for a particular solution formed by the superposition of an exponential and an affine function. The affine function must then be a solution of Eq. (53) with H_{v}(0) = 0(55)and the exponential function must be solution of Eq. (53), keeping only the exponential part on the righthand side (56)Applying the boundary condition defined by Eq. (48) to the full solution J_{γ3} = J_{γ3P1} + J_{γ3P2} + J_{γ3H}, we find C_{4} = 0. The full solution of Eq. (44)is therefore given by (57)We can get an expression for the source function by replacing J_{γ} and J_{γ3} in the radiative equilibrium equation (Eq. (45)): (58)To get the complete solution of the problem, we need to determine the two remaining integration constants C_{1} and C_{3} using the boundary condition (49). For that we need to calculate J_{1}(0), J_{2}(0), H_{1}(0) and H_{2}(0). The first two quantities can be evaluated by inserting the values of J_{γ}(0) and J_{γ3}(0) from Eqs. (51) and (57) into system (42): Noting that
we can evaluate H_{1}(0) and H_{2}(0) by inserting Eqs. (58) into (50). Then, Eq. (49) is a linear system of two equations with two unknowns. After some calculations we get the expressions for C_{1} and C_{3}where we have where we defined
3.4. Atmospheric temperature profile
Using the relations B = σT^{4}/π, , and and Eq. (58), we can derive the equation for the temperature at any optical depth (76)with
3.5. Grey limit
In the grey limit, γ_{P} → 1 (as γ_{1} and γ _{2} ) and we obtain If we also assume that we obtain and E → −1/ and the solution converges towards that of Guillot (2010) (see Eq. (30)). For other values of our model differs from the solutions of Guillot (2010, see also Hansen 2008) because of the different boundary conditions used in the two models (see Sect. 2.3.2). However, calculations show that the value of C obtained here differs from the same coefficient extracted from Eq. (30) by at most 12% and that the two solutions also converge for . As seen in Fig. 5, in the semigrey limit, and when calculating the full temperature profile, our model differs by at most 2% from the Guillot (2010) model. The difference between the various solutions must be attributed to the Eddington approximation.
3.6. Using the model
The temperature vs. optical depth profile for our irradiated picketfence model is given by Eq. (76). The profile has been derived using the Rosseland optical depth as vertical coordinate. It is therefore valid for any functional form of the Rosseland opacities. Equation (4)allows us to switch from τ to P as the vertical coordinate. Although, for convenience, this expression contains four different variables, γ_{1},γ_{2},γ_{P}, and τ_{lim}, it must be kept in mind that, besides the Rosseland mean opacity, there are only two independent variables in the problem. The variables β and R ≡ γ_{1}/γ_{2} = κ_{1}/κ_{2} are the ones to consider to control the shape of the thermal opacities. The variables γ_{P} and τ_{lim} are the ones to consider to control the profile itself. The variable γ_{P} is directly related to the skin temperature of the planet (see Sect. 4.3) whereas τ_{lim} is the optical depth at which the irradiated picketfence model differs from the semigrey model. The steps to use our model are as follow:

1)
choose the pair of variables suitable for the problem: (R, β) or (γ_{P}, τ_{lim}), for example;

2)
using Eqs. (87) to (95), calculate the values of γ_{P}, γ_{1}, γ_{2}, and τ_{lim};

3)
using Eqs (77) to (81) and Eqs. (66) to (75), calculate the coefficients A, B, C, D, and E.;

4)
using Eq. (76), calculate the temperature/optical depth profile;

5)
using Eq. (4), calculate the pressure/optical depth relationship and therefore the pressure/temperature profile.
For Rosseland opacities depending on the temperature, step 5) can be iterated until convergence. Given the apparent complexity of the solution, we provide a readytouse code^{3} in different languages that gives the temperature/optical depth profile (steps 1 to 4) or the temperature/pressure profile given a Rosseland mean opacity.
The relationship between the different variables are listed below:
3.7. About averaging
Equation (76) can thus be considered to depend on κ_{R}, γ_{P} ≡ κ_{P}/κ_{R}, and β. While κ_{R} can be considered a function of pressure and temperature (e.g. extracted from a known Rosseland opacity table) when deriving the atmospheric temperature profile, it is important to realize that the analytical solution remains valid only if γ_{P} and β are held constant. This analytical solution therefore cannot accommodate consistent Rosseland and Planck opacities as a function of depth. A solution consisting of atmospheric slices with different values of γ_{P} have been derived by Chandrasekhar (1935) for the nonirradiated case, but it becomes too complex to be handled easily.
Furthermore, the solution is provided only for one fixed direction of the incoming irradiation. When considering the case of a nonresolved planet around a star, any information acquired on its atmosphere will have been averaged over at least a fraction of its surface. Solving this problem for the particular case of Eq. (76) goes beyond the scope of the present work, but it can be approximated relatively well on the basis of the study by Guillot (2010). This work shows that given an irradiation flux at the substellar point , where T_{∗} is the star’s effective temperature, R_{∗} its radius, and D the starplanet distance, the average temperature profile of the planet will be very close to that obtained from the onedimensional solution with an average angle and an average irradiation temperature T_{irr} = (1 − A)^{1/4}f^{1/4}T_{sub}, where A is the (assumed) Bond albedo of the atmosphere and f is a correction factor, equal to 1/4 when averaging on the entire surface of the planet and equal to 1/2 when averaging on the dayside only. This corresponds to the socalled isotropic approximation. In the semigrey case, it is found to be within 2% of the actual average for a typical hotJupiter (see Fig. 2 of Guillot 2010).
For the interpretation of spectroscopic and photometric data of secondary eclipses, the dayside average is often used (f = 1/2). For the calculation of evolution models, the global average is the correct physical quantity to be used when the composition and opacity variations in latitude and longitude are not precisely known (see Guillot 2010). In that case, f = 1/4 which is equivalent to setting the irradiation temperature equal to the usual equilibrium temperature defined as T_{eq} ≡ T_{∗}(R_{∗}/2D)^{1/2} (Saumon et al. 1996). Obviously detailed interpretations must use an approach mixing threedimensional dynamical and radiative transfer models (see Guillot 2010; Heng et al. 2012; Showman et al. 2009).
3.8. Adding several bands in the visible
Although, for the simplicity of the derivation, our model used only one spectral band in the visible channel, it can be easily extended to n visible bands. The most important point is that our equations, and in particular Eq. (43), are linear in the visible. Thus, the equations can be solved for any linear combination of visible bands. In that case the first momentum of the visible intensity (see Eq. (40)) writes (97)where β_{vi} is the relative spectral extent of the ith band and γ_{vi} = κ_{vi}/κ_{R} with κ_{vi} the opacity in the ith visible band. Equation (76)then becomes (98)where C_{i}, D_{i}, and E_{i} are the coefficients C, D, and E given by Eqs. (79)to (81)where have been replaced by .
4. Comparisons
4.1. Comparison of nonirradiated solutions
Figure 4 shows a comparison between our results and the solutions of King (1955) and Chandrasekhar (1935). The solutions are extremely close, the temperatures being always less than a few percentage points of each other. Our solution is almost identical to that of Chandrasekhar (1935), a consequence of using the Eddington approximation and similar boundary conditions. The difference of these with the exact solution from King (1955) can be attributed to the Eddington approximation.
The nongrey effects lead to colder temperatures at small optical depths. When β is close to unity, a blanketing effect leads to a heating of the deeper layers too. All solutions have the correct behaviour (see also Sect. 5).
Fig. 4 Comparison of the nonirradiated solutions of the radiative transfer problem within the socalled picketfence model approximation (see text). The left panel shows temperature (in T_{eff} units) vs. optical depth. The right panel shows the relative temperature difference between our model and other works. The models shown correspond to the solutions of King (1955) (blue lines), Chandrasekhar (1935) (green lines), and this work (red lines). Different models correspond to the grey case (plain), i.e. R = 1, and 2 nongrey cases: β = 0.01, R = 10^{3} (dashed) and β = 0.7, R = 10^{3} (dotted) where R ≡ κ_{1}/κ_{2}. The red and green lines are so similar that they are almost indistinguishable in the left panel. 

Open with DEXTER 
Fig. 5 Comparison between our model in the semigrey limit and Guillot (2010). We used γ_{v} = 0.25 (plain line) and γ_{v} = 10 (dashed line). For the Guillot (2010) model we show the curves for two different boundary conditions:f_{H} = 1/2 (blue) and f_{H} = 1/ (green). We used μ_{∗} = 1/. 

Open with DEXTER 
4.2. Comparison of irradiated solutions
The solutions presented in this work for the irradiated semigrey case (i.e. R ≡ κ_{1}/κ_{2} = 1) are very similar to those of Guillot (2010). As seen in Fig. 5, the solutions obtained either with f_{H} = 1/2, or have relative differences of up to 2% with those of this work. These differences are of the same kind as those arising from the use of the Eddington approximation compared to exact solutions discussed previously. They are inherent to the approximation made on the angular dependency of the radiation field and implicitly linked to the choice of the different boundary conditions discussed in Sect. 2.3.2.
4.3. Comparison of skin temperatures
Fig. 6 Skin temperature of the planet given by our irradiated picketfence model for different values of γ_{v} and in the nonirradiated case. Curves for β = 0.01 (plain lines) and β = 0.5 (dashed lines) are shown. Skin temperature from Chandrasekhar (1935) and King (1956) are also shown. For the irradiated case we used and T_{int} = 0. The γ_{v} = 0.1, the nonirradiated, and the Chandrasekhar (1935) curves are closely packed. 

Open with DEXTER 
As discussed previously, the skin temperature (temperature at the limit of zero optical depth) is an important outcome of radiative transfer and in the case of nonirradiated models, an exact solution is available. We compare our results to other analytical results in Fig. 6. In the limit of a nonirradiated planet and in the limit , our skin temperature converges to the one derived by Chandrasekhar (1935). This is an important test for the model, as for low values of γ_{v}, most of the stellar flux is absorbed in the deep layers of the planet and the model is expected to behave as a nonirradiated model with the same effective temperature. Moreover, we note that for low values of γ_{v}, the skin temperature is affected only by γ_{P} as was already claimed by King (1956) and Chandrasekhar (1935). This conclusion no longer applies for higher values of γ_{v} for which the skin temperature also depends on β. This can be seen by comparing the dotted lines and plain lines of the same colour in Fig. 6. At a given value of γ_{P}, a higher value of β corresponds to a smaller κ_{2}/κ_{1}. Depending on the value of β, the stellar irradiation can be absorbed in a region which can be optically thick to the two thermal bands, only one, or none, leading to different behaviour for the skin temperature.
5. Consequences of nongrey effects
In this section we study the physical processes that shape our nongrey temperature profile. To overcome the apparent complexity of our solution, we first derive an approximate expression for the thermal fluxes at the top of the atmosphere. We then obtain a much simpler expression for the skin and the deep temperatures. Comparing these expressions with their semigrey equivalent, we get physical insights into the processes that shape the temperature profile.
5.1. Estimation of the fluxes in the different bands
In steady state, all the energy that penetrates the atmosphere must be radiated away. Thus, the radiative equilibrium at the top of the atmosphere is of great importance to understand how the nongrey effects shape the temperature profile. In particular, whether the thermal fluxes are transported by the channel of highest opacity (channel 1) or the channel of lowest opacity (channel 2) is of particular importance.
As seen in Eq. (76), the contribution to the final temperature of the internal luminosity and of the external irradiation are independent. Thus, the thermal fluxes can be split into two independent contributions that can be studied separately: (99)Figure 7 shows which thermal band actually carries the thermal flux H_{irr}(0) out of the atmosphere (the flux F_{irr}(0) is equal to 4πH_{irr}(0)). This depends strongly on whether the stellar irradiation is absorbed in the upper or in the deep atmosphere. If it is deposited in the deep layers of the planet (i.e. γ_{v} ≪ 1), most of the flux is transported by the second thermal channel whatever the width of the second channel. Conversely, when the stellar irradiation is deposited in the upper atmosphere, most of the flux is carried by the first thermal channel whatever the width of the first channel. The tipping point, i.e. when each channel carries half of the flux, is reached when . Figure 8 shows the variations of τ_{lim} with the width and the strength of the two thermal opacity bands; τ_{lim} increases with β but decreases with R ≡ κ_{1}/κ_{2}. It always corresponds to an optical depth where the first channel is optically thick and the second is optically thin.
For high values of γ_{P} (i.e. γ_{P} > 2), we can approximate the ratio of the thermal fluxes related to the irradiation by a much simpler expression: (100)As shown in Fig. 7 this expression correctly matches the expression of the analytical model. Depending on the value of , the expression reduces to (101a) (101b)We now look for a similar expression for the thermal fluxes resulting from the internal luminosity (H_{int}). Because the internal luminosity irradiates the atmosphere from below, the resulting thermal fluxes behave similarly to the irradiated when γ_{v} → 0, thus we have (102)As γ_{P} is always greater than one and β is always lower than one, the internal luminosity is always transported by channel 2, the channel of lowest opacity.
5.2. The skin temperature
The skin temperature reveals the behaviour of the atmosphere at low optical depths. This is the part of the atmosphere probed during the transit of an exoplanet in front of its host star and is therefore of particular importance to interpret the observations. Figure 9 shows that in the irradiated case nongrey effects always tend to lower the skin temperature compared to the semigrey case. This upper atmospheric cooling is already significant (>10%) for slightly nongrey opacities (i.e. γ_{P} ≈ 2). For higher values of γ_{P} the cooling is stronger, reaching 50% for γ_{P} ≈ 10−1000. Conversely to the nonirradiated case, the skin temperature is not only a function of γ_{P} but also depends on β, i.e. not only are the mean opacities relevant, but also their actual shape. For high values of β, when the stellar irradiation is absorbed in the upper layers of the atmosphere (e.g. γ_{v} = 100) the cooling is more efficient than when the stellar irradiation is absorbed in the deep layers (e.g. γ_{v} = 0.01), whereas for low values of β the cooling is independent of γ_{v}.
The skin temperature results directly from the radiative equilibrium of the upper atmosphere. Using the boundary condition (49)in the radiative equlibrium Eq. (23)evaluated at τ = 0 we can write (103)where the skin temperature is given by . The skin temperature, depends on the values of H_{1}(0) and H_{2}(0) and thus on whether the stellar irradiation is absorbed in the deep atmosphere or in the upper atmosphere.
5.2.1. Case of deep absorption of the irradiation flux
Fig. 7 Ratio of the total flux in the two thermal bands in function of for different β and for γ_{P} = 100 given by our analytical model (plain) and by Eq. (100) (dashed). We used μ_{∗} = 1/ and T_{int} = 0. 

Open with DEXTER 
Fig. 8 Value of τ_{lim} in function of the width of the lines β and their strength κ_{1}/κ_{2}. The xaxis is in logit scale, where the function logit is defined as logit(x) = log (x/(1 − x)). 

Open with DEXTER 
When , the stellar irradiation is absorbed in the deep layers of the atmosphere, where the second thermal band, the band of lowest opacity, is optically thick. Thus, most of the flux is transported by the second thermal band and we have H_{2}(0) = H_{∞} − H_{v}(0). For high values of γ_{P}, using Eq. (101)and Eq. (102)we get which is always larger than one. Thus, although most of the flux is in the second thermal band, it is the first band, the band of highest opacity, that sets the radiative equilibrium. Neglecting the second term in Eq. (103)and calculating H_{1}(0) with Eqs. (101)and (102)we obtain: (104)Noting that for high values of γ_{P}, , and γ_{1} ≈ γ_{P}/β, the equation becomes (105)Replacing the fluxes by their equivalent temperature we get an expression for the skin temperature valid for and γ_{P} > 2: (106)When , the first term dominates and the expression differs by a factor of 1/ from the semigrey case (Eq. (31)). Because γ_{P} > 1 for nongrey opacities, the skin temperature is always smaller in the nongrey case than in the grey case, as shown in Fig. 9. Physical interpretation. When most of the irradiation is absorbed where both thermal channels are optically thick. The flux is mainly transported by the channel of lowest opacity κ_{2} but only the residual flux transported by the channel of highest opacity κ_{1} contributes to the radiative equilibrium at the top of the atmosphere. Because it represents only a small part of the total flux, the upper atmosphere does not need to radiate a lot of energy and thus the upper atmospheric temperatures are smaller than in the semigrey case. The larger the departure from the semigrey opacities, the cooler the skin temperature, without lower bounds.
5.2.2. Case of shallow absorption of the irradiation flux
When , most of the stellar irradiation is absorbed in the upper atmosphere, where only the first thermal band is optically thick. According to Eq. (101), most of the flux originating from the irradiation H_{irr}(0) is carried by the first thermal band, the band of highest opacity. Conversely, following Eq. (102), the internal luminosity is still transported by the second thermal channel, as in the γ_{v}τ_{lim} < 1 case. As γ_{1} > γ_{2}, the radiative equilibrium of the upper atmosphere is still determined by the channel of highest opacity, channel 1, and the second term of Eq. (103)can be neglected. Conversely to the case γ_{v}τ_{lim} < 1, the top boundary condition now reads H_{1}(0) ≈ H_{1, int} − H_{v}(0). Using Eq. (102)to calculate H_{1, int} and noting that for high values of γ_{P}, γ_{P} ≈ βγ_{1}, the radiative equilibrium becomes (107)Replacing the fluxes by their equivalent temperatures we get an expression for the skin temperature valid for and γ_{P} > 2(108)This relation differs from the case as the factor 1/ before the irradiation temperatures is replaced by a factor 1/β. Thus, the skin temperature no longer becomes arbitrarily low. However, for high values of γ_{v}, the second term in the parenthesis dominates and the skin temperature decreases proportionally to 1/, which is faster than in the case . As an example, in Fig. 9, for γ_{v} = 100, the skin temperature decreases much faster when γ_{P} increases for large values of β (i.e. when ). Physical interpretation. When , most of the incident irradiation is absorbed in the upper atmosphere, where the second channel is optically thin. Therefore it is mainly transported by the channel of highest opacity: channel 1. Similarly to the case , the radiative equilibrium at the top of the atmosphere is set by the channel of highest opacity, the one that carries most of the thermal flux. Therefore all the flux from the irradiation contributes to the radiative equilibrium of the upper layers and the skin temperature cannot cool as much as in the case, its lowest value being . However, for high values of and as long as , decreases faster than in the case . This confines the stratosphere due to (i.e. the atmospheric levels with a temperature inversion) around the τ = τ_{lim} level whereas it extends up to τ = 0 in the semigrey case (see Figs. 12−14 hereafter).
Fig. 9 Contours of the relative difference between the skin temperature in the nongrey model and in the semigrey model for different values of γ_{v} in function of the width of the lines and their strength. The nongrey atmosphere is 10% (resp. 50%) cooler than the semigrey atmosphere above the blue (resp. green) lines. The dashed lines are contours of γ_{P}. We used μ_{∗} = 1/. 

Open with DEXTER 
Fig. 10 Contours of the relative difference between the deep temperature in the nongrey model and in the semigrey model for different values of γ_{v} in function of the width of the lines and their strength. The nongrey atmosphere is 10% hotter (resp. cooler) than the semigrey atmosphere inside the red (resp. blue) contours. The dashed lines are contours of γ_{P}. We used μ_{∗} = 1/. 

Open with DEXTER 
5.3. The deep temperature
The temperature of the deep atmosphere is a fundamental outcome from radiative transfer models as it reveals the energy exchange between the planet and its surroundings. Therefore, it is often used as a boundary condition of planetary interior models. We define the deep temperature as (109)Thus, the temperature of the deep atmosphere can be approximated as between the τ ≈ 1 level and the radiative/convective boundary. For irradiated planets, the deep temperature corresponds to the isothermal zone around τ ≈ 1 and is close to the temperature at 10 bar often used as a boundary condition for interior models (e.g., Burrows et al. 1997; Guillot & Showman 2002). As seen in Fig. 10, the deep temperature has a complex behaviour. For low values of γ_{v}, whenever β becomes large enough, the deep temperature increases compared to the semigrey case, an effect known as the line blanketing effect in the stellar literature (e.g. Milne 1921; Chandrasekhar 1935; Hubeny & Lanz 1995). This effect is always maximum when (see Fig. 8). Conversely, for high values of (i.e. ), the deep atmosphere warms up only for high values of γ_{P} () whereas it becomes cooler than in the semigrey case for lower values of γ_{P}, a behaviour that was not spotted in previous analytical models.
The deep atmospheric temperature is directly set by the boundary condition at the top of the atmosphere. From Eq. (58), we see that when τ → ∞, (110)where C_{1} is set by the top boundary condition (49)applied on J_{γ}(0) (see Eq. (52)): (111)Similarly to the skin temperature, the deep temperature depends on H_{1}(0) and H_{2}(0), and also depends on whether the thermal flux is transported by the first or by the second thermal channel, i.e. whether γ_{v}τ_{lim} is larger or smaller than one.
5.3.1. Case of deep absorption of the irradiation flux
In the case , most of the thermal flux is transported by the second thermal channel and because γ_{1} ≫ γ_{2} we can write (112)applying the radiative equilibrium at the top of the atmosphere, and considering that most of the flux is carried by the second thermal channel, we get (113)Thus, we can calculate C_{1} and obtain (114)For high values of γ_{P}, γ_{2} ≈ (1 − β). Replacing the fluxes by their equivalent temperatures we get an expression for T_{deep} valid for and γ_{P} > 2: (115)This expression differs from the semigrey value of Guillot (2010) by a factor 1/(1 − β) multiplying the first term. Thus, when β → 1, the temperature becomes warmer than in the semigrey case, as seen for the low values of γ_{v} in Fig. 10. Physical interpretation. When , most of the flux from the star is absorbed in the deep atmosphere and is principally transported by the channel of lowest opacity (channel 2), even when the width of this channel is smaller than the width of the first thermal channel. Whenever β → 1, the width of the second channel decreases. In order to keep transporting most of the thermal flux, the flux per wavelength in the second channel must increase. This increases the temperature where the second channel is optically thick, i.e. in the deep atmosphere. This is equivalent to the line blanketing effect that has been well studied in stars (see Milne 1921; Chandrasekhar 1935; Hubeny & Lanz 1995, for example).
5.3.2. Case of shallow absorption of the irradiation flux
When most of the irradiation flux is absorbed in a region where the second thermal channel is optically thin. Thus the flux is carried by the first thermal channel and we have H_{1, irr}(0) ≫ H_{2, irr}(0). However, because γ_{2} ≪ γ_{1}, we can use Eq. (101b) to show that , which is larger than 1. Thus Eq. (112)remains valid. However, conversely to the case , the top boundary condition now reads (116)where H_{2, irr} is given by Eq. (101b) and H_{1, int} by Eq. (102). This leads to: (117)Again, for large values of γ_{P}, γ_{2} → 1 − β and replacing the fluxes by their equivalent temperatures we get an expression for T_{deep} valid for and γ_{P} > 2: (118)When , the contribution to the deep temperature of the irradiation temperature becomes inversely proportional to . As γ_{P} > 1, the deep temperature is smaller in the nongrey case than in the semigrey case. This is illustrated by the cases γ_{v} = 10 and γ_{v} = 100 in Fig. 10. When , the term in becomes very small compared to the term in 1/(1 − β) and the expression converges toward equation Eq. (115), valid for . Physical interpretation. When , the incident irradiation is absorbed in the upper atmosphere, where only the channel of highest opacity is optically thick. Thus, the channel of highest opacity κ_{1} transports all the energy and radiates it directly to space. The incident irradiation is no more transported to the deep atmosphere, leading to a cooler deep atmosphere.
Fig. 11 Ratio of the monochromatic flux in the two bands F_{ν2}/F_{ν1} = βH_{2}(0)/(1 − β)H_{1}(0) in function of the opacity ratio κ_{1}/κ_{2} for different bandwidths β and visible to infrared opacities γ_{v}. We used . 

Open with DEXTER 
5.4. Outgoing flux
During secondary eclipse observations, the flux emitted by the planet can be observed in different bands (e.g. Seager & Deming 2010). The detection of molecular species in the emission spectrum of an exoplanet depends strongly on the flux contrast between the continuum and the molecular band considered, which in turn depends on the temperature profile. Figure 11 shows the flux per wavelength emitted in the first band (F_{ν1} = 4πH_{1}(0)/β) over the flux per wavelength emitted in the second band (F_{ν2} = 4πH_{2}(0)/(1 − β)). This would be the expected contrast in the emission spectrum of the planet between the spectral features and the continuum. For a nonirradiated atmosphere and for low values of γ_{v} this is a monotonic function of the opacity ratio κ_{1}/κ_{2}. The flux in the band of lowest opacity is always bigger than the flux in the band of highest opacity, i.e. we see absorption bands. For large values of γ_{v}, whenever a strong temperature inversion happens the absorption bands turn into emission bands. Those different behaviours are captured by the simple expression (100). Note that in all cases, for large values of κ_{1}/κ_{2} we have: (119)
6. Resulting temperature profiles
No matter how strong the nongreyness of the opacities is, there is always a region, at high enough optical depth, where the nongrey solution converges toward the grey solution (see e.g. Fig. 4). The transition between a regime where the grey model is accurate to a regime where the nongrey effects are of prime importance is set by the parameter τ_{lim}. For optical depths lower than τ_{lim}, nongrey effects are always important, whereas for optical depths higher than τ_{lim}, nongrey effects are present only if γ_{v}τ_{lim} < 1 and β → 1. Three distinct situations can be observed in Fig. 8. For narrow lines (β < 0.1), τ_{lim} is always smaller than one, for larger lines or molecular bands (0.1 < β < 0.9), τ_{lim} is close to one, whereas for inverted lines (0.9 < β < 1), τ_{lim} can reach much higher values. Thus, when γ_{v} ≫ 1, few nongrey effects are expected in the deep atmosphere, contrary to the cases γ_{v} ≈ 1 and γ_{v} ≪ 1.
In the case of narrow lines (β < 0.1), only the nongrey cooling of the upper atmosphere is effective. As shown in Fig. 12, the profile remains close to the semigrey model at large optical depths. However, at low optical depths, for τ < τ_{lim}, the atmosphere can be much cooler than in the semigrey case (case R = 1). In particular, in the γ_{v} = 10 case, the nongrey cooling localizes the temperature inversion to a specific layer, contrary to the semigrey case where it extends to the top of the atmosphere. The envelope of all the profiles (shaded area) is much wider than in the semigrey case (see Fig. 1).
Fig. 12 Pressure/temperature profiles for an irradiated planet (T_{int} = T_{irr}/10 and μ_{∗} = 1/). The shaded area shows the full range of parameters 10^{3} < β < 10^{1}, 1 < R < 10^{4} and 0.01 < γ_{v} < 100. The lines are profiles obtained for β = 0.1; for R = 1 (plain lines), R = 100 (dashed lines), and R = 10^{4} (dotted lines); and for γ_{v} = 0.1 (blue), γ_{v} = 1 (green), and γ_{v} = 10 (red). 

Open with DEXTER 
In the case of large lines or molecular bands (0.1 < β < 0.9) shown in Fig. 13, both the nongrey cooling of the upper atmosphere and the blanketing effect are important. Whereas the upper atmosphere undergoes an efficient cooling, the lower atmosphere (τ > 1) can experience a significant warming via the blanketing effect. Lowering the ability of the deep atmosphere to cool down efficiently can significantly affect the evolution of the planet (Parmentier & Guillot 2011; Budaj et al. 2012; Spiegel & Burrows 2013; Rauscher & Showman 2013) and could contribute to the radius anomaly of hotJupiters (e.g. Guillot & Showman 2002; Laughlin et al. 2011). Whenever , the stellar irradiation is deposited at a level where nongrey effects lower the ability of the atmosphere to cool down efficiently. This leads to an efficient and localized warming causing a temperature inversion in the profile at , even when none is expected from the semigrey model (i.e. even when ). This happens, for example, when β ≈ 0.5 for γ_{v} = 10, when β ≈ 0.9 for γ_{v} = 1, and for β ≈ 0.99 for γ_{v} = 0.1 (see Fig. 14).
In the case of inverted lines (β > 0.9), shown in Fig. 14, both the upper temperature and the deep temperature are affected by the nongrey effects. The upper atmosphere cools significantly compared to the semigrey case. The deep atmosphere can either warm up because of the blanketing effect but, for high values of γ_{v} it can also become cooler than in the semigrey case (see the case γ_{v} = 10 and R = 100 in Fig. 14). Temperatures as cool as 0.5T_{eff, μ∗} can be reached. This is fundamentally different from the semigrey case where the deep temperature is always larger than 2^{1/4}T_{eff, μ∗} (see Fig. 3).
As β increases, τ_{lim} increases and the blanketing effect disappears. Eventually, when β → 1, the opacities, and thus the profile, become semigrey again.
In summary, our irradiated picketfence model can reach the whole temperature range span by the numerical models (see the shaded area in Figs. 12 to 14). Our model should therefore be preferred to classical semigrey models as an approximate solution for the temperature profile of irradiated atmospheres.
Fig. 13 Pressure/temperature profiles for an irradiated planet (T_{int} = T_{irr}/10 and μ_{∗} = 1/). The shaded area shows the full range of parameters 0.1 < β < 0.9, 1 < R < 10^{4} and 0.01 < γ_{v} < 100. The lines are profiles obtained for β = 0.5; for R = 1 (plain lines), R = 100 (dashed lines), and R = 10^{4} (dotted lines); and for γ_{v} = 0.1 (blue), γ_{v} = 1 (green), and γ_{v} = 10 (red). 

Open with DEXTER 
Fig. 14 Pressure/temperature profiles for an irradiated planet (T_{int} = T_{irr}/10 and μ_{∗} = 1/). The shaded area shows the full range of parameters 10^{3} < 1 − β < 10^{1}, 1 < R < 10^{4} and 0.01 < γ_{v} < 100. The lines are profiles obtained for β = 0.9; for R = 1 (plain lines), R = 100 (dashed lines), and R = 10^{4} (dotted lines); and for γ_{v} = 0.1 (blue), γ_{v} = 1 (green), and γ_{v} = 10 (red). 

Open with DEXTER 
Main quantities used in this paper.
7. Conclusion
We derived an analytic nongrey model to approximate the structure of a planeparallel irradiated planetary atmosphere. Our model includes both thermal and visible nongrey opacities. The thermal and visible opacities are in the form of a two different picketfence opacity functions. the thermal opacities are parametrized by the ratio of the visible to the infrared Rosseland mean opacities (γ_{v}), the ratio of the Planck to the Rosseland mean thermal opacities (γ_{P}), and the spectral width of the lines (β). The model is valid for any functional form of the Rosseland mean opacities, the ones obtained from an opacity table for example. However, it cannot account for both realistic Rosseland mean and Planck mean opacities. Their ratio, γ_{P} and the width of the lines, β, must be held constant through the atmosphere. Although the model is limited to two thermal opacity bands, it can take into account any number of visible opacity bands, each band adding two new parameters, the strength of the band γ_{vi} and its width β_{vi} Our model solves the inability of previous analytical models to reach temperatures as cold as predicted by the numerical calculations. For opacities dominated by strong and narrow lines (β < 0.1), nongrey opacities lead to a colder upper atmosphere, but converges toward the grey model at optical depth greater than τ_{lim} (see Fig. 8). For opacities dominated by wide lines, or molecular bands (β ≈ 0.5), nongrey opacities still allow the upper atmosphere to cool down more efficiently, but also inhibit the cooling of the deep atmosphere. In that case, a significant warming of the deep atmosphere can happen, down to optical depths much greater than τ_{lim}. This planetary blanketing effect could contribute to the radius anomaly of hot Jupiters.
Temperature inversions that were not predicted by previous analytical models occur whenever because of the interaction between the incoming stellar irradiation and the nongrey thermal opacities. These could have interesting observational consequences.
We show that the internal flux is always transported by the spectral channel of lowest opacity. Conversely, the absorbed irradiation flux is transported by the spectral channel of lowest opacity only when . For values of larger than , it is transported by the spectral channel of highest opacity.We provide simple analytical expressions for the outgoing thermal flux in the different spectral bands.
Finally, our model allows for a much greater range of temperature profiles than other analytical and semianalytical solutions of the radiative transfer equations for irradiated atmospheres. We encourage the community to use it when fast calculations of atmospheric temperature profiles are needed. Given the apparent complexity of the solution, a code is available at the CDS or at www.oca.eu/parmentier/nongrey.
In stellar physics the effective temperature is usually what we call the internal temperature. In both planetary and stellar fields, the effective temperature aims at representing the total energy budget of the atmosphere. Although in stellar physics most of the flux comes from the deep interior, this is not true in irradiated atmospheres, and is a better representation of the total energy budget of the studied slice of atmosphere. The energy budget of the whole atmosphere is therefore .
We noticed that the values of L_{1} and L_{3} in Chandrasekhar’s book were inverted and corrected this here.
Acknowledgments
This work was performed in part thanks to a joint Fulbright Fellowship to V.P. and T.G. The whole project would not have been possible without the help and support of Douglas Lin. We also acknowledge Jonathan Fortney and Mark Marley for many useful discussions, and the University of California Santa Cruz for hosting us while this work was carried out.
References
 Avrett, E. H., & Krook, M. 1963, ApJ, 137, 874 [NASA ADS] [CrossRef] (In the text)
 Budaj, J., Hubeny, I., & Burrows, A. 2012, A&A, 537, A115 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Burrows, A., Marley, M., Hubbard, W. B., et al. 1997, ApJ, 491, 856 [NASA ADS] [CrossRef] (In the text)
 Burrows, A., Hubeny, I., Budaj, J., Knutson, H. A., & Charbonneau, D. 2007, ApJ, 668, L171 [NASA ADS] [CrossRef] [MathSciNet] (In the text)
 Chandrasekhar, S. 1935, MNRAS, 96, 21 [NASA ADS] (In the text)
 Chandrasekhar, S. 1960, Radiative transfer (New york: Dover Publications) (In the text)
 Chevallier, L., Pelkowski, J., & Rutily, B. 2007, J. Quant. Spectr. Radiat. Transf., 104, 357 [NASA ADS] [CrossRef] (In the text)
 Eddington, A. S. 1916, MNRAS, 77, 16 [NASA ADS] (In the text)
 Fortney, J. J., Lodders, K., Marley, M. S., & Freedman, R. S. 2008, ApJ, 678, 1419 [NASA ADS] [CrossRef] (In the text)
 Guillot, T. 2010, A&A, 520, A27 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Guillot, T., & Havel, M. 2011, A&A, 527, A20 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Guillot, T., & Showman, A. P. 2002, A&A, 385, 156 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Hansen, B. M. S. 2008, ApJS, 179, 484 [NASA ADS] [CrossRef] (In the text)
 Heng, K., Menou, K., & Phillipps, P. J. 2011, MNRAS, 413, 2380 [NASA ADS] [CrossRef] (In the text)
 Heng, K., Hayek, W., Pont, F., & Sing, D. K. 2012, MNRAS, 420, 20 [NASA ADS] [CrossRef] (In the text)
 Hubeny, I., & Lanz, T. 1995, ApJ, 439, 875 [NASA ADS] [CrossRef] (In the text)
 Hubeny, I., Burrows, A., & Sudarsky, D. 2003, ApJ, 594, 1011 [NASA ADS] [CrossRef] (In the text)
 King, I. J. I. F. 1955, ApJ, 121, 711 [NASA ADS] [CrossRef] (In the text)
 King, J. I. F. 1956, ApJ, 124, 272 [NASA ADS] [CrossRef] (In the text)
 Krook, M. 1963, ApJ, 137, 863 [NASA ADS] [CrossRef] (In the text)
 Laughlin, G., Crismani, M., & Adams, F. C. 2011, ApJ, 729, L7 [NASA ADS] [CrossRef] (In the text)
 Line, M. R., Zhang, X., Vasisht, G., et al. 2012, ApJ, 749, 93 [NASA ADS] [CrossRef] (In the text)
 Matsui, T., & Abe, Y. 1986, Nature, 322, 526 [NASA ADS] [CrossRef] (In the text)
 Meador, W. E., & Weaver, W. R. 1980, Journal of Atmospheric Sciences, 37, 630 [NASA ADS] [CrossRef] (In the text)
 Mihalas, D. 1978, Stellar atmospheres, 2nd edn. (W.H. Freeman and Co.) (In the text)
 Mihalas, D., & Mihalas, B. W. 1984, Foundations of radiation hydrodynamics (New york: Oxford University Press) (In the text)
 MillerRicci, E., & Fortney, J. J. 2010, ApJ, 716, L74 [NASA ADS] [CrossRef] (In the text)
 Milne, E. A. 1921, MNRAS, 81, 510 [NASA ADS] (In the text)
 Mordasini, C., Alibert, Y., Georgy, C., et al. 2012a, A&A, 547, A112 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Mordasini, C., Alibert, Y., Klahr, H., & Henning, T. 2012b, A&A, 547, A111 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Münch, G. 1946, ApJ, 104, 87 [NASA ADS] [CrossRef] (In the text)
 Parmentier, V., & Guillot, T. 2011, in EPSCDPS Joint Meeting 2011, 1367 (In the text)
 Pierrehumbert, R. T. 2010, Principles of Planetary Climate (Cambridge: University Press) (In the text)
 Pujol, T., & North, G. R. 2003, Tellus A, 55, 328 [NASA ADS] [CrossRef] (In the text)
 Rauscher, E., & Menou, K. 2013, ApJ, 764, 103 [NASA ADS] [CrossRef] (In the text)
 Rauscher, E., & Showman, A. P. 2013, ApJ, submitted [arXiv:1309.7052] (In the text)
 Robinson, T. D., & Catling, D. C. 2012, ApJ, 757, 104 [NASA ADS] [CrossRef] (In the text)
 Rutily, B., Chevallier, L., Pelkowski, J., & Bergeat, J. 2008, J. Quant. Spectr. Radiat. Transf., 109, 28 [NASA ADS] [CrossRef] (In the text)
 Saumon, D., Hubbard, W. B., Burrows, A., et al. 1996, ApJ, 460, 993 [NASA ADS] [CrossRef] (In the text)
 Seager, S., & Deming, D. 2010, ARA&A, 48, 631 [NASA ADS] [CrossRef] (In the text)
 Shaviv, N. J., Shaviv, G., & Wehrse, R. 2011, Icarus, 216, 403 [NASA ADS] [CrossRef] (In the text)
 Showman, A. P., Fortney, J. J., Lian, Y., et al. 2009, ApJ, 699, 564 [NASA ADS] [CrossRef] (In the text)
 Spiegel, D. S., & Burrows, A. 2013, ApJ, 772, 76 [NASA ADS] [CrossRef] (In the text)
 Unno, W., & Yamashita, Y. 1960, PASJ, 12, 157 [NASA ADS] (In the text)
 Weaver, C. P., & Ramanathan, V. 1995, J. Geophys. Res., 100, 11585 [NASA ADS] [CrossRef] (In the text)
All Tables
All Figures
Fig. 1 Optical depth vs. atmospheric temperature in units of the effective temperature. A numerical solution obtained from Fortney et al. (2008) (thick black line) is compared to the semigrey analytical solutions of Guillot (2010) for values of the greenhouse factor ranging from 0.01 to 100 (black to red lines). Low values of γ_{v} are redder. We used , T_{irr} = 1250K corresponding to the dayside average profile of a planet at 0.05 au from a sunlike star, and accounting for the albedo obtained in the numerical model. The internal temperature is T_{int} = 125 K and gravity 25ms^{2}. The effective temperature of the studied slice of atmosphere^{1} is defined by . For the numerical solution, the relation between pressure and optical depth was calculated using Rosseland mean opacities; TiO and VO opacities were not included. 

Open with DEXTER  
In the text 
Fig. 2 Simplified thermal opacities for the picketfence model. β = δν/Δν is the equivalent bandwidth (see text). 

Open with DEXTER  
In the text 
Fig. 3 Minimum temperature of the semigrey model in terms of the effective temperature as a function of γ_{v}/μ_{∗}. 

Open with DEXTER  
In the text 
Fig. 4 Comparison of the nonirradiated solutions of the radiative transfer problem within the socalled picketfence model approximation (see text). The left panel shows temperature (in T_{eff} units) vs. optical depth. The right panel shows the relative temperature difference between our model and other works. The models shown correspond to the solutions of King (1955) (blue lines), Chandrasekhar (1935) (green lines), and this work (red lines). Different models correspond to the grey case (plain), i.e. R = 1, and 2 nongrey cases: β = 0.01, R = 10^{3} (dashed) and β = 0.7, R = 10^{3} (dotted) where R ≡ κ_{1}/κ_{2}. The red and green lines are so similar that they are almost indistinguishable in the left panel. 

Open with DEXTER  
In the text 
Fig. 5 Comparison between our model in the semigrey limit and Guillot (2010). We used γ_{v} = 0.25 (plain line) and γ_{v} = 10 (dashed line). For the Guillot (2010) model we show the curves for two different boundary conditions:f_{H} = 1/2 (blue) and f_{H} = 1/ (green). We used μ_{∗} = 1/. 

Open with DEXTER  
In the text 
Fig. 6 Skin temperature of the planet given by our irradiated picketfence model for different values of γ_{v} and in the nonirradiated case. Curves for β = 0.01 (plain lines) and β = 0.5 (dashed lines) are shown. Skin temperature from Chandrasekhar (1935) and King (1956) are also shown. For the irradiated case we used and T_{int} = 0. The γ_{v} = 0.1, the nonirradiated, and the Chandrasekhar (1935) curves are closely packed. 

Open with DEXTER  
In the text 
Fig. 7 Ratio of the total flux in the two thermal bands in function of for different β and for γ_{P} = 100 given by our analytical model (plain) and by Eq. (100) (dashed). We used μ_{∗} = 1/ and T_{int} = 0. 

Open with DEXTER  
In the text 
Fig. 8 Value of τ_{lim} in function of the width of the lines β and their strength κ_{1}/κ_{2}. The xaxis is in logit scale, where the function logit is defined as logit(x) = log (x/(1 − x)). 

Open with DEXTER  
In the text 
Fig. 9 Contours of the relative difference between the skin temperature in the nongrey model and in the semigrey model for different values of γ_{v} in function of the width of the lines and their strength. The nongrey atmosphere is 10% (resp. 50%) cooler than the semigrey atmosphere above the blue (resp. green) lines. The dashed lines are contours of γ_{P}. We used μ_{∗} = 1/. 

Open with DEXTER  
In the text 
Fig. 10 Contours of the relative difference between the deep temperature in the nongrey model and in the semigrey model for different values of γ_{v} in function of the width of the lines and their strength. The nongrey atmosphere is 10% hotter (resp. cooler) than the semigrey atmosphere inside the red (resp. blue) contours. The dashed lines are contours of γ_{P}. We used μ_{∗} = 1/. 

Open with DEXTER  
In the text 
Fig. 11 Ratio of the monochromatic flux in the two bands F_{ν2}/F_{ν1} = βH_{2}(0)/(1 − β)H_{1}(0) in function of the opacity ratio κ_{1}/κ_{2} for different bandwidths β and visible to infrared opacities γ_{v}. We used . 

Open with DEXTER  
In the text 
Fig. 12 Pressure/temperature profiles for an irradiated planet (T_{int} = T_{irr}/10 and μ_{∗} = 1/). The shaded area shows the full range of parameters 10^{3} < β < 10^{1}, 1 < R < 10^{4} and 0.01 < γ_{v} < 100. The lines are profiles obtained for β = 0.1; for R = 1 (plain lines), R = 100 (dashed lines), and R = 10^{4} (dotted lines); and for γ_{v} = 0.1 (blue), γ_{v} = 1 (green), and γ_{v} = 10 (red). 

Open with DEXTER  
In the text 
Fig. 13 Pressure/temperature profiles for an irradiated planet (T_{int} = T_{irr}/10 and μ_{∗} = 1/). The shaded area shows the full range of parameters 0.1 < β < 0.9, 1 < R < 10^{4} and 0.01 < γ_{v} < 100. The lines are profiles obtained for β = 0.5; for R = 1 (plain lines), R = 100 (dashed lines), and R = 10^{4} (dotted lines); and for γ_{v} = 0.1 (blue), γ_{v} = 1 (green), and γ_{v} = 10 (red). 

Open with DEXTER  
In the text 
Fig. 14 Pressure/temperature profiles for an irradiated planet (T_{int} = T_{irr}/10 and μ_{∗} = 1/). The shaded area shows the full range of parameters 10^{3} < 1 − β < 10^{1}, 1 < R < 10^{4} and 0.01 < γ_{v} < 100. The lines are profiles obtained for β = 0.9; for R = 1 (plain lines), R = 100 (dashed lines), and R = 10^{4} (dotted lines); and for γ_{v} = 0.1 (blue), γ_{v} = 1 (green), and γ_{v} = 10 (red). 

Open with DEXTER  
In the text 