© ESO, 2010

1. Introduction

Since the discovery of the first extrasolar planets around pulsars (Wolszczan & Frail 1992) and main-sequence stars (Mayor & Queloz 1995), more than 400 planets and planet candidates have been detected, most of which are Jupiter-mass objects. However, in recent years, detections of smaller planets have been announced, and today 21 planets with (minimum) masses below 10 Earth masses (m_⊕) are known (e.g. Rivera et al. 2005; Howard et al. 2009). Two of these planets, CoRoT-7 b (Léger et al. 2009) and GL 1214 b (Charbonneau et al. 2009), have been discovered by the transit method and produce a measurable radial velocity variation, hence their true mass and mean density could be derived. They are the first two examples of so-called “Super-Earths” (1   m_⊕ < m ≲ 10   m_⊕).

Modeling studies of the atmospheres of terrestrial (exo)planets (e.g. Kasting et al. 1993; Joshi et al. 1997; Forget & Pierrehumbert 1997; Pierrehumbert & Erlick 1998; Selsis et al. 2002; Segura et al. 2003, 2005; Grenfell et al. 2007a,b; Kitzmann et al. 2010) have aimed at characterizing the response of the atmospheric system to changes in planetary (e.g. atmospheric composition, presence of clouds) and stellar parameters (e.g. central star type, orbital distance). These studies mostly assessed the potential surface habitability of such planets and predicted their spectral appearance. Surface habitability is commonly related to the presence of liquid water on the surface, hence temperatures above 273 K, because liquid water seems to be the fundamental requirement for life as we know it on Earth.

In view of these modeling activities, the system Gliese 581 (Gl 581) with four planets (Bonfils et al. 2005; Udry et al. 2007; Mayor et al. 2009) is particularly interesting. It hosts at least three potentially low-mass, therefore possibly terrestrial planets. The one closest to the central star, Gl 581e, was announced by Mayor et al. (2009) and has a minimum mass of 1.94 m_E with an orbital distance of 0.03 AU to the star. The two outer low-mass planets, Gl 581c and d, were discovered by Udry et al. (2007). Mayor et al. (2009) refined the orbital distances and minimum masses of these two planets, obtaining 5.36 m_E (Gl 581c at a distance of 0.07 AU to the star) and 7.09 m_E (Gl 581d at a distance of 0.22 AU to the star), respectively.

Taking into account photometric (López-Morales et al. 2006; Mayor et al. 2009) and dynamical (Beust et al. 2008; Mayor et al. 2009) constraints, the inclination of the Gl 581 system most probably lies between 40° and 85°. This implies masses no higher than 1.56 times the minimum masses, so all three low-mass planets are likely to be Super-Earths.

The habitability of Gl 581c and d has been investigated by Selsis et al. (2007) and von Bloh et al. (2007a), based on the discovery data published by Udry et al. (2007). In their study, von Bloh et al. (2007a) applied a box-model of the Earth, which incorporates the carbonate-silicate cycle, requiring a surface reservoir of liquid water and a tectonically active planet. It includes the exchange of CO₂ between the mantle and crust of the planet and its atmosphere by assuming parameterizations for continental growth and spreading rate. The CO₂ partial pressure and the stellar luminosity were related to the surface temperature through a simple energy balance equation between stellar and emitted thermal fluxes (see also von Bloh et al. 2007b for more details). Selsis et al. (2007), on the other hand, used previous model results for the habitable zone (HZ) from the seminal study of Kasting et al. (1993) to estimate the boundaries of the HZ in the Gl 581 system. The results of Kasting et al. (1993) were obtained with a 1D radiative-convective model. Furthermore, Selsis et al. (2007) used model results for early Mars from Mischna et al. (2000) to illustrate the uncertainty in the outer limit of the HZ because of the (possible) presence of CO₂ clouds.

Both Selsis et al. (2007) and von Bloh et al. (2007a) concluded that the inner planet Gl 581c is unlikely to be habitable, because it is closer to the star than the inner boundary of the HZ, whereas the outer planet, Gl 581d, just might be habitable. Based on the calculations by Selsis et al. (2007), Mayor et al. (2009) concluded that Gl 581d lies in the HZ of its central star, considering that the refined orbit means it receives more than 30% more stellar energy than previously thought.

Very recently, Wordsworth et al. (2010) presented 1D modeling studies of high CO₂ atmospheres of Gl 581d, varying the CO₂ pressure and other parameters such as relative humidity and surface albedo. They found habitable surface conditions with CO₂ partial pressures of  ≥ 5 bar. In a similar approach to Wordsworth et al. (2010), we applied a different 1D radiative-convective model to the atmosphere of Gl 581d considering different appropriate planetary scenarios. By calculating atmospheric pressure, temperature, and water profiles, a reasonable range of planetary surface conditions has been investigated with regard to potential habitability.

The paper is organized as follows: Sect. 2 describes the properties of the planetary system Gl 581. The model used is described in Sect. 3. A description of the model input parameters and of the runs performed is given in Sect. 4. Results are described and discussed in Sect. 5. We give our conclusions in Sect. 6.

2. The planetary system Gl 581

2.1. Main properties of the star Gl 581

Gl 581 is a quiet M 3 star (Udry et al. 2007). Its main parameters are summarized in Table 1.

In order to perform proper atmospheric modeling studies, we need a spectrum of Gl 581. Therefore, we derived a high-resolution spectrum (Fig. 1) from two sources, i.e. an UV spectrum measured by the IUE (International Ultraviolet Explorer) satellite in 1987 and a synthetic Nextgen model spectrum (Hauschildt et al. 1999).

The model parameters of the available NextGen model spectrum used here are T_eff = 3200 K, log g = 4.5,  [ Fe / H ]  = 0.0 and R_Gl   581 = 0.38   R _⊙ . The synthetic spectrum was calculated for wavelengths from 0.334–971 μm. It was merged with the IUE spectrum from 0.185–0.334 μm.

Based on the 6.27 pc distance of Gl 581 (Butler et al. 2006), the spectrum was scaled to the orbital distance of Gl 581d (0.22 AU, Mayor et al. 2009) and binned to the spectral intervals required for the model code (see Sect. 3).

The orbit of Gl 581d could be highly eccentric (Table 2). The mean stellar flux received over an eccentric orbit is given by

(1)

where F(a) is the flux received by the planet at the distance of its semi-major axis a. As shown by Williams & Pollard (2002) with a 3D model approach, the overall climate of a planet on an eccentric orbit behaves almost as if the planet were constantly receiving this averaged flux . For Gl 581d, the application of Eq. (1) results in an increase of the stellar flux of about 8% compared to the circular case.

Assuming that the atmospheric response timescales to changes in stellar flux (which are on the order of weeks/months) are comparable to the orbital period of the planet (about teo months), was used for the model calculations. Note, though, that this approach depends on the properties of the central star and the planetary system. For larger orbital distances, for instance, flux variations caused by orbital eccentricities should be taken into account.

2.2. Properties of the planet Gl 581d

The semi-major axis is known from the measured orbital periods and the mass of the central star via Kepler’s third law. This we took from Mayor et al. (2009) to be a_Gl   581d = 0.22 AU. The best fit of the radial velocity data was obtained for a highly eccentric orbit with e_d = 0.38 (see Mayor et al. 2009).

In order to be consistent with the studies of Selsis et al. (2007) and von Bloh et al. (2007a), the stated minimum mass was assumed to be the true mass m of the planet. The planetary radius r is then taken from a theoretical mass-radius relationship for terrestrial planets (Sotin et al. 2007), finally yielding the surface gravity.

Usually, planets in the HZ of M stars (or closer to the star) are assumed to be tidally locked in a synchronous rotation, i.e. they rotate with the same rate ω as their orbital period P (e.g. Lammer et al. 2007; Scalo et al. 2007). However, as already noted by Goldreich & Peale (1968), Grießmeier et al. (2005), Levrard et al. (2007) and Correia et al. (2008), planets that show significant eccentricities – possibly like Gl 581d – are unlikely to become locked in a 1:1-resonance. Therefore, Gl 581d is not considered to rotate synchronously in this study.

An important model parameter for 1D radiative-convective models is the surface albedo. One approach used in these models is to adjust the surface albedo in the model (e.g., Kasting et al. 1993; Grenfell et al. 2007b; Goldblatt et al. 2009a) so that the models reproduce prescribed reference scenarios (e.g. modern Earth with a surface temperature of 288 K). Therefore, the values of the model surface albedos largely depend on the impact of clouds on the surface temperature. This impact, however, is unknown for Super-Earths. Hence, instead of an arbitrary choice of the surface albedo used in the tuning of 1D models, we adopted the measured surface albedo of Earth here. Thus, the unconstrained effect of clouds in Super-Earth atmospheres is explicitly excluded from our cloud-free simulations. The surface albedo value used here is A_surf = 0.13 (Earth mean, see Kitzmann et al. 2010; and Rossow & Schiffer 1999).

Table 2 summarizes the assumed planetary parameters of Gl 581d.

3. Atmospheric model

3.1. General model description

A 1D, cloud-free radiative-convective column model was used for the calculation of the atmospheric structure, i.e. the temperature and pressure profiles for different Gl 581d scenarios.

The model is originally based on the one described by Kasting et al. (1984a) and Kasting et al. (1984b). Further developments are described by e.g. Kasting (1988), Kasting (1991), Kasting et al. (1993), Mischna et al. (2000) and Pavlov et al. (2000). Additional updates of the thermal radiation scheme of the model have been introduced by Segura et al. (2003). The model version used here is based on the version of von Paris et al. (2008) where a more detailed model description is given. The current model uses the radiative transfer scheme MRAC in the IR (see the introduction of this scheme in von Paris et al. 2008). The water profile in the model is calculated based on the relative humidity distribution of Manabe & Wetherald (1967). Above the cold trap, the water profile is set to an isoprofile at the cold trap value.

The model considers N₂, H₂O, and CO₂ as atmospheric gases. Other radiative trace species might be present in the atmospheres of exoplanets (e.g. SO₂, O₃, or CH₄), which could alter the radiative budget significantly. But because the presence of these gases highly depends on the planetary scenario (e.g., outgassing, volcanism, formation, biosphere), our model atmospheres are restricted to the two most important greenhouse gases of the Earth’s atmosphere (H₂O and CO₂), using N₂ as additional background gas. This is based on the observation that in the solar system all terrestrial atmospheres (Venus, Earth, Mars, Titan) contain significant amounts of N₂.

Temperature profiles from the surface up to the mid-mesosphere are calculated by solving the equation of radiative transfer and performing convective adjustment, if necessary. The convective lapse rate is assumed to be adiabatic. Water profiles are calculated by assuming a fixed relative humidity profile (Manabe & Wetherald 1967) through the troposphere.

As a globally averaged 1D column model, the model does not incorporate dynamical processes such as winds and latitudinal energy transport. Thus potentially important 3D effects are not included. For example, the effect of atmospheric collapse and global energy re-distribution on tidally locked planets was investigated by Joshi et al. (1997) and Joshi (2003). The study of Spiegel et al. (2008) showed that slowly rotating planets feature a very efficient energy redistribution because of the weaker Coriolis force (compared to Earth). Hence a 1D model could give a reasonably accurate picture of global mean temperatures (and therefore, habitability) on planets orbiting M stars. In addition to rotation, also atmospheric optical depth and mass of a planet are important for assessing temperature gradients (Joshi et al. 1997). Therefore, 3D atmospheric modeling of Gl 581d is needed in the future to investigate its potential habitability in more detail.

3.2. Computational details

In this work, we employed an improved version of the model of von Paris et al. (2008):

• 1.

  in the present study, CO₂ condensation is included for the calculation of the adiabatic lapse rate in the convective regime according to the treatment by Kasting (1991) and Kasting et al. (1993). The saturation vapor pressure of CO₂ is taken from Ambrose (1956). CO₂ condensation is generally assumed to occur when the atmosphere is (super-)saturated with respect to CO₂. This criterion is expressed by the super saturation ratio S_sat, which should be higher than unity:

  (2)

  where p_CO2 is the partial CO₂ pressure and P_sat,CO2 the saturation vapor pressure of CO₂. The numerical value of S_sat used here is taken to be 1.34, based on measurements by Glandorf et al. (2002). This implies that the atmosphere has super-saturated with respect to CO₂ for CO₂ condensation to occur. Note that Kasting (1991) and Kasting et al. (1993) assumed S_sat = 1, i.e. the model atmosphere was assumed to be saturated;

• 2.

  the calculations now also include Rayleigh scattering by water vapor (σ_ray,H2O, in [cm^-2]) calculated via

  (3)

  where D is the depolarization ratio, r the refractivity and λ the wavelength in μm. The numerical factor 4.577 × 10^-21 is taken from Allen (1973). This work assumes D = 0.17 from Marshall & Smith (1990). The refractivity is calculated as r = 0.85·r_dryair (Edlén 1966). The refractivity of dry air is obtained from an approximate formula given by Bucholtz (1995);

• 3.

  the heat capacity of water vapor has been included in the calculations based on the so-called Shomate equation (Parks & Shomate 1940) with parameters from Chase (1998);

• 4.

  one critical issue for the IR radiative transfer scheme used here is the so-called binary species parameter η (Mlawer et al. 1997; Colaprete & Toon 2003; von Paris et al. 2008). In previous studies, the variation of absorption cross sections of different mixtures was assumed to be linear according to the relative concentrations of the major absorbers (here, H₂O and CO₂). For atmospheres where the relative concentrations of water and carbon dioxide vary by several orders of magnitude, this approach depends on the assumed reference ratio. Hence, a logarithmic η grid has been introduced where relative concentrations are allowed to vary by 15 orders of magnitude;

• 5.

  the IR H₂O continuum calculation used in von Paris et al. (2008) has been replaced by the CKD continuum formulation (Clough et al. 1989; Schreier & Böttger 2003). In this way, the continuum is no longer restricted to the atmospheric window region (between 8 and 12 μm) and incorporates contributions from H₂O self as well as foreign continua and, for the first time in the context of exoplanet 1D models, the CO₂ foreign continuum in the IR;

• 6.

  the temperature-pressure grid for the absorption coefficient calculations in the IR has been extended to T = 100 K and p = 1000 bar. This was done as described by von Paris et al. (2008), using line parameters from Rothman et al. (1995);

• 7.

  the model has been further modified to account for changing pressure levels during the calculations. Under some conditions described here, the atmospheric content of water or carbon dioxide is controlled by the surface temperature caused by evaporation or condensation. The surface pressure is thus determined by the surface temperature and must be re-calculated after each model time step to adjust to the newly calculated value of the surface temperature. The model pressure grid is therefore no longer fixed, but dynamically determined between two consecutive iterations in the model.

The IR radiative transfer updates mentioned above were tested against high-resolution line-by-line (lbl) calculations with the MIRART radiative transfer code (Schreier & Böttger 2003), which was found to compare well with other lbl codes (e.g., Melsheimer et al. 2005).

The total IR flux was within 5% of the lbl values for all atmospheric levels. Fluxes for each spectral interval of the radiative transfer code were usually within 10% of the lbl values, except for a few near-IR bands in the middle to upper atmosphere. However, these bands do not contribute strongly to the cooling rate (less than 10^-6) or the stratospheric energy budget, hence these deviations are insignificant for the calculation of the temperature profiles. Further comparative studies with the improved model are presented in Appendix A.

4. Model input for Gl 581 and performed runs

Important properties characterizing planetary atmospheres are the total surface pressure and the composition of the atmosphere. These determine the greenhouse effect and the volatile reservoir, hence the habitability of a planet. Atmospheric pressure and composition are determined by the accretion and outgassing history, whether from the interior or from impacts (e.g. Pepin 1991), atmospheric loss because of escape to space (e.g. Kulikov et al. 2007; Kulikov et al. 2006) or incorporation into planetary surface material. But these processes are not known for the Gl 581 system. Hence, in order to study the influence of atmospheric properties on the potential habitability of Gl 581d, we performed a parameter study considering different atmospheric scenarios.

We assumed that the initial water delivered to the planet via outgassing and impacts has been retained, hence a reservoir of liquid water is present on the planetary surface. This assumption was also made by Selsis et al. (2007), von Bloh et al. (2007a) and Wordsworth et al. (2010) for the Gl 581 planetary system.

The total surface pressure (1, 2, 5, 10, 20 bar) and CO₂ concentration (0.95, 0.05 and 355 ppm vmr, respectively) were varied. N₂ is assumed to act as a filling background gas. The range of surface pressures was chosen to represent scenarios adopted in the literature for early Earth and early Mars in terms of atmospheric column density. The values of CO₂ concentrations are chosen to represent modern Earth (355 ppm CO₂ concentration), early Earth (5% CO₂) and Mars or Venus (95% CO₂) conditions, i.e. representative values of terrestrial atmospheres within the solar system. Table 3 summarizes the set of runs for Gl 581d performed here.

5. Results and discussion

5.1. Temperature profiles

Figure 2 shows the temperature-pressure profiles for the low CO₂ set of runs (G1) of Table 3. For all runs the surface temperatures are far lower than the freezing point of water (273 K, indicated by the dotted vertical line in Fig. 2), which implies that liquid water is not present on the surface of the planet. Upon increasing the surface pressure from 1 bar to 20 bar, the surface temperatures increase by 4.2 K (from 203.6 K to 208.9 K).

Also shown in Fig. 2 is the equilibrium temperature T_eq of the planet. It is about 197 K, hence habitable conditions would require a greenhouse effect (GHE) of 76 K. As can be inferred from Fig. 2, the actual GHE in the model atmospheres is about 6 to 12 K, which is much less than the approximately 30 K GHE provided by the modern Earth atmosphere.

An interesting feature of the temperature profiles of the low CO₂ case is the absence of a convective troposphere for all G1 runs. These atmospheres are all in radiative equilibrium.

Figure 3 shows the temperature-pressure profiles for the medium CO₂ set of runs (G2) of Table 3. For the 20 bar run, the value of the surface temperature was 313.3 K, hence exceeded the freezing point of water. Thus, in this particular scenario Gl 581d might be habitable. The other runs with 1, 2, 5, and 10 bar though all showed surface temperatures below 273 K. Note, however, that for the 10 bar run the obtained surface temperature is 268.1 K, which is fairly close to the freezing point of water.

In contrast to the low CO₂ case, the increase of surface pressure from 1 to 20 bar has a huge effect on surface temperature, which increased by about 105 K. This is caused by an enhanced absorption of stellar radiation by CO₂ and H₂O near-IR absorption bands and a massive greenhouse effect, which results in the onset of convection for runs with higher surface pressures.

Figure 4 illustrates this effect, showing net (i.e., frequency-integrated) stellar and thermal downward (F_d) and upward (F_u) fluxes for the medium CO₂ (G2) 20 bar run. Much of the incoming stellar radiation (~70%) is absorbed by CO₂ and water in the lower atmosphere, as illustrated by the left panel in Fig. 4. The difference C_r = F_u − F_d for the thermal fluxes is the radiative cooling. A strong GHE is indicated by a low value of C_r. In the lower atmosphere of the 20 bar run, both thermal components balance each other, as can be seen in Fig. 4, meaning that the thermal radiation is efficiently trapped in the atmosphere. The value of C_r in the bottom atmosphere layer is about 5 Wm^-2, which corresponds to about 1% of the surface emission (σ_B Wm^-2, σ_B Stefan’s constant, T_surf = 313 K surface temperature). On Earth, the value of C_r is about 70 Wm^-2, about 20% of the surface emission.

The atmospheric structure of the medium CO₂ runs in Fig. 3 is very different from the low CO₂ runs in Fig. 2. Firstly, the temperature inversion in the upper atmosphere of the G2 runs at pressures below 10 mbar is much more pronounced (~30 K) than for the low CO₂ runs. This is owing to the strong absorption of stellar radiation by CO₂ and H₂O in the near-IR bands (at 2, 2.7 and 4.3 μm).

Secondly, convective tropospheres develop for the 5, 10, and 20 bar runs. These extend to an altitude of about 0.5–1 surface scale heights (about 2–5 bar) above the surface, which is comparable to the troposphere extension on present Earth. Hence, the lapse rate is much steeper in the medium CO₂ cases than in the low CO₂ cases. For example, the convective lapse rate in the 20 bar medium CO₂ case is about 22 K km^-1 near the surface (the dry adiabatic lapse rate varies linearly with gravity), whereas the radiative lapse rate for the 20 bar low CO₂ run is only about 9 K km^-1. The value of 22 K km^-1 is very close to the dry adiabatic lapse rate of 23 K km^-1 in the medium CO₂ case. Despite the high surface temperature of 313 K and a corresponding partial pressure of water of about 70 mbar, the water concentrations are only of the order on 10^-3 vmr near the surface. Hence, the lapse rate is close to the dry adiabat (see discussion in Ingersoll 1969), even if appreciable amounts of water are present in the atmosphere.

In Fig. 5, the temperature-pressure profiles for the high CO₂ set of runs (G3) are shown. Except for the low-pressure runs with 1 and 2 bar surface pressure, all scenarios showed surface temperatures above 273 K. Hence, the results suggest that CO₂-rich atmospheres may imply habitable conditions on Gl 581d.

The atmospheric structure in the high-CO₂ case is different from those in the low and medium CO₂ cases. Even the 1 and 2 bar runs now show convective tropospheres, albeit not very extended ones. More massive tropospheres develop in the 5 bar, 10 bar, and 20 bar runs owing to the onset of CO₂ condensation in the middle atmosphere, as indicated by the CO₂ condensation curve in Fig. 5. The respective tropopauses are located 3–5 surface scale heights above the planetary surface. On present Earth, this would correspond to tropopause heights of about 20–40 km compared to the roughly 10 km tropopause altitude observed today.

The tropospheres in the 5, 10, and 20 bar high CO₂ runs are divided into two regimes, an upper troposphere with CO₂ condensation (see CO₂ condensation curve in Fig. 5) and a lower troposphere with H₂O condensation. This temperature structure is comparable to atmospheric structures calculated for models of the early Mars atmosphere (Kasting 1991; Mischna et al. 2000; Colaprete & Toon 2003).

From Fig. 5 it is also evident that temperature profiles are similar to each other for the 5, 10, and 20 bar runs. The temperature profiles of these runs differ by about 2–5 K at equal pressures, as is the case for the 10 and 20 bar medium CO₂ runs (Fig. 3). This is because of two reasons. Firstly, the lower atmospheres of the 5, 10, and 20 bar high-CO₂ runs become optically thick for thermal radiation, as is observed for the medium CO₂ 10 and 20 bar runs. Secondly, an increase in surface pressure does not result in an increase in planetary albedo because Gl 581 emits only a negligible amount of radiation in the spectral range where the Rayleigh cross section is largest. Therefore, the global energy balance is no longer affected by the lower part of the model atmosphere, hence temperature profiles do not differ significantly.

A widely used assumption for the assessment of the HZ is the so-called maximum greenhouse, introduced by Kasting et al. (1993) for an Earth-like planet around the Sun. With increasing amounts of CO₂ in the atmosphere of a planet located near the outer boundary of the HZ, the GHE will become more and more saturated, i.e. the optical depth is near or higher than unity for all CO₂ bands. Then, any further increase of CO₂ will only increase Rayleigh scattering, hence increase the planetary albedo. Thus, for increasing CO₂ partial pressure, surface temperatures will show a maximum and then decrease with increasing CO₂ pressure. The same behavior was found for early Mars (Mischna et al. 2000). The different behavior for our high CO₂ runs compared to the results of Kasting et al. (1993) and Mischna et al. (2000) is illustrated in Fig. 6. They are summarized in terms of atmospheric column density and normalized surface temperature T_surf / T_eq where T_surf is the surface temperature of the planet. As is obvious from this figure, the high-CO₂ runs from this work do not feature a maximum greenhouse effect. However, a maximum greenhouse effect will probably occur at still higher pressures, hence column densities (not investigated here). The absence of a maximum greenhouse effect for the runs of this work is caused by four factors.

Firstly, Gl 581 emits much more radiation in the near- to mid-IR, and less in the visible than the Sun. Consequently, the contribution of Rayleigh scattering to the planetary albedo is much less efficient for planets around Gl 581 than around the Sun because of the λ^-4-dependence of the Rayleigh scattering cross section. Secondly, the stronger near-IR emission of Gl 581 leads to more heating by near-IR absorption bands of H₂O and CO₂. Thirdly, the simulations of Kasting et al. (1993) were done at constant surface temperatures of 273 K (hence, constant partial pressure of 6.5 mbar water), which neglects the positive feedback provided by increased water vapor at higher surface temperatures. Fourthly, because of the higher gravity of Gl 581d compared to Earth, the same column amount of CO₂ (i.e., x coordinate in Fig. 6) is reached at much higher pressures, e.g. 20 bar on Gl 581d compared to 8 bar on Earth. The pressure broadening of absorption lines then leads to an enhanced absorption in the line wings. Because the line centers are usually optically thick, a higher absorption coefficient in the line wings can significantly increase the overall absorption of radiation. It was recently suggested that this kind of behavior could have helped to warm the early Earth by invoking higher N₂ partial pressures than today (Goldblatt et al. 2009a). It is also mainly responsible for the fact that surface temperatures for the 1 bar high CO₂ case are about 85 K lower than for the 20 bar medium CO₂ case, although both atmospheres contain the same amount of CO₂. When interpreting this result, one has to bear in mind the potential uncertainties related to the radiative transfer in high-pressure atmospheres (see sensitivity study in Appendix B).

5.2. Comparison with other studies of Gl 581d

For several model scenarios (5, 10, and 20 bar high CO₂ and 20 bar medium CO₂ runs) of this work, surface temperatures were found to be above 273 K, i.e. these results imply habitable surface conditions on Gl 581d. In all other scenarios, Gl 581d was found to be uninhabitable with surface temperatures below 273 K.

The results are summarized in Fig. 7 in the considered parameter space, i.e. the surface pressure – CO₂ concentration plane.

Yet with regard to the uncertainties in the radiative transfer associated with CO₂ continuum absorption (see the sensitivity study for surface temperature in Appendix B), the modeling results regarding the habitability must of course be treated with caution.

Still, our results imply that, given plausible Venus- or early Earth-like atmospheric compositions and taking into account reasonable uncertainties in the radiative transfer formulations, Gl 581d can be classified as a habitable planet.

The study of von Bloh et al. (2007a) concluded that Gl 581d represents a habitable planet even for relatively low CO₂ partial pressures of about 4–5 bar. Taking into account the more than 30% increase in insolation for Gl 581d owing to the revision of orbital parameters by Mayor et al. (2009), the results of von Bloh et al. (2007a) imply that Gl 581d could indeed be habitable for even less dense CO₂ atmospheres. The findings of our study do not support these implications of von Bloh et al. (2007a) though.

This disagreement is partly owed to the parameterizations and empirical criteria employed. For example, the parametrization of the relation between CO₂ partial pressure, planetary albedo, and surface temperature used by von Bloh et al. (2007a) is based on work by Williams & Kasting (1997), Budyko (1982) and Chamberlain (1980). These studies were done for Earth, orbiting around the Sun, and thus did not account for different central stars with a different spectral energy distribution or more massive terrestrial planets. In contrast, the radiative-convective model used here takes these effects into account.

The study of Selsis et al. (2007) concluded that a pure CO₂ atmosphere without additional greenhouse gases is unlikely to provide habitable conditions on Gl 581d. Like in the study of von Bloh et al. (2007a), they did this for the then known orbital parameters of Gl 581d. Selsis et al. (2007) also provide a discussion of the uncertainties of the limits of the outer HZ in view of early Mars. From this discussion, Gl 581d could still be in the HZ when CO₂ clouds or additional greenhouse gases are taken into account. Note that, as already stated by Selsis et al. (2007), the use of parameterizations of the outer boundary of the HZ provided by Kasting et al. (1993) is uncertain for low-mass stars below about 3700 K effective temperature.

Based on the calculations by Selsis et al. (2007), Mayor et al. (2009) concluded that Gl 581d is a habitable planet, considering that the revised orbit means that it receives more than 30% more stellar energy than previously thought. Because according to our study the low CO₂ and most of the medium CO₂ scenarios are uninhabitable (see Figs. 2 and 3), this disagrees with their general implication. This can be motivated in more detail as follows: given that the equilibrium temperature of a planet only increases as the fourth root of stellar energy input, the increase of received stellar energy would only lead to an increase of about 10 K in equilibrium temperature (from about 185 K to 195 K). Thus, only an atmosphere providing a large greenhouse effect would be able to warm the planet above 273 K, hence result in habitable conditions (see for example Fig. 5).

The recent atmospheric modeling studies by Wordsworth et al. (2010) qualitatively agree with our simulations. Wordsworth et al. (2010) use a 1D radiative-convective model like ours. They also incorporate a correlated-k approach for the radiative transfer, similar to what what we did. However, they used a different parameterization of the continuum absorption of CO₂ than the one used here (for a sensitivity study of CO₂ continuum absorption, again see Appendix B). Their calculated surface temperatures for the high-CO₂ cloud-free cases (set G3 in Table 3) are comparable with the results obtained here (about 310 K for a 10 bar atmosphere, about 350 K for a 20 bar atmosphere).

6. Conclusions

Detailed model calculations of possible atmospheres for the low-mass extrasolar planet Gl 581d have been presented in this explorative modeling study. Using an improved 1D radiative-convective climate model, several key atmospheric parameters (e.g. surface pressure, atmospheric composition) were varied to investigate their influence on resulting surface conditions. The planetary scenarios investigated here are consistent with assumptions made in the literature regarding surface pressures and atmospheric compositions of terrestrial planets.

Gl 581d is a potentially habitable planet, because for massive CO₂ atmospheres (5 or more bar surface pressure with CO₂ concentrations of 95%, 20 bar with 5% CO₂), the surface temperatures exceeded 273 K, i.e. the freezing point of water. For these massive CO₂ atmospheres, the surface temperatures could be as high as 357 K.

The results of our high CO₂ scenarios are confirmed by the very recent model results of Wordsworth et al. (2010). But our model calculations imply that also atmospheres that are not CO₂-dominated (i.e., 5% vmr instead of 95% vmr) could result in habitable conditions for Gl 581d.

For atmospheric scenarios with less CO₂, the planet was found to be uninhabitable in our calculations.

Nevertheless, Gl 581d is the first extrasolar (potentially terrestrial) planet where habitable conditions are at least conceivable within a reasonable range of surface pressures and CO₂ concentrations.

Acknowledgments

We thank the anonymous referee for his/her fast and detailed comments, which helped improve the paper. P.v.P. thanks J.W. Stock for useful discussions and comments on the text. This research has been supported by the Helmholtz Gemeinschaft through the research alliance “Planetary Evolution and Life”. Pascal Hedelt acknowledges support from the European Research Council (Starting Grant 209622: E₃ARTHs).

Appendix A: Comparative studies used for model validations

A.1. Runs

Several sets of comparison runs were performed to compare the model results of the improved model with other published work and benchmark calculations.

For all comparative runs, the central star was the Sun. The solar input spectrum is based on the high-resolution spectrum provided by Gueymard (2004) (e.g., Kitzmann et al. 2010).

The first set of comparison runs tested the response of the model to doubling/quadrupling the CO₂ content within the modern Earth atmosphere to compare the results with investigations regarding the anthropogenic greenhouse effect (runs W1-W8 in Table A.1).

The second set of comparison runs that tested the model sensitivity to changes in CO₂ content was performed with a reduced solar luminosity S of 0.8 times the present value to compare to published calculations of different models with the same assumptions (runs S1-S6 in Table A.1).

The simulation to which the runs from Table A.1 were compared included an adjustment of the model surface albedo in a way that modern Earth reference calculations yield a surface temperature of 288 K. For the model used here, however, a reference surface temperature of 284.5 K was chosen to account for the missing greenhouse effect of methane, ozone, and nitrous oxide, which are not included in the present model. This value of the surface temperature is calculated with a modern Earth model (Grenfell et al. 2007b, when using the RRTM radiative transfer scheme developed by Mlawer et al. 1997) upon excluding the greenhouse effect provided by methane, ozone, and nitrous oxide. The model surface albedo is 0.24, which is slightly higher than the value of 0.21 from von Paris et al. (2008), but still compatible with surface albedos of other cloud-free models in the literature (e.g., Goldblatt et al. 2009a; Haqq-Misra et al. 2008, both studies using 0.23 as their surface albedo value). Note that the actual global value for Earth is approximately 0.13 (Kitzmann et al. 2010; Rossow & Schiffer 1999), as stated above.

The third set of comparison runs followed a proposed evolutionary sequence of the Earth’s atmosphere based on Hart (1978) in terms of time t_b before present (runs H1-H8 in Table A.2). To this end, the model surface albedo was adjusted to yield a surface temperature of 288 K for present luminosity to be consistent with model studies using this sequence, which also introduced this particular adjustment (Kasting & Ackerman 1986; Kiehl & Dickinson 1987). The resulting surface albedo was 0.185, in contrast to the assumed 0.24 for all other runs. The assumed N₂ partial pressure was 770 mb.

A.2. Results from comparative studies

Figure A.1 summarizes the results of the runs of set W (Table A.1), i.e. tests for the assumed anthropogenic greenhouse effect. Shown is the increase ΔT of surface temperature with respect to run W1, which is a measure of the strength of the greenhouse effect caused by CO₂.

Figure A.1 also shows values for the same numerical experiment performed with two other models. The first model is actually a version of the model used here, but it uses RRTM (Mlawer et al. 1997) instead of MRAC (see model description above) in the IR radiative transfer. The second model is from Goldblatt et al. (2009b) (also used in Goldblatt et al. 2009a), which is also a 1D radiative-convective model. It differs from the model used here in that the correlated-k radiative transfer in the IR is based on Hitran 1992. Furthermore, the radiative transfer uses a 2-stream method instead of the diffusivity approximation used by MRAC. The main difference though is the numerical scheme employed to reach the steady state atmosphere. Goldblatt et al. (2009b) use a Newton-Raphson method, not a time-stepping algorithm. Figure A.1 shows that the agreement between the three models is indeed very good.

Figure A.2 compares the calculations for set S from Table A.1 with calculations done by Haqq-Misra et al. (2008) (similar to runs S3 and S4). They used an updated version (e.g., new absorption coefficients) of the 1D radiative-convective model of Kasting & Ackerman (1986). The difference between their model and ours is that the IR radiative transfer is performed by exponential-sum fitting of transmission instead of correlated-k. Figure A.2 implies that the results from Haqq-Misra et al. (2008) and this work for the two specific runs S3 and S4 agree very well.

The results of the model calculations for set H from Table A.2 were compared with the model results of Kiehl & Dickinson (1987) and Kasting & Ackerman (1986) (see Figs. 7 and 8 in Kiehl & Dickinson 1987; Fig. 2 in Kasting 1989) and are shown in Fig. A.3. All models included in this test are 1D radiative-convective models.

The model by Kasting & Ackerman (1986) uses band models and exponential sums in the radiative transfer instead of the correlated-k. Also, they calculated the relative humidity in the atmosphere by a slightly different approach than the one used by us.

The study by Kiehl & Dickinson (1987) uses a relatively high-resolution radiative transfer code (500 intervals in the thermal region) based on band models. A second difference is the convective adjustment scheme that is not based on a simple formulation of the adiabatic lapse rate, the approach used by Kasting & Ackerman (1986) and here.

Calculated surface temperatures are within 2 K of the results of the two other studies. These differences are quite small given that the models differ with respect to e.g., radiative transfer and the treatment of convection, as stated above.

A.3. Conclusions from comparative studies

The purpose of this section was to test the updated and improved model against other published model results. It was shown that the model agrees well with recent models of von Paris et al. (2008) and Haqq-Misra et al. (2008) as well as earlier model calculations (Kiehl & Dickinson 1987; Kasting & Ackerman 1986). Furthermore, by comparing with benchmark studies regarding the anthropogenic greenhouse effect, the model was found to agree also very well with other model results (e.g., Goldblatt et al. 2009b).

Appendix B: Sensitivity studies for the effect of radiative transfer approximations on surface temperature

In order to assess the sensitivity of the surface temperatures calculated for the atmospheric scenarios of Table 3 to details in the radiative transfer (such as continuum absorption or line mixing), an additional sensitivity study was performed. Based on this sensitivity study, the results regarding surface temperatures and habitability can be interpreted with more confidence. Also, it serves as an illustration of the uncertainties associated with atmospheric modeling in general.

The sensitivity runs have been performed for the simulations with medium and high CO₂ concentrations (sets G2–G3 in Table 3), because these runs resulted in habitable surface conditions.

B.1. Performed runs

Both H₂O and CO₂ show significant collision-induced continuum absorption in the mid- to far-infrared. Especially the formulations of the foreign continua of H₂O and CO₂ could be relatively uncertain, because they are calculated for N₂-O₂ background atmospheres (see, e.g. Halevy et al. 2009).

The first set of sensitivity runs was performed in order to assess the influence of the H₂O continuum on the surface temperature (runs CH1-CH4 in Table B.1). These tests were done for the 5, 10, and 20 bar atmosphere of set G3 as well as the 20 bar run from set G2 (Table 3) by multiplying the H₂O continuum by arbitrary factors f_H2O of 0.2 and 0.5.

The second set of sensitivity runs was performed to assess the influence of the CO₂ self-continuum on the surface temperature (runs CS1-CS4 in Table B.1). The CO₂ self-continuum was multiplied by arbitrary factors f_CO2 of 0.2 and 0.5, respectively. To test the influence of CO₂ foreign-continuum absorption, the continuum was removed for the same runs, i.e. multiplied by 0 (runs CF1-CF4 in Table B.1). These tests were performed for the same runs as the H₂O sensitivity runs.

An additional challenge in the radiative transfer arises from so-called line-mixing. Assuming a Voigt line profile (i.e. a convolution of Lorentz and Gauss profiles) is no longer justified. Comparisons of computer simulations with experimental data by Rodrigues et al. (1999) have shown that including line-mixing into the calculations can result in a decrease of absorption coefficients by up to a factor of 2. Line-mixing parameters for CO₂ are included in the HITRAN database (Rothman et al. 2005), but not in the HiTemp database (Rothman et al. 1995) used for the absorption coefficients in the model. Hence, a Voigt profile was assumed for all lines at all pressures during the calculations of the absorption coefficients (von Paris et al. 2008). In order to test the sensitivity of the results to line-mixing, a third set of sensitivity runs was performed (runs LM1-LM4 in Table B.1). Here, the optical depth in the main atmospheric band of CO₂ (15 μm) was arbitrarily decreased in the lower troposphere (i.e. at pressures higher than 100 mb) by factors f_LM of 0.2 and 0.5, respectively. Again, these tests were performed for the same runs as the H₂O sensitivity runs.

Table B.1 summarizes the sensitivity runs.

B.2. Effect of H₂O continuum

The runs of sets CH1-CH4 from Table B.1 were done in order to test the influence of the H₂O continuum on calculated surface temperatures. The obtained results did not differ significantly from the control runs. For a decrease of H₂O continuum of a factor of 5, the corresponding decrease of surface temperature was less than 2 K in each case.

Thus, these results indicate that the formulation of the H₂O continuum is not critical for the calculations presented here.

B.3. Effect of CO₂ continuum

Sensitivity tests were performed to assess the influence of the CO₂ continuum on surface temperature (sets CS1-CS4 and CF1-CF4 from Table B.1). For the CO₂ self-continuum, the effect was relatively large. The strong mid-IR band around 7 μm dominates absorption and is the main opacity source in this spectral region. Also, the far-IR bands show significant contribution to the overall absorption. Hence, the effect of a decrease of the CO₂ continuum can be expected to be very strong.

For set CS1, the surface temperature decreased from 357 K (control run) to 346 K and 329 K on decreasing the continuum absorption by a factor of 2 and 5, respectively.

Table B.2 summarizes the effect on surface temperatures.

From the results in Table B.2, one infers that the uncertainties in the CO₂ self continuum opacity will not alter the principle conclusions on surface habitability provided that they do not exceed a factor of 2–5. The results of the sensitivity tests imply, however, that more detailed modeling and measurements of the CO₂ self-continuum absorption are needed in the future.

Upon removing the CO₂ foreign continuum (sets CF1-CF4 in Table B.1), surface temperatures decrease by 1–3 K for the high CO₂ runs (CF1-CF3). For the 20 bar medium CO₂ run, however, excluding the CO₂ foreign continuum decreased the surface temperature by 33 K to 280 K. This implies that the foreign continuum is an important opacity source and should be included in all future simulations.

B.4. Effect of line mixing

The effect of line mixing on surface temperatures was investigated with the sets LM1-LM4 from Table B.1.

For all cases, the surface temperatures were almost unaffected (less than 1 K decrease). These effects are important for interpreting measurements on Mars or Venus, but are not likely to significantly affect results regarding exoplanets, where only first-order estimates can be done so far.

B.5. Summary of sensitivity studies

The purpose of this section was to investigate the influence of the radiative transfer approximations on the calculated surface temperatures.

It was shown that the model results do depend to a certain extent on the details of the radiative transfer. However, only the effect of the CO₂ continuum is likely to be important for our study.

References

All Tables

All Figures