© ESO, 2017

1. Introduction

Clouds can have an important impact on the climate of terrestrial planets by either trapping the infrared radiation in the lower atmosphere (greenhouse effect) or by scattering incident stellar radiation back to space (albedo effect). The position and extension of the habitable zone (HZ) around different types of stars thus depends on the presence of clouds (Marley et al. 2013). Especially, the outer HZ boundary might be influenced by the formation of CO₂ ice clouds and their corresponding climatic impact. If a planet is located farther away from its host star, its atmospheric and surface temperatures become cooler due to the decrease in stellar insolation. To sustain liquid water on the surface, a thick atmosphere composed of a greenhouse gas, such as CO₂, is required. If the terrestrial planet is still geologically active, CO₂ can accumulate in the atmosphere by volcanic outgassing (see Pierrehumbert 2010, for an overview). With decreasing atmospheric temperatures, carbon dioxide will condense at some point to form clouds composed of CO₂ ice crystals.

In contrast to other types of condensates important for habitable, terrestrial planets – such as liquid H₂O or water ice – dry ice is more or less transparent in the infrared except within a few strong absorption bands (Hansen 1997, 2005; see the corresponding refractive index in Fig. 1 in Kitzmann et al. 2013). Thus, as argued by Kasting et al. (1993) or Forget & Pierrehumbert (1997), a classical greenhouse effect by absorption and re-emission of thermal radiation is unlikely to occur for CO₂ ice clouds. However, as pointed out by Forget & Pierrehumbert (1997) and Pierrehumbert & Erlick (1998), CO₂ ice particles can efficiently scatter thermal radiation back to the planetary surface thereby creating a scattering greenhouse effect. Depending on the cloud properties a scattering greenhouse effect can outweigh the cloud’s albedo effect and can, in principle, increase the surface temperature above the freezing point of water (see also Mischna et al. 2000; Colaprete et al. 2003, for studies of CO₂ clouds in the early Martian atmosphere).

Most atmospheric modelling studies on the climatic effects of CO₂ clouds so far have been limited to the early Martian atmosphere. For a fully cloud-covered early Mars with a thick CO₂ dominated atmosphere and CO₂ clouds composed of spherical CO₂ ice particles, Forget & Pierrehumbert (1997) determined that in contrast to the cloud-free distance of 1.67 au by Kasting et al. (1993), the outer boundary of the HZ should be located at 2.4 au. This value has been further used by Selsis et al. (2007) to extrapolate the effects of CO₂ clouds on the outer HZ boundary towards other main sequence central stars, assuming that the radiative effects of CO₂ clouds are not a function of the incident stellar radiation or the properties of the planet and its atmosphere. The cloud-free outer HZ has recently been slightly revised by Kopparapu et al. (2013) who used an updated version of the Kasting et al. (1993) model. The corresponding results of these earlier studies on the outer boundary of the classical habitable zone are summarised in Fig. 1.

Fig. 1 Estimates for the outer boundary of the classical habitable zone from previous studies. The red and blue curves show the cloud-free boundaries from Kasting et al. (1993) and Kopparapu et al. (2013), respectively. The results for an atmosphere with 100% CO2 clouds from Selsis et al. (2007) is shown in yellow. This estimate is based on a calculation by Forget & Pierrehumbert (1997), which is additionally marked in the figure.

For the super-Earth Gliese 581 d orbiting the M3V dwarf Gliese 581, a simplified description of CO₂ ice particle formation was included in the one-dimensional (1D) and three-dimensional (3D) atmospheric models by Wordsworth et al. (2010b, 2011). These studies stated that the CO₂ clouds contribute to the greenhouse effect by increasing the planet’s surface temperature. However, they also employed a simplified treatment of the radiative transfer by using two-stream methods.

The first extensive study on the radiative effect of CO₂ clouds in atmospheres of terrestrial planets around main-sequence dwarf stars was done by Kitzmann et al. (2013). In that radiative transfer study, the properties of the CO₂ ice particles (particle sizes, optical depths) were varied over a large parameter range to calculate their radiative effects using a high-order discrete ordinate radiative transfer method. The results of this study suggest, that the simplified two-stream radiative transfer methods employed in previous model studies strongly overestimated the scattering greenhouse effect of CO₂ clouds, concluding that more accurate radiative transfer schemes are absolutely required for atmospheric models to predict the correct climatic impact of CO₂ ice particle clouds.

Following this numerical radiative transfer study, Kitzmann (2016) reinvestigated the impact of CO₂ ice clouds in the atmosphere of early Mars by using a radiative-convective atmospheric model with an accurate radiative transfer. While the results in Kitzmann (2016) suggest, that carbon dioxide clouds still yield a net greenhouse effect under certain conditions, the impact on the surface temperature is much less pronounced than found in the previous studies on this topic. This reduced heating effect will also affect the previous estimates on the extension of the HZ due to the presence of CO₂ clouds.

In this work, I study the climatic effect of CO₂ ice clouds in CO₂-dominated atmospheres of terrestrial planets around different types of main-sequence dwarf stars by using a 1D radiative-convective atmospheric model with an accurate multi-stream radiative transfer. In particular, the impact of carbon dioxide clouds on the position of the outer boundary of the classical habitable zone is investigated. Section 2 gives an overview of the atmospheric model used in this study, as well as the CO₂ cloud description. The climatic impact of CO₂ clouds is discussed in Sect. 3, while their effect on the outer HZ boundary is studied in Sect. 4. Concluding remarks and a summary are given in Sect. 5.

2. Model description

For the calculations in this publication, I use a 1D radiative-convective atmospheric model, previously used to study the climatic impact of CO₂ clouds in the atmosphere of early Mars. A detailed model description is presented in the following.

The model features a state-of-the-art radiative transfer treatment based on opacity sampling and a general discrete ordinate method, able to accurately treat anisotropic scattering. The model currently considers N₂, CO₂, and H₂O as atmospheric species.

In the dry atmosphere, N₂ and CO₂ are considered to be well mixed, that is, they have a constant mixing ratio throughout the atmosphere. For water, the relative humidity profile of Manabe & Wetherald (1967) is used, with a fixed relative humidity at the surface of 77%. Godolt et al. (2016) showed that the results using this relative humidity profile are a good approximation in comparison to full 3D studies of cold atmospheres.

The atmospheric model is stationary, that is, it doesn’t contain an explicit time dependence. It does, however, use the usual approach of time stepping to calculate the atmospheric temperatures (e.g. Kasting et al. 1993). The radiative equilibrium temperatures are obtained via $\frac{\mathrm{d}}{\mathrm{d}\mathit{t}}\mathit{T}\mathrm{\left(}\mathit{z}\mathrm{\right)}\mathrm{=}\mathrm{-}\frac{\mathit{g}}{{\mathit{c}}_{\mathit{p}}\mathrm{\left(}\mathit{z}\mathrm{\right)}}\frac{\mathrm{d}\mathit{F}\mathrm{\left(}\mathit{z}\mathrm{\right)}}{\mathrm{d}\mathit{p}\mathrm{\left(}\mathit{z}\mathrm{\right)}}\mathit{,}$(1)where g is the gravitational acceleration, c_p the heat capacity at constant pressure and F the wavelength-integrated radiation flux. It should be noted that the time t in Eq. (1)needs not to be treated as a real time variable. It is, rather, just an iteration parameter to drive the atmosphere into equilibrium. For the present model, t is not a global constant, but is allowed to change as a function of grid point. This allows, for example, to use relatively large values of t at grid points with very high thermal inertia (usually at high pressures), whereas it can be smaller at the top of the atmosphere to stabilise the convergence in this region.

In regions where the atmosphere is found to be convectively unstable, convective adjustment is performed. The convective lapse rate for regions with CO₂ and H₂O condensation is modelled after Wordsworth & Pierrehumbert (2013).

Approximately one hundred grid points are used to discretise the vertical extension of the atmosphere. For the surface albedo, I use the mean Earth-like value of 0.13 (Kitzmann et al. 2010).

2.1. Radiative transfer

A single radiative transfer scheme with a multi-stream discrete ordinate method is used for the entire wavelength range from 0.1 μm to 500 μm. This method solves the transfer equation $\mathit{\mu }\frac{\mathrm{d}{\mathit{I}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{\tau }}_{\mathit{\lambda }}\mathit{,\mu }\mathrm{\right)}}{\mathrm{d}{\mathit{\tau }}_{\mathit{\lambda }}}\mathrm{=}{\mathit{I}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{\tau }}_{\mathit{\lambda }}\mathit{,\mu }\mathrm{\right)}\mathrm{-}{\mathit{S}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{\tau }}_{\mathit{\lambda }}\mathit{,\mu }\mathrm{\right)}$(2)at several distinct values of the angular variable μ (streams) and afterwards computes the averaged radiation field quantities (radiation flux, mean intensity) by angular integration. The general source function S_λ(τ_λ) takes the incident radiation from the central star (S_{λ, ∗}), the local thermal emission, described by the Planck function B_λ, and the contributions due to scattering into account. It is given by $\begin{array}{ccc}\hspace{0.17em}\begin{array}{c}\begin{array}{c}\end{array}{\mathit{S}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{\tau }}_{\mathit{\lambda }}\mathrm{\right)}\mathrm{=}\hspace{0.17em}\hspace{0.17em}{\mathit{S}}_{\mathit{\lambda ,}\mathrm{\ast }}\mathrm{\left(}{\mathit{\tau }}_{\mathit{\lambda }}\mathrm{\right)}\mathrm{+}\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{\omega }}_{\mathit{\lambda }}\mathrm{)}\mathit{B}\lambda \begin{array}{c}\end{array}\\ \begin{array}{c}\end{array}\mathrm{+}\frac{{\mathit{\omega }}_{\mathit{\lambda }}}{\mathrm{2}}{\mathrm{\int }}_{-1}^{\mathrm{+}\mathrm{1}}{\mathit{p}}_{\mathit{\lambda }}^{\mathrm{0}}\mathrm{\left(}\mathit{\mu ,}{\mathit{\mu }}^{\mathrm{\prime }}\mathrm{\right)}{\mathit{I}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{\mu }}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{d}{\mathit{\mu }}^{\mathrm{\prime }}\mathit{,}\end{array}& & \end{array}$(3)where ω_λ is the single scattering albedo and ${\mathit{p}}_{\mathit{\lambda }}^{\mathrm{0}}\mathrm{\left(}\mathit{\mu ,}{\mathit{\mu }}^{\mathrm{\prime }}\mathrm{\right)}$ is the azimuthally-averaged scattering phase function.

Here, the scattering phase function is represented as a series of Legendre polynomials (Chandrasekhar 1960) ${\mathit{p}}_{\mathit{\lambda }}^{\mathrm{0}}\mathrm{\left(}\mathit{\mu ,}{\mathit{\mu }}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{=}\sum _{\mathit{n}\mathrm{=}\mathrm{0}}^{\mathrm{\infty }}\mathrm{(}\mathrm{2}\mathit{n}\mathrm{+}\mathrm{1}\mathrm{)}{\mathit{P}}_{\mathit{n}}\mathrm{\left(}\mathit{\mu }\mathrm{\right)}{\mathit{P}}_{\mathit{n}}\mathrm{\left(}{\mathit{\mu }}^{\mathrm{\prime }}\mathrm{\right)}{\mathit{\chi }}_{\mathit{\lambda ,n}}$(4)with the Legendre polynomials P_n(μ) and the phase function moments χ_λ,n. The moments are defined by the integration of the full phase function p_λ(α) with respect to the scattering angle α, weighted by the Legendre polynomials P_n(μ)${\mathit{\chi }}_{\mathit{\lambda ,n}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}{\mathrm{\int }}_{-1}^{\mathrm{+}\mathrm{1}}{\mathit{p}}_{\mathit{\lambda }}\mathrm{\left(}\mathit{\alpha }\mathrm{\right)}{\mathit{P}}_{\mathit{n}}\mathrm{\left(}\cos \mathit{\alpha }\mathrm{\right)}\mathrm{d}\cos \mathit{\alpha }\mathit{.}$(5)For most phase functions, this series is infinite but is, in practice, truncated at a certain n = N_max. As mentioned by Chandrasekhar (1960), the value of N_max is equal to the number of streams used to solve the radiative transfer equation.

To treat the strong forward scattering peak of the phase function in case of large size parameters, the δ-M scaling method from Wiscombe (1977) is used. This method approximates the forward scattering peak by a δ distribution and removes it from the phase function series. The δ-M scaling, however, produces errors in the computed intensities which are therefore corrected by using the method presented in Nakajima & Tanaka (1988).

For the numerical implementation of the discrete ordinate method, I use the C-DISORT radiative transfer code (Hamre et al. 2013). Throughout this study, eight computational streams are used in all calculations. Tests by doubling the number of streams showed no changes in the resulting atmospheric and surface temperatures or radiation fluxes.

2.1.1. Opacity sampling

For the description of wavelength-dependent transport coefficients, I adopt the opacity sampling method. This method has been introduced in the context of the cool atmospheres of late-type stars which are dominated by molecular absorption (Sneden et al. 1976). It has the advantage of operating in the normal wavelength/wavenumber space, which allows the transport coefficients of different species to be directly added (Mihalas 1978). The rationale for employing opacity sampling is based on the fact, that wavelength-averaged quantities, such as the total radiation flux or the mean intensity, converge well before all spectral lines are fully resolved.

In the opacity sampling approach, the equation of radiative transfer is solved at distinct wavelength points, at which its results are identical to a line-by-line radiative transfer method. There are a number of strategies on how to distribute these distinct wavelengths, such that, for example, the wavelength-integrated flux converges for a small number of points.

As mentioned in Kitzmann (2016), the distribution of wavelengths at which the equation of radiative transfer is solved, is treated separately in three different wavelength regions. In the infrared, the points are sampled along the Planck black body curves for different temperatures, adopted from the method published by Helling & Jorgensen (1998). Here, the points are sampled for 30 temperatures between 100 K and 400 K, covering the entire range of atmospheric temperatures encountered in the atmospheric scenarios of this study. Between 0.3 μm and 5 μm, 20 000 wavelength points are distributed logarithmic equidistantly in wavenumbers. For smaller wavelengths, approximately 100 points are used to treat the smooth Rayleigh slope.

Fig. 2 Histogram of the distribution of opacity sampling points as a function of wavelength. The histogram shows the number of wavelength points in 100 bins, equidistantly distributed in log space. The total number of sampling points is approximately 40 000.

Figure 2 shows a histogram of the distribution of opacity sampling points as a function of wavelength. The figure nicely illustrates the high sampling rate in the thermal infrared where most of the atmospheric infrared radiation is transported and the drop-off at the flanks of the high and low temperature Planck curves.

In total approximately 40 000 distinct wavelengths are used in this study. Further increasing the spectral resolution has no impact on the atmospheric fluxes and temperature profiles. In principle, the number of points could also be reduced significantly without any large impact on the temperatures. For example, reducing the wavelengths in the visible and near infrared by 10 000 yields changes in the surface temperatures of only 0.1 K. In fact, reasonable values for the surface temperatures can already be obtained by using only a couple of hundred wavelengths.

2.1.2. Absorption coefficients

The molecular absorption coefficients used in this study are calculated with the open source Kspectrum code (version 1.2.0). It should be noted that the currently available version of Kspectrum contains a bug in the calculation of the sub-Lorentzian line profiles for CO₂. A fixed version of the code has been forwarded to the code’s author but has not been made publicly available so far. Using the HITRAN 2012 database (Rothman et al. 2013), absorption cross-sections of CO₂ and H₂O are obtained for pressures between 10^-6 bar and 300 bar and temperatures between 100 K and 640 K. The cross-sections are tabulated in equidistant wavenumber steps of 0.01 cm^-1. The opacity sampling points are always chosen to coincide with one of the tabulated wavenumbers, thus avoiding spectral interpolation of the cross-sections.

In the case of CO₂, the sub-Lorentzian line profiles of Perrin & Hartmann (1989) are employed. The line profiles are truncated at a distance of 500 cm^-1 from the line centre. The continuum contribution and dimer absorption between 1100 and 2000 cm^-1 is taken from Baranov et al. (2003) while the descriptions in Gruszka & Borysow (1997) are used to treat collision-induced absorption for wavenumbers smaller than 250 cm^-1. For the self and foreign continuum absorption of H₂O, the MT-CKD formulation (Mlawer et al. 2012) is used. Following the requirements of the MT-CKD model, the line wings of H₂O are truncated at distances of 25 cm^-1.

Recently, Ozak et al. (2016) showed the importance of taking into account the line-coupling of the CO₂ absorption lines in the infrared for the climate of the early Mars. The use of purely Lorentzian far wing line shapes can lead to an overestimation of the CO₂ greenhouse effect in the infrared. The effect of line mixing is partly included in the model by using the sub-Lorentzian line profiles, though it may still overestimate the greenhouse warming by the CO₂ gas under certain conditions. Thus, the HZ boundaries presented here could be an upper limit.

2.1.3. Molecular Rayleigh scattering

Molecular Rayleigh scattering (Lord Rayleigh 1899) is included for CO₂, H₂O, and N₂. The corresponding cross-section is computed via ${\mathit{\sigma }}_{\mathrm{Rayleigh}\mathit{,\nu }}\mathrm{=}\frac{\mathrm{24}{\mathit{\pi }}^{\mathrm{3}}{\mathit{\nu }}^{\mathrm{4}}}{{\mathit{n}}_{\mathrm{ref}}^{\mathrm{2}}}\mathrm{·}{\left(\frac{\mathit{n}\mathrm{(}\mathit{\nu }{\mathrm{)}}^{\mathrm{2}}\mathrm{-}\mathrm{1}}{\mathit{n}\mathrm{(}\mathit{\nu }{\mathrm{)}}^{\mathrm{2}}\mathrm{+}\mathrm{2}}\right)}^{\mathrm{2}}\mathrm{·}\mathit{K}\mathrm{\left(}\mathit{\nu }\mathrm{\right)}\mathit{,}$(6)where ν is the wavenumber, n the refractive index, n_ref a reference particle number density, and K the King factor. The King factor describes a correction to account for anisotropic molecules. It can also be written as a function of the depolarization factor D$\mathit{K}\mathrm{\left(}\mathit{\nu }\mathrm{\right)}\mathrm{=}\frac{\mathrm{6}\mathrm{+}\mathrm{3}\mathit{D}\mathrm{\left(}\mathit{\nu }\mathrm{\right)}}{\mathrm{6}\mathrm{-}\mathrm{7}\mathit{D}\mathrm{\left(}\mathit{\nu }\mathrm{\right)}}\mathrm{·}$(7)For water, the refractive index from Wagner & Kretzschmar (2008) and the depolarisation factor of 3 × 10^-4 from Murphy (1977) are adopted. The refractive indices and King factors for CO₂ and N₂ are taken from Sneep & Ubachs (2005). Note that Eq. (13) in Sneep & Ubachs (2005) for the refractive index of CO₂ contains a typographical error. The nominator of the last term should read 0.1218145 × 10^-6 instead of the 0.1218145 × 10^-4 factor stated in their Eq. (13).

2.2. Cloud description

The cloud description is identical to the one also used by Forget & Pierrehumbert (1997), Kitzmann et al. (2013), or Kitzmann (2016). The size distribution of the CO₂ ice particles is described by a modified gamma distribution $\mathit{f}\mathrm{\left(}\mathit{a}\mathrm{\right)}\mathrm{=}\frac{{\mathrm{(}{\mathit{a}}_{\mathrm{eff}}{\mathit{\sigma }}^{\mathrm{)}}}^{\mathrm{2}\mathrm{-}\mathrm{1}\mathit{/}\mathit{\sigma }}}{\mathrm{\Gamma }\mathrm{(}{\frac{\mathrm{1}\mathrm{-}\mathrm{2}\mathit{\sigma }}{\mathit{\sigma }}}^{\mathrm{)}}}\mathit{a}\left(\frac{\mathrm{1}}{\mathit{\sigma }}\mathrm{-}\mathrm{3}\right){\mathrm{e}}^{\mathrm{-}\frac{\mathit{a}}{{\mathit{a}}_{\mathrm{eff}}\mathit{\sigma }}}$(8)where a_eff is the effective particle radius and σ the effective variance, for which a value of 0.1 is used (Forget & Pierrehumbert 1997).

The optical properties of the ice particles are obtained via Mie theory, thereby assuming spherical particles (Mie 1908). The refractive index for dry ice is taken from Hansen (1997, 2005). A plot with the optical properties for some selected values of the effective radius and an optical depth of one can be found in Kitzmann (2016). Throughout this study, the cloud optical depth τ refers to the particular wavelength of λ = 0.1 μm.

The scattering phase function of the CO₂ ice cloud particles is approximated by the analytical Henyey-Greenstein function (Henyey & Greenstein 1941) ${\mathit{p}}_{\mathrm{HG}\mathit{,\lambda }}\mathrm{\left(}{\mathit{g}}_{\mathit{\lambda }}\mathit{,\alpha }\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}\mathrm{-}{\mathit{g}}_{\mathit{\lambda }}^{\mathrm{2}}}{{\mathrm{(}\mathrm{1}\mathrm{+}{\mathit{g}}_{\mathit{\lambda }}^{\mathrm{2}}\mathrm{-}\mathrm{2}{\mathit{g}}_{\mathit{\lambda }}\cos {\mathit{\alpha }}^{\mathrm{)}}}^{\mathrm{3}\mathit{/}\mathrm{2}}}\mathit{,}$(9)where g_λ is the asymmetry parameter obtained from Mie theory.

This function lacks the complicated structure and detailed features of the full Mie phase function, though preserves its average quantities, such as, most notably, the asymmetry parameter g_λ, that is ${\mathit{g}}_{\mathrm{HG}\mathit{,\lambda }}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}{\mathrm{\int }}_{-1}^{\mathrm{+}\mathrm{1}}{\mathit{p}}_{\mathrm{HG}\mathit{,\lambda }}\mathrm{\left(}{\mathit{g}}_{\mathit{\lambda }}\mathit{,\alpha }\mathrm{\right)}\cos \mathit{\alpha }\hspace{0.17em}\mathrm{d}\cos \mathit{\alpha }\mathrm{=}{\mathit{g}}_{\mathit{\lambda }}\mathit{.}$(10)In case of the Henyey-Greenstein function, the phase function moments required for the Legendre series (4)have the simple form ${\mathit{\chi }}_{\mathit{\lambda ,n}}\mathrm{=}{\mathit{g}}_{\mathit{\lambda }}^{\mathit{n}}\mathit{.}$(11)

2.3. Stellar spectra

Stellar spectra for four main-sequence stars with different stellar effective temperatures T_eff are used in this study. This includes σ Boo (F2V, T_eff = 6722 K), the Sun (G2V, T_eff = 5777 K), the young K-dwarf ϵ Eri (K2V, T_eff = 5072 K), and AD Leo (M3.5eV, T_eff = 3400 K). The spectra are a composite of stellar atmosphere models and measurements. Details on the spectra can be found in Kitzmann et al. (2010).

3. Climatic impact of CO₂ ice clouds

In this section, I first investigate the climatic effect of the CO₂ clouds for planets orbiting different central stars. The aim is to test the impact of the stars’ spectral energy distribution on the efficiency of the net greenhouse effect. The results presented in Kitzmann et al. (2013) suggest that the climatic impact is directly affected by the spectral distribution of the incident stellar radiation. M-stars, for example, seem to results in only a very small net greenhouse effect of the CO₂ ice particles. This is caused by the shift of the stellar spectrum more towards the near- and mid-infrared compared to solar-type stars, such as our Sun. The shift of the spectrum to longer wavelengths makes Rayleigh scattering very inefficient and, thus, leads to a reduced net greenhouse effect by the cloud particles. However, Kitzmann et al. (2013) were not able to quantify the impact on the surface temperature.

In order to be comparable, a common scenario is used for all central stars. In each case, the atmosphere is composed of six bar CO₂ gas. Neither N₂ nor H₂O are considered for these calculations. The planets are placed at orbital distances, where a surface temperature of 273.15 K (freezing point of water) is obtained for each central star. For these cases, the corresponding temperature-pressure profiles are shown in Fig. 3, along with the saturation vapour pressure curve of CO₂.

In all cases, the troposphere is separated into two different convective regimes: A dry adiabatic region in the lower troposphere and a moist CO₂ adiabatic temperature profile in the upper part. This upper part indicates the atmospheric region, where a CO₂ ice cloud could potentially form.

Fig. 3 Temperature-pressure profiles for planets with cloud-free atmospheres around four different main-sequence stars. The corresponding atmospheres are composed of six bar of CO2 gas. The black, dashed-dotted line shows the saturation vapour pressure curve of CO2. The areas coloured in blue indicate the positions where CO2 cloud layers are placed in the subsequent calculation.

To study the impact of the CO₂ clouds, a cloud layer is inserted into these saturated regions in the following. For the solar-type stars, the cloud layer is centred around 0.1 bar. This corresponds roughly to the position of the cloud layer in the studies of Forget & Pierrehumbert (1997) and Mischna et al. (2000). In the case of the M-type dwarf, the cloud layer has to be located deeper in the atmosphere (0.4 bar) because the strong absorption of stellar radiation by the atmosphere leads to a warmer middle atmosphere and, thus, to a reduced tropopause height. The positions of the cloud layers is marked in Fig. 3.

A supersaturated CO₂-dominated atmosphere would provide a large amount of condensible material. Thus, the CO₂ ice particles could potentially reach large particle sizes. The only detailed microphysical study by Colaprete et al. (2003) obtained particle sizes of up to 1000 μm for an early Martian atmosphere composed of two bar of CO₂. On the other hand, the more simple cloud schemes in the 3D GCM by Forget et al. (2013) resulted in particle sizes almost two orders of magnitudes smaller. One probable cause for this discrepancy could be the high, critical supersaturation (1.35, experimentally derived by Glandorf et al. 2002) used in Colaprete et al. (2003) to initiate the heterogeneous nucleation. Additionally, the clouds condense from the major atmospheric constituent meaning that a substantial amount of condensible material is available at this high supersaturation. This can result in the growth of very large ice particles.

Fig. 4 Impact of CO2 ice clouds on the surface temperature in atmospheres of planets around different main-sequence stars. Surface temperatures are shown as a function of optical depth (at a wavelength of λ = 0.1 μm) and effective particle size of the CO2 ice clouds. The black contour lines are given in steps of 5 K. The white contour line indicates the clear-sky case.

To cover this parameter space, the effective particle sizes are varied between 10 μm and 500 μm in this study. The resulting surface temperatures for optical depths of the cloud layer between 0.1 and 20 are shown in Fig. 4 for all four central stars. The white contour lines indicate the clear-sky temperatures of 273.15 K.

The results presented in Fig. 4 confirm the findings from Kitzmann et al. (2013) and Kitzmann (2016) with respect to the fact, that CO₂ clouds are only effective in producing a net greenhouse effect within a small parameter range. The extent of this parameter range, on the other hand, depends crucially on the spectral type of central star; this is relatively large for F-type stars and very small for M-dwarfs.

On the other hand, the most effective particle radius for a net greenhouse effect is almost independent of the stellar type and is approximately 25 μm. As discussed in Kitzmann (2016), this is not too surprising because the scattering greenhouse effect occurs in the infrared wavelength region. In order to be effective scatterers in this wavelength region, the particle size must be comparable to the IR wavelengths.

The ability of the cloud layer to limit the loss of shortwave radiation to space due to molecular Rayleigh scattering, however, crucially depends on the incident stellar spectrum. The increased scattering of shortwave radiation between the lower atmosphere and the cloud base allows more radiation to be kept within the atmosphere than in the clear-sky case. Since the molecular Rayleigh scattering is proportional to ∝ λ^-4, it is much less effective for cooler stars because the stellar spectrum is more shifted to near-infrared wavelengths. Thus, the stellar radiation of late-type stars is predominantly absorbed in the upper planetary atmosphere rather than being scattered. Consequently, the net heating effect of the CO₂ ice clouds is considerably smaller in these cases.

Figure 5 shows the resulting surface temperature as a function of cloud layer’s optical depth but for a fixed effective particle radius of 25 μm. As already discussed in Kitzmann (2016) for the conditions of the early Martian atmosphere, the CO₂ clouds are only effective at optical depths of approximately 5 to 8. For larger optical depths, the net heating effect becomes increasingly smaller, up to the point where the presence of the cloud layer results in atmospheric cooling. As already apparent from Fig. 4, the irradiance from an F-type star produces the largest net greenhouse effect, with a temperature of up to 30 K higher than the clear-sky case. The planet around the M-dwarf, on the other hand, doesn’t benefit strongly from the cloud’s greenhouse effect. Here, the largest temperature increases are at most 6 K. On the contrary, for optical depths larger than 10, the CO₂ ice cloud produces a net cooling effect that reaches temperature decreases of more than 20 K at an optical depth of 20.

These results will, of course, have a direct impact on the outer boundaries of the habitable zone around these different central stars. This is investigated further in the following section.

4. Effect of CO₂ clouds on the outer boundary of the habitable zone

4.1. Cloud-free limit

Before the impact of clouds is evaluated, I first perform calculations without the presence of clouds to obtain the cloud-free outer boundary of the habitable zone, as done by, for example, Kasting et al. (1993) or Kopparapu et al. (2013). These publications used an inverse modelling approach. This approach makes assumptions on the temperature profiles in the radiative equilibrium part of the atmospheres, which allows them to obtain the outer HZ boundary without iterating the model into thermal equilibrium. Here, in contrast to that, I use full atmospheric calculations, that is, no a priori assumptions of the temperature profile are made and all results represent fully converged model calculations.

Additionally, Kasting et al. (1993) and Kopparapu et al. (2013) used an increased surface albedo to account for the net cooling effect of water droplet and ice clouds. The surface albedo was tuned to yield a mean Earth surface temperature for an Earth-like planet around the Sun. Assuming that the impact of clouds neither depends on the central star nor on the atmospheric temperatures, this surface albedo was kept constant in all cases. Kitzmann et al. (2010), on the other hand, showed, that the greenhouse effect of the water ice clouds, for example, depends on the atmospheric temperatures. For lower temperatures, this greenhouse effect has an increased efficiency, such that the tuned surface albedo would need to be decreased to capture this effect.

Fig. 5 Surface temperatures as a function of the cloud’s optical depth. Results are shown for the four different central stars and an effective cloud particle radius of 25 μm. The horizontal dotted line denotes the clear-sky case (273.15 K surface temperature).

However, since the outer boundary is determined for a thick CO₂ atmosphere, the exact value of the surface albedo is not of great importance because most of the incident stellar radiation will be absorbed or scattered by the CO₂ molecules before it reaches the surface. Only at low CO₂ surface pressures should deviations be expected to occur.

Following the assumptions of Kasting et al. (1993) and Kopparapu et al. (2013), I use an atmosphere composed of N₂, CO₂, and H₂O. The amount of N₂ is fixed at one bar, while the vertical distribution of H₂O is given by the relative humidity parametrisation of Manabe & Wetherald (1967). The atmospheric CO₂ content, on the other hand, is a free parameter. By varying the amount of stellar insolation, I obtain the orbital distances, where the assumed planet possesses a surface temperature of 273.15 K for a given surface pressure of CO₂. The corresponding results for all four central stars are presented in Fig. 6.

The results clearly show the well-known maximum greenhouse effect of CO₂ (Kasting et al. 1993). Carbon dioxide is only an efficient greenhouse gas below a certain maximum partial pressure. For surface pressures higher than this value, the greenhouse effect is offset by the molecular Rayleigh scattering, such that the net effect of the CO₂ molecules would be a dominating cooling effect. This maximum greenhouse effect is a function of the spectral distribution of the incident stellar radiation. It is more dominant for an F-type star and less effective for late-type stars because their spectra are more shifted towards the near infrared, which makes Rayleigh scattering rather inefficient. The results of Fig. 6 agree overall with the values of Kopparapu et al. (2013). Most deviations can be found at lower CO₂ surface pressures, where the different value of the surface albedo influences the results. Figure 6 allows us to obtain the cloud-free outer boundary of the habitable zone, through use of the minima of the curves, that is, the smallest possible insolations.

Fig. 6 Maximum greenhouse effect of CO2 dominated atmospheres for planets orbiting different main-sequence stars. The figure shows the stellar insolation required to obtain a temperature of 273.15 K as a function of the CO2 surface pressure. The stellar flux at the top of the atmosphere F∗ is depicted in units of the solar constant S⊙. The minima of the curves correspond to the maximum greenhouse effect and indicate the outer boundaries of the cloud-free habitable zone.

These critical insolations for the maximum greenhouse effect are shown in Fig. 7 (upper panel) as a function of the stellar effective temperature. Additionally, the figure depicts the corresponding results from Kopparapu et al. (2013). It should be noted that Kopparapu et al. (2013) use a value of 1360 W m^-2 for S_⊙, whereas here a solar constant of 1366 W m^-2 is employed.

Fig. 7 Cloud-free outer boundary of the habitable zone. The critical insolations F∗ in units of the present day solar constant S⊙ for the maximum greenhouse effect (upper panel) and the corresponding orbital distances (lower panel) of the outer HZ boundary are shown as a function of the stellar effective temperature. The results of Kopparapu et al. (2013) are additionally presented for comparison (red curves).

The deviations found in terms of the stellar insolation at the outer boundary as a function of T_eff are rather small. These minor differences can be explained by the use of different surface albedos and the employed radiative transfer schemes.

The outer boundary position can also be expressed in terms of orbital distance d by the usual relation $\mathit{d}\mathrm{=}\sqrt{{\mathit{S}}_{\mathrm{\odot }}\mathit{/}{\mathit{F}}_{\mathrm{\ast }}}$, where F_∗ is the stellar insolation at the top of the atmosphere and [d] = au. The resulting cloud-free distances are shown in the lower panel of Fig. 7, again with the outer HZ boundary determined by Kopparapu et al. (2013) for comparison. Figure 7 suggests, that the differences in terms of orbital distances between the two different studies are relatively small and almost negligible, considering the different modelling approaches. As noted in the introduction, the revised orbital distances by Kopparapu et al. (2013) also differ slightly from the original values obtained by Kasting et al. (1993).

Following the approach of Selsis et al. (2007), I perform a polynomial fit of the resulting orbital distances as a function of the stellar effective temperature. The results for the outer HZ boundary locations in astronomical units are expressed by a 4th-order polynomial of the form: $\mathit{d}\mathrm{=}{\mathit{c}}_{\mathrm{1}}{\mathit{T}}_{\mathrm{eff}}\mathrm{+}{\mathit{c}}_{\mathrm{2}}{\mathit{T}}_{\mathrm{eff}}^{\mathrm{2}}\mathrm{+}{\mathit{c}}_{\mathrm{3}}{\mathit{T}}_{\mathrm{eff}}^{\mathrm{3}}\mathrm{+}{\mathit{c}}_{\mathrm{4}}{\mathit{T}}_{\mathrm{eff}}^{\mathrm{4}}\mathit{,} \mathrm{with}\mathrm{\left[}\mathit{d}\mathrm{\right]}\mathrm{=}\mathrm{au}\mathit{.}$(12)The corresponding parameters c_i are given in Table 1.

4.2. Impact of CO₂ clouds

In this subsection the impact of the CO₂ ice clouds on the position of the outer HZ boundary is calculated. The calculations are restricted to the most optimum cases, that is, cloud optical depth and particle radii as well as the CO₂ surface pressure are chosen such that the highest possible net heating effect of the CO₂ ice clouds is obtained. Additionally, in accordance with Sect. 3, the cloud coverage is set to 100%. Thus, the values for the outer boundary presented in this section represent the upper limits of the impact of CO₂ ice cloud on the outer HZ location.

The effective radius of particle size distribution is chosen based on the results from the previous subsection. The results suggest, that the most effective particle size is independent from the central type and given by approximately a_eff ≈ 25 μm. This particle size will therefore be used in the following for all calculations.

For the most optimum cases, the atmospheric CO₂ partial pressure and the optical depths of the clouds are slightly higher than what has been presented in the previous section. For example, the maximum greenhouse effect of the clear-sky atmosphere in the G-star case occurs around a CO₂ partial pressure of approximately 5.7 bar (see Fig. 6), while the most effective optical depth of the cloud layer in the 6 bar pure CO₂ atmosphere is approximately 6.5. When varying both, the CO₂ content and the cloud optical thickness, the highest impact on the surface temperature occurs for a carbon dioxide surface pressure of 8 bar and an optical depth of 8. However, the difference between these two cases in terms of the stellar insolation required to obtain a surface temperature of 273.15 K is only about 0.008 S_⊙ and, thus, has no impact on the position of the outer HZ boundary.

Fig. 8 Impact of CO2 ice clouds on the position of the outer HZ boundary. The solid red curve shows the cloud-free boundary, while the solid blue curve shows the effect of clouds on the habitable zone. The results from Selsis et al. (2007; dashed, blue curve) and Forget & Pierrehumbert (1997) for the impact of clouds on the HZ boundary are additionally shown for comparison.

The orbital distances of the HZ boundary for all four central stars are again fitted with the polynomial in Eq. (12), with the resulting parameters given in Table 1. The results of this polynomial fit are shown in Fig. 8 along with the results from Selsis et al. (2007) and Forget & Pierrehumbert (1997) for comparison.

The increase of the orbital distances compared to the cloud-free case is clearly a function of the stellar effective temperature. For a planet around cool M-dwarfs, the increase in distance is relatively small. In the best-case scenario presented in this section, the increase is only approximately 0.05 au. Given the fact that, in reality, the atmospheric and CO₂ cloud properties will be less than optimal, one shouldn’t expect a strong positive effect of carbon dioxide clouds in these cases. On the contrary, for less than optimal cloud properties, the outer boundary could even be shifted to smaller orbital distances than the clear-sky case for cooler central stars. For stars with higher effective temperatures, the increase in orbital distance can be up to 0.5 au. For the F-type star, for example the presence of CO₂ ice clouds would allow the planet to be located at 2.54 au, which is 0.44 au farther than the clear-sky value of 2.1 au.

As mentioned in the introduction, the only atmospheric modelling calculations for the impact of the carbon dioxide clouds on the position of the outer HZ boundary was done by Forget & Pierrehumbert (1997), who state a distance of 2.4 au for the Sun. This result has been used by Selsis et al. (2007) to scale the cloud-free distances from Kasting et al. (1993). It should be noted though, that the absorption coefficients for CO₂ molecules in Kasting et al. (1993) and Forget & Pierrehumbert (1997) were overestimated because the continuum absorption was accounted for twice in their atmospheric model (see Wordsworth et al. 2010a, for details). The less effective greenhouse effect by the CO₂ ice clouds in this study results in a much smaller value for the orbital distance in this case. According to Fig. 8, the orbital distance for a G2V star is about 2 au and, thus, 0.4 au smaller than what has previously been estimated. Compared to the results from Selsis et al. (2007), the new estimates provided by this work shift the boundary towards smaller orbital distances by approximately 0.5 au for F-type stars and approximately 0.2 au for M-dwarfs.

5. Summary

In this study the climatic effects of CO₂ ice clouds in CO₂-dominated atmospheres of terrestrial planets around different types of main-sequence dwarf stars and the impact of CO₂ clouds on the location of the outer boundary of the classical habitable zone are investigated. A radiative-convective atmospheric model employing an accurate discrete ordinate radiative transfer was used to calculate surface temperatures for a broad range of particle sizes and optical depths of CO₂ ice particles.

As a first step, the climatic impact of CO₂ clouds in a thick, CO₂-dominated atmosphere was evaluated for different stellar spectral energy distributions and as a function of the effective particle radius and the cloud’s optical depth. As already pointed out by Kitzmann et al. (2013), the heating and cooling effect should be a direct function of the incident stellar spectrum. The results from this work suggest that for these thick, CO₂-dominated atmospheres the cloud’s radiative forcing yields, at most, temperature increases smaller than 6 K for cool M-stars. On the other hand, central stars with higher effective temperatures resulted in surfaces temperatures approximately 30 K higher than the corresponding clear-sky case. Furthermore, the parameter range where a CO₂ ice cloud can have a net positive effect is rather limited for late-type stars. A heating effect is only obtained in a small parameter range with respect to optical depths and effective particle sizes. Outside this parameter range, the clouds would have a dominating albedo effect, thus cooling the surface. The most effective particle size for the scattering greenhouse effect was found to be approximately 25 μm, independent from the central star type.

Following this, the impact of dry ice clouds on the position of the classical HZ’s outer boundary was investigated. Atmospheric calculations without clouds were performed first to obtain the cloud-free habitable zone boundary. The results of this cloud-free limit are very similar to the orbital distances published by Kopparapu et al. (2013).

The outer boundary influenced by clouds has been calculated for the most optimum scenarios and, thus, represents an upper limit. The additional greenhouse effect allows the orbital distance of the outer HZ boundary to be increased by up to 0.5 au for F-stars. Stars with lower effective temperatures, on the other hand, don’t benefit strongly from the CO₂ ice clouds. The less effective net greenhouse effect results in an extension of only 0.05 au.

Compared to the orbital distance of the outer HZ boundary of 2.4 au for the Sun published by Forget & Pierrehumbert (1997), which has been obtained using a simplified two-stream method, the revised distance obtained by the more accurate radiative transfer treatment in this study is 0.4 au smaller. Thus, all parametrisations of the outer boundary which rely on that value (e.g. Selsis et al. 2007 or Kaltenegger & Sasselov 2011) should be revisited.

It should be noted that the outer HZ boundary calculated in this work only considers the classical CO₂-dominated atmospheres. The presence of other greenhouse gases, such as CH₄, might also lead to an additional extension of the HZ. Other possible mechanisms also include, for example, the greenhouse effect provided by a H₂-dominated atmosphere as studied by e.g. Pierrehumbert & Gaidos (2011), for example.

There are also still open questions remaining regarding the climatic impact of CO₂ clouds. Models including a treatment for CO₂ ice cloud formation have so far strongly disagreed on the resulting particle sizes (e.g. Forget et al. 2013; Colaprete et al. 2003), such that additional research in this area is warranted. Additionally, all models have, so far, assumed spherical particles and used Mie theory to obtain the optical properties of the cloud particles. However, as shown by laboratory measurements of Behnken (1912), Wahl (1913), or Wergin et al. (1997) CO₂ ice crystals can have cubic or octahedral shapes. Combinations of both (cuboctahedra) or more complicated shapes such as rhobic-dodecahedral crystals can also occur. Since the greenhouse effect of carbon dioxide clouds is determined by the scattering properties of the ice particles, their shapes should have a large impact on the resulting optical properties (e.g. scattering phase function) and, thus, also on the climatic impact.

Acknowledgments

D.K. gratefully acknowledges the support of the Center for Space and Habitability of the University of Bern and the MERAC Foundation for partial financial assistance. This work has been carried out within the frame of the National Centre for Competence in Research PlanetS supported by the Swiss National Science Foundation. D.K. acknowledges the financial support of the SNSF.

References

All Tables

All Figures

Fig. 1 Estimates for the outer boundary of the classical habitable zone from previous studies. The red and blue curves show the cloud-free boundaries from Kasting et al. (1993) and Kopparapu et al. (2013), respectively. The results for an atmosphere with 100% CO2 clouds from Selsis et al. (2007) is shown in yellow. This estimate is based on a calculation by Forget & Pierrehumbert (1997), which is additionally marked in the figure.

In the text

	Fig. 2 Histogram of the distribution of opacity sampling points as a function of wavelength. The histogram shows the number of wavelength points in 100 bins, equidistantly distributed in log space. The total number of sampling points is approximately 40 000.
In the text

	Fig. 3 Temperature-pressure profiles for planets with cloud-free atmospheres around four different main-sequence stars. The corresponding atmospheres are composed of six bar of CO2 gas. The black, dashed-dotted line shows the saturation vapour pressure curve of CO2. The areas coloured in blue indicate the positions where CO2 cloud layers are placed in the subsequent calculation.
In the text

	Fig. 4 Impact of CO2 ice clouds on the surface temperature in atmospheres of planets around different main-sequence stars. Surface temperatures are shown as a function of optical depth (at a wavelength of λ = 0.1 μm) and effective particle size of the CO2 ice clouds. The black contour lines are given in steps of 5 K. The white contour line indicates the clear-sky case.
In the text

	Fig. 5 Surface temperatures as a function of the cloud’s optical depth. Results are shown for the four different central stars and an effective cloud particle radius of 25 μm. The horizontal dotted line denotes the clear-sky case (273.15 K surface temperature).
In the text

Fig. 6 Maximum greenhouse effect of CO2 dominated atmospheres for planets orbiting different main-sequence stars. The figure shows the stellar insolation required to obtain a temperature of 273.15 K as a function of the CO2 surface pressure. The stellar flux at the top of the atmosphere F∗ is depicted in units of the solar constant S⊙. The minima of the curves correspond to the maximum greenhouse effect and indicate the outer boundaries of the cloud-free habitable zone.

In the text

Fig. 7 Cloud-free outer boundary of the habitable zone. The critical insolations F∗ in units of the present day solar constant S⊙ for the maximum greenhouse effect (upper panel) and the corresponding orbital distances (lower panel) of the outer HZ boundary are shown as a function of the stellar effective temperature. The results of Kopparapu et al. (2013) are additionally presented for comparison (red curves).

In the text

	Fig. 8 Impact of CO2 ice clouds on the position of the outer HZ boundary. The solid red curve shows the cloud-free boundary, while the solid blue curve shows the effect of clouds on the habitable zone. The results from Selsis et al. (2007; dashed, blue curve) and Forget & Pierrehumbert (1997) for the impact of clouds on the HZ boundary are additionally shown for comparison.
In the text