Microphysical Modeling of Mineral Clouds in GJ1214 b and GJ436 b: Predicting Upper Limits on the Cloud-top Height

We extend the microphysical model originally developed by Ohno & Okuzumi (2017) to predict the vertical distributions of a cloud in the atmosphere. The cloud model of Ohno & Okuzumi (2017) adopts a 1D Eulerian framework, and provides the vertical distributions of number (n_c) and mass (ρ_c) densities of cloud particles by taking into account the vertical transport of cloud particles due to the updraft motion and gravitational settling, and the particle growth via condensation and coalescence (see Sections 2.3 and 2.5). In this study, we take into account the vertical transport of cloud particles via eddy diffusion (e.g., Ackerman & Marley 2001).

Following previous studies, such as Charnay et al. (2015a) and Morley et al. (2013, 2015), we consider the clouds composed of solid KCl particles formed through the condensation of KCl vapor. The initial cloud particles are assumed to form at the cloud base through the condensation of vapor onto the small nuclei that already exist in the atmosphere, the process so called heterogeneous nucleation. On Earth, such small nuclei, called the cloud condensation nuclei (CCNs), include sea salt, volcano ash, and dust from the land (Rogers & Yau 1989). The amount of CCNs on exoplanets is still highly uncertain and their composition is also uncertain, and therefore we take the number density of CCNs as a free parameter. The height of the cloud base is determined from the comparison between the atmospheric temperature and condensation temperature (see Section 2.2). The condensation temperature is defined as the temperature at which the partial pressure of a volatile is equal to its saturation vapor pressure.

Following Ohno & Okuzumi (2017), we assume that the cloud particles have the characteristic radius r_c and corresponding mass m_c = (4π/3)ρ_intr_c³, where ρ_int is the internal density of the particles. The internal density can vary significantly if the particles grow into porous aggregates (Kataoka et al. 2013). In this study we simply assume ρ_int = ρ_p, where ρ_p is the material density of the condensate, but we will discuss the influences of varying the internal density in Section 5.3. Assuming the mass distribution is narrowly peaked at m ≈ m_c, the number and mass densities are related by ρ_c = m_cn_c. Such frameworks are called the double-moment bulk schemes in meteorology (e.g., Ziegler 1985; Ferrier 1994) and the characteristic size method in planetary formation community (e.g., Birnstiel et al. 2012; Ormel 2014; Sato et al. 2016). This method allows us to derive the physical understanding from calculations more clearly, and to perform the calculations with very little computational time compared with spectral bin schemes (e.g., Brauer et al. 2008) that solve the evolution of the full size distribution.

We investigate the influences of atmospheric metallicity on the vertical profiles of clouds in super-Earths. In this paper, the atmospheric metallicity refers to the ratio of atmospheric heavy element abundance to that of the solar atmosphere (i.e., ${({N}_{{\rm{Z}}}/({N}_{{\rm{H}}}+{N}_{\mathrm{He}}))/({N}_{{\rm{Z}}}/({N}_{{\rm{H}}}+{N}_{\mathrm{He}}))}_{\mathrm{solar}}$). Recent theoretical studies suggest that the atmospheres of super-Earths potentially have metallicities higher than solar, and even higher than 100× solar, depending on the properties of the building blocks of the respective planets (Fortney et al. 2013; Venturini et al. 2016). Also, the interior modeling showed that GJ1214 b might have a steam atmosphere mainly composed of water vapor (Rogers & Seager 2010; Valencia et al. 2013). Therefore, we take the atmospheric metallicity as a free parameter widely ranging from the metallicity of 1× solar to water vapor atmosphere. The metallicity difference provides the different pressure–temperature structure, total cloud mass, and eddy diffusion coefficient.

To determine the location of the cloud base, we construct the pressure–temperature structure using the analytical model of radiative atmosphere described by Guillot (2010) under the assumption of hydrostatic equilibrium. Guillot (2010) derived the analytical solution of global mean thermal profiles that gives good agreement with the predictions from sophisticated simulations. The stellar effective temperature, radii, semimajor axis, and planetary radii of GJ1214 b and GJ436 b are taken from the Exoplanet.eu catalog.² Following Guillot (2010), the temperature in each atmospheric layer is given by

$$\begin{eqnarray}{T}^{4} & = & \displaystyle \frac{3{T}_{\mathrm{int}}^{4}}{4}\left[\displaystyle \frac{2}{3}+\tau \right]+\displaystyle \frac{3{T}_{\mathrm{irr}}^{4}}{4}f\\ & & \times \left[\displaystyle \frac{2}{3}+\displaystyle \frac{1}{\gamma \sqrt{3}}+\left(\displaystyle \frac{\gamma }{\sqrt{3}}-\displaystyle \frac{1}{\gamma \sqrt{3}}\exp (-\gamma \sqrt{3}\tau )\right)\right],\end{eqnarray} \tag{ 1 }$$

where τ is the vertical infrared optical depth τ given by

$$\begin{eqnarray}&&\tau (z)={\int }_{z}^{\infty }{\rho }_{{\rm{g}}}{\kappa }_{\mathrm{th}}{dz}^{\prime} ,\end{eqnarray} \tag{ 2 }$$

where κ_th is the atmospheric infrared opacity. The f = 1/4 is the heat redistribution factor under the assumption of the radiation redistributed around the entire planet, T_int is the intrinsic effective temperature, T_irr is the irradiation effective temperature, and the γ = κ_v/κ_th is the ratio of the visible to infrared opacities, respectively. For GJ1214 b, we take T_int = 60 K (Rogers & Seager 2010) and γ = 0.038, so that reproduces the P–T structure predicted by radiative transfer models of Miller-Ricci & Fortney (2010). For GJ436 b, we take T_int = 300 K (Morley et al. 2017) and γ = 0.05 that is in a good agreement with the retrieved P–T structure (Miguel et al. 2015).

We calculate τ using the fitting formula of Rosseland mean opacity of a cloud-flee atmosphere described by Freedman et al. (2014). This fitting formula is a function of atmospheric metallicity, pressure, and temperature, and valid for P = 10⁻⁶–3 × 10² bar and T = 75–4000 K. Although the opacity table for higher metallicity (>50× solar) is not available so far, the fitting formula can provide the qualitative results for such high metallicity atmospheres. For water vapor atmospheres, we use the opacity of 50× solar metallicity that yields the similar P–T structure to that for a water vapor (Miller-Ricci & Fortney 2010) for simplicity. We also neglect the opacity of cloud particles that might change the location of the cloud base, but we plan to investigate these impacts in future study.

Figure 1 shows the vertical P–T structures of GJ1214 b for 1×  and 100× solar metallicity and the condensation temperature at each atmospheric layer. Here we predict the condensation temperature for each volatile using the saturation vapor pressure described in Rogers & Yau (1989), Ackerman & Marley (2001), and Morley et al. (2012). The vapor species has a solid phase if the atmospheric temperature is lower than its condensation temperature. Therefore, for each volatile species, the cloud base is expected to be placed at the location where the P–T curve intersects the curve of condensation temperature of the species. Figure 1 indicates that the KCl, ZnS, and Na₂S are condensible for 1× solar metallicity case, and KCl and ZnS are condensible for 100× solar metallicity. Since the abundance of KCl vapor is higher than that of ZnS vapor for solar like atmosphere (Morley et al. 2012), we focus on the mineral clouds of KCl in this study. The cloud base for KCl is placed at ∼0.4 bar for 1× solar metallicity, ∼0.1 bar for 10× solar metallicity, and ∼0.07 bar for 100× solar metallicity, respectively, which is in good agreement with the prediction of previous studies (Miller-Ricci Kempton et al. 2012; Morley et al. 2013; Charnay et al. 2015a).

Zoom In Zoom Out Reset image size

We calculate the vertical distributions of the number and mass densities of cloud particles by taking into account their growth and vertical transport. The microphysics of cloud formation is complex (see, e.g., Rossow 1978; Rogers & Yau 1989; Pruppacher & Klett 1997; Seinfeld & Pandis 2006). However, Ohno & Okuzumi (2017) showed that inclusion of condensation and collisional growth is enough to approximately reproduce the observations of terrestrial water clouds and Jovian ammonia clouds. Therefore we take into account the condensation and collisional growth in this study.

Following Charnay et al. (2015a), we consider the clouds formed through the large-scale atmospheric motion driven by the intense day–night heating contrast. Previous studies showed that the global averaged distributions of such clouds can be approximately reproduced by a 1D advection-diffusion model with an empirical parameterization of the eddy diffusion coefficient K_z (Parmentier et al. 2013; Charnay et al. 2015a). Hence the master equations used here are constructed by adding the source terms expressing particle growth to the 1D advection-diffusion model:

$$\begin{eqnarray}&&\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}=\displaystyle \frac{\partial }{\partial z}\left[{n}_{{\rm{g}}}{K}_{{\rm{z}}}\displaystyle \frac{\partial }{\partial z}\left(\displaystyle \frac{{n}_{{\rm{c}}}}{{n}_{{\rm{g}}}}\right)+{v}_{{\rm{t}}}(r){n}_{{\rm{c}}}\right]-{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coll}},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\rho }_{{\rm{c}}}}{\partial t}=\displaystyle \frac{\partial }{\partial z}\left[{\rho }_{{\rm{g}}}{K}_{{\rm{z}}}\displaystyle \frac{\partial }{\partial z}\left(\displaystyle \frac{{\rho }_{{\rm{c}}}}{{\rho }_{{\rm{g}}}}\right)+{v}_{{\rm{t}}}(r){\rho }_{{\rm{c}}}\right]+{\left(\displaystyle \frac{\partial {\rho }_{{\rm{c}}}}{\partial t}\right)}_{\mathrm{cond}},\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\rho }_{{\rm{v}}}}{\partial t}=\displaystyle \frac{\partial }{\partial z}\left[{\rho }_{{\rm{g}}}{K}_{{\rm{z}}}\displaystyle \frac{\partial }{\partial z}\left(\displaystyle \frac{{\rho }_{{\rm{v}}}}{{\rho }_{{\rm{g}}}}\right)\right]-{\left(\displaystyle \frac{\partial {\rho }_{{\rm{c}}}}{\partial t}\right)}_{\mathrm{cond}},\end{eqnarray} \tag{ 5 }$$

where v_t is the terminal velocity of cloud particles, K_z is the eddy diffusion coefficient, and ρ_v is the vapor mass density. The terminal velocity depends on the particle size and atmospheric density, as introduced in Section 2.4. Each source term, introduced in Section 2.5, expresses the particle growth via condensation and collision of each particle. Without these terms, the Equations (3)–(4) are reduced to the 1D transport model for fixed size particles used by Parmentier et al. (2013) and Charnay et al. (2015a).

The eddy diffusion coefficient K_z represents the strength of effective vertical mixing for cloud particles. In this study, we adopt the empirical formula of K_z proposed by Charnay et al. (2015a),

$$\begin{eqnarray}&&{K}_{{\rm{z}}}={K}_{0}{\left(\displaystyle \frac{P}{{P}_{0}}\right)}^{-2/5},\end{eqnarray} \tag{ 6 }$$

where K₀ is the value of K_z at a reference pressure P₀. Charnay et al. (2015a) derived this formula from 3D GCM simulations that takes into account the transport of fixed size particles. Since they suggested that K_z is almost independent of particle size (see the Figure 14 of Charnay et al. 2015a), we use Equation (6) for all ranges of particle size in our calculations. The exponent of −2/5 is similar to the K_z ∝ P^−1/3 predicted by mixing theory (Ackerman & Marley 2001) and K_z ∝ P^−1/2 predicted by other GCM simulations for hot Jupiter (Parmentier et al. 2013). According to Charnay et al. (2015a), we choose the reference pressure of P₀ = 1 bar and take the values of K₀ as summarized in Tables 1 and 2. For GJ1214 b, we use the values of K₀ derived from the power-law fitting to the GCM data (see Figure 14 of Charnay et al. 2015a), which are metallicity-dependent. For GJ436 b, we adopt K₀ = 2.5 × 10³ m² s⁻¹ independent of the metallicity. The adopted value is similar to the value K_z ∼ 10²–10³ m² s⁻¹ suggested by Madhusudhan & Seager (2011) to explain the disequilibrium chemistry seen in GJ436 b's emission spectra. Our adopted value is two orders of magnitude lower than the earlier prediction by Lewis et al. (2010) based on the rms velocity from 3D GCM calculations. However, Parmentier et al. (2013) and Charnay et al. (2015a) recently pointed out that the values of the eddy diffusion coefficient predicted in this way tend to be one or two orders of magnitude higher than those directly determined from vertical particle distribution. If we take this into account, our choice of K₀ is consistent with the GCM results by Lewis et al. (2010).

A terminal velocity v_t is determined by the balance between gravitational force and gas frictional force. The gas frictional force depends on the behavior of the gas flow around the settling particles, and varies with the particle size, settling velocity, and the mean free path of gas particles (e.g., Rossow 1978; Woitke & Helling 2003). In this study, we adopt the following formula of the terminal velocity,

$$\begin{eqnarray}&&{v}_{{\rm{t}}}({r}_{{\rm{c}}})=\displaystyle \frac{2\beta {{gr}}_{{\rm{c}}}^{2}{\rho }_{{\rm{p}}}}{9\eta }{\left[1+{\left(\displaystyle \frac{0.45{{gr}}_{{\rm{c}}}^{3}{\rho }_{{\rm{g}}}{\rho }_{{\rm{p}}}}{54{\eta }^{2}}\right)}^{2/5}\right]}^{-5/4},\end{eqnarray} \tag{ 7 }$$

where η is the dynamic viscosity of the atmosphere and β is the slip correction factor. β accounts for the transition of gas drag behavior from viscous flow (Stokes's law) to free molecular flow (Epstein's law) around the particle, given by (Davies 1945)

$$\begin{eqnarray}&&\beta =1+{\mathrm{Kn}}_{{\rm{g}}}[1.257+0.4\exp (-1.1/{\mathrm{Kn}}_{{\rm{g}}})],\end{eqnarray} \tag{ 8 }$$

where Kn_g = l/r_c is the gas Knudsen number and l is the gas mean free path given in the Appendix. Equation (7) without β is same as the Equation (23) in Ohno & Okuzumi (2017) that asymptotically reaches the Stokes's law for a laminar flow limit, Newton's law for a turbulent flow limit, and reproduces the intermediate regime predicted by the experiment (see the Figure 7 in Ohno & Okuzumi 2017). The top panel of Figure 2 shows the terminal velocity as a function of particle size and atmospheric pressure. Figure 2 shows the terminal velocity increases with height in the upper atmosphere because of Epstein's law arisen from the low atmospheric density.

Zoom In Zoom Out Reset image size

We also show the ratio of the mixing timescale τ_mix to the falling timescale τ_fall in the bottom panel of Figure 2. Each timescale is defined as

$$\begin{eqnarray}&&{\tau }_{\mathrm{mix}}=\displaystyle \frac{{H}^{2}}{{K}_{{\rm{z}}}}\end{eqnarray} \tag{ 9 }$$

and

$$\begin{eqnarray}&&{\tau }_{\mathrm{fall}}=\displaystyle \frac{H}{{v}_{{\rm{t}}}},\end{eqnarray} \tag{ 10 }$$

where H = k_BT/m_gg is the pressure scale height, respectively. Here we assume T = 500 K, K₀ = 10³ m² s⁻¹, and a solar composition atmosphere. Cloud particles ascend if τ_mix ≪ τ_fall, and fall if τ_mix ≫ τ_fall. Figure 2 indicates that the cloud particles are required to maintain their size r_c ≲ 5 μm to ascend above 10⁻³ bar suggested for GJ436 b (Knutson et al. 2014a), and r ≲ 0.3 μm to ascend above 10⁻⁵ bar under the assumed parameters.

The cloud particles ascend from the cloud base while growing through condensation and collision with each other. Condensation dominates the growth of small particles due to the relatively short timescale. The growth rate of ρ_c via condensation depends on the behavior of vapor molecule motion, and is expressed by (Rogers & Yau 1989; Woitke & Helling 2003)

$$\begin{eqnarray}{\left(\displaystyle \frac{\partial {\rho }_{{\rm{c}}}}{\partial t}\right)}_{\mathrm{cond}} & = & 4\pi {r}_{{\rm{c}}}^{2}{n}_{{\rm{c}}}({\rho }_{{\rm{v}}}-{\rho }_{{\rm{s}}})\\ & & \times \min \left[{C}_{\mathrm{re}},\displaystyle \frac{D}{{r}_{{\rm{c}}}}{\left(1+\left(\displaystyle \frac{{m}_{{\rm{v}}}L}{{k}_{{\rm{B}}}T}-1\right)\displaystyle \frac{{LD}{\rho }_{{\rm{s}}}}{{KT}}\right)}^{-1}\right],\end{eqnarray} \tag{ 11 }$$

where ρ_v is the vapor mass density, ρ_s is the saturation vapor density, ${C}_{\mathrm{re}}=\sqrt{{k}_{{\rm{B}}}T/2\pi {m}_{{\rm{v}}}}$ is the relative velocity of vapor molecules, m_v is the mass of the vapor molecules, L is the specific latent heat of condensation, and D is the molecular diffusion coefficient of vapor in ambient air, respectively. The first formula in the bracket corresponds to the free molecular flow regime (Woitke & Helling 2003) in which the vapor molecules are freely impinging onto the particles. The second formula corresponds to the diffusive regime (Rogers & Yau 1989) in which the vapor molecules behave as continuum.

Collisional growth is induced by the relative velocity arisen from both gravitational settling and Brownian motion of particles. In this paper, we refer to the collisional growth by gravitational settling as coalescence and that by Brownian motion as coagulation. Then the decrease in number density via collisional growth is expressed by

$$\begin{eqnarray}&&{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coll}}={\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coag}}+{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coal}},\end{eqnarray} \tag{ 12 }$$

where $| \partial {n}_{{\rm{c}}}/\partial t{| }_{\mathrm{coag}}$ is the decrease in number density for coagulation and $| \partial {n}_{{\rm{c}}}/\partial t{| }_{\mathrm{coal}}$ is that for coalescence. The expression of $| \partial {n}_{{\rm{c}}}/\partial t{| }_{\mathrm{coag}}$ depends on the particle Knudsen number Kn_p defined as

$$\begin{eqnarray}&&{\mathrm{Kn}}_{p}=\displaystyle \frac{\beta }{6\eta {r}_{{\rm{c}}}^{2}}\sqrt{\displaystyle \frac{{m}_{{\rm{c}}}{k}_{{\rm{B}}}T}{2\pi }}.\end{eqnarray} \tag{ 13 }$$

The Brownian motion of particles is diffusive for Kn_p ≪ 1 and ballistic for Kn_p ≫ 1. The rate of decrease of particle number density via coagulation is given by (Seinfeld & Pandis 2006)

$$\begin{eqnarray}{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coag}}=\left\{\begin{array}{ll}8\sqrt{\displaystyle \frac{\pi {k}_{{\rm{B}}}T}{{m}_{{\rm{c}}}}}{r}_{{\rm{c}}}^{2}{n}_{{\rm{c}}}^{2} & ({\mathrm{Kn}}_{p}\gt 1/\sqrt{2})\\ \displaystyle \frac{4{k}_{{\rm{B}}}T\beta }{3\eta }{n}_{{\rm{c}}}^{2} & ({\mathrm{Kn}}_{p}\lt 1/\sqrt{2}).\end{array}\right.\end{eqnarray} \tag{ 14 }$$

The transition takes place at r_c ≈ 0.07 μm under the assumptions of T = 1000 K, P = 0.1 bar, and m_g = 2 amu, which are equivalent to the parameters for the cloud base.

For the coalescence growth, the rate of decrease of number density $| \partial {n}_{{\rm{c}}}/\partial t{| }_{\mathrm{coal}}$ is given by (Rossow 1978)

$$\begin{eqnarray}&&{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coal}}\approx 2\pi {r}_{{\rm{c}}}^{2}{n}_{{\rm{c}}}^{2}{\rm{\Delta }}{vE},\end{eqnarray} \tag{ 15 }$$

where Δv is the relative velocity induced by the gravitational settling, and E is the collection efficiency defined as the ratio of the effective collisional cross section to the geometric cross section (e.g., Pruppacher & Klett 1997). For the relative velocity, Sato et al. (2016) and Krijt et al. (2016) showed that the characteristic size approach with Δv = 0.5v_t(r) is in good agreement with the results of spectral bin schemes, and therefore we assume Δv = 0.5v_t(r_c). The collection efficiency E accounts for the effect of the gas flow around the particle moving relative to the background gas, and is expressed in terms of Stokes number

$$\begin{eqnarray}&&\mathrm{Stk}=\displaystyle \frac{{v}_{{\rm{t}}}({r}_{{\rm{c}}}){\rm{\Delta }}v}{{{gr}}_{{\rm{c}}}},\end{eqnarray} \tag{ 16 }$$

which is defined as the ratio of the stopping time = v_t(r_c)/g to the crossing time ∼r_c/Δv. When Stk ≪ 1, the particles are strongly coupled to the gas flow around the other particles, and hence E behaves as E ≈ 0 (Rossow 1978). We evaluate E using a smoother analytic function of Guillot et al. (2014) given by

$$\begin{eqnarray}&&E=\max [0,1-0.42{\mathrm{Stk}}^{-0.75}],\end{eqnarray} \tag{ 17 }$$

which vanishes at Stk ≲ 0.3 and approaches unity at Stk ≫ 1. If Kn_g > 1, we assume E = 1 because the influence of the gas on the particle trajectory should be weak in that region (Rossow 1978).

We numerically solve the Equations (3)–(5) until the system reaches to the steady-state profiles. The initial number density of the cloud particles at the cloud base is parameterized by the CCN number density n_CCN. We take the n_CCN as a free parameter widely ranging as 10⁶–10¹⁵ m⁻³. Since the composition of the CCNs in close-in super-Earths is unknown, we assume the bulk density of the CCN as that of KCl. This assumption does not affect the calculated cloud vertical profiles as long as the mass fraction of the CCNs in the cloud particles is very small. Therefore, we choose the upper limit of n_CCN so that the total mass of CCNs does not exceed that of KCl vapor at the cloud base. We set the radii of CCNs as r_CCN = 0.001 μm, and then n_CCN ≈ 10¹⁵ m⁻³ corresponds to the upper limit for our calculations.

We choose the flux of a lower boundary condition so that n_c/n_g, ρ_c/ρ_g, and ρ_v/ρ_g keep the values of the cloud base. We adopt the zero-flux boundary condition at the top of the computational domain, which is located at P = 10⁻⁸ bar. The vertical coordinate z is discretized into linearly spaced bins. We use the different grid width for different atmospheric metallicity, as summarized in Table 1. The time increment Δt is chosen at every time step so that the fractional decreases in n, ρ_c, and ρ_v do not exceed 0.5 (i.e., Δt ≤ −0.5 × min[(∂ln n/∂t)⁻¹, (∂ln ρ_c/∂t)⁻¹, (∂ln ρ_v/∂t)⁻¹]). However, this expression yields very small Δt because the time increment determined by condensation is much shorter than that for collisional growth and vertical transport. To avoid this, we adjust the time increment as Δt ≤ −0.5 × (∂ ln n/∂ t)⁻¹ if ${(\partial \mathrm{ln}{\rho }_{{\rm{c}}}/\partial t)}^{-1}\lt 0.1\times {(\partial \mathrm{ln}n/\partial t)}^{-1}$. In this case, we convert the all excess/lack of vapor from saturation value into cloud particles.

We calculate the mean molecular weight of the atmosphere assuming a hydrogen–helium–water mixture in accordance with Fortney et al. (2013). Elemental abundances are taken from Lodders (2003). The mixing ratio of KCl vapor q_KCl below the cloud base is calculated assuming the number of KCl molecules is equal to that of K. For the steam atmosphere and the metallicity of 1000× solar, we evaluate the mean molecular weight as that of water, and q_KCl as a ratio of K to O because the atmosphere is dominated by water rather than hydrogen for extremely metal-enriched cases. The Appendix summarizes the evaluation of other physical parameters (e.g., viscosity) required for our calculations. Tables 1 and 2 show the values of q_KCl, K₀, and H at the upper isothermal region for GJ1214 b and GJ436 b, respectively.

In this section, we particularly focus on the physical mechanisms that control the vertical distributions of the cloud particle size. Figure 3 shows the calculated vertical profiles of mineral clouds in GJ1214 b. We find that the particles grow mainly near the cloud base (left column in Figure 3) and stop growing in the upper atmosphere where P ≤ 10⁻³ bar. This occurs because the mixing timescale ${\tau }_{\mathrm{mix}}\propto {K}_{{\rm{z}}}^{-1}\propto {P}^{2/5}$ decreases with height, and eventually becomes shorter than the timescales of condensation, coagulation, and coalescence. This trend is also seen in the results of a recent cloud model that takes into account the evolution of particle size distribution (Gao et al. 2018). The final particle radius ranges from 1 to 2 μm for the metallicity of 1× solar, 0.9–4 μm for 10× solar, 1.5–10 μm for 100× solar, and 5–30 μm for water vapor atmosphere, respectively. Figure 3 indicates that the final particle size decreases with the n_CCN and approaches a minimum value in the limit of high n_CCN. In Section 3.2 we explain how the final particle size is determined. We also find that a higher metallicity leads to a larger final size, although its effect is small compared with that of CCN number density.

Zoom In Zoom Out Reset image size

The cloud mass mixing ratio, defined as ρ_c/ρ_g, steeply decreases with height above the height where τ_fall < τ_mix. This can be understood from the transport equations. In the upper atmosphere, the source terms expressing the particle growth are negligible, as mentioned above. Therefore, in a steady state, the vertical mixing of particles should balance with sedimentation,

$$\begin{eqnarray}&&-{\rho }_{{\rm{g}}}{K}_{{\rm{z}}}\displaystyle \frac{\partial }{\partial z}\left(\displaystyle \frac{{\rho }_{{\rm{c}}}}{{\rho }_{{\rm{g}}}}\right)-{v}_{{\rm{t}}}{\rho }_{{\rm{c}}}=0.\end{eqnarray} \tag{ 18 }$$

When the τ_mix ≪ τ_fall, Equation (18) indicates that ρ_c/ρ_g is nearly constant for height, which is seen in the lower region of Figure 3. When the τ_mix ≫ τ_fall, the mass mixing ratio decreases with height due to the particle sedimentation.

The vertical distribution of the cloud mass density also depends on the CCN number density and atmospheric metallicity (the right column of Figure 3). A larger CCN number density leads to a higher mass density at high altitude, because the final particle size decreases with n_CCN as mentioned before. We also find that a higher metallicity yields a higher cloud mass density at high altitude. This metallicity dependence arises because the final particle size is insensitive to the metallicity, while the cloud mass density at the cloud base is approximately proportional to the metallicity. Furthermore, the dependence of τ_mix ∝ H² also yields the higher cloud mass at high altitude for higher metallicity cases, because the H decreases with increasing atmospheric metallicity.

The final particle size determines how high the cloud particles can ascend. Here we discuss the mechanisms that control the final particle size. Figure 4 shows the vertical distributions of the timescales of vertical mixing, condensation, coagulation, and coalescence for three cases: depleted CCNs (n_CCN = 10⁶ m⁻³), enriched CCNs (n_CCN = 10¹² m⁻³), and enriched vapor (steam atmospheres). The timescales of condensation, coagulation, and coalescence are defined as

$$\begin{eqnarray}&&{\tau }_{\mathrm{cond}}={\rho }_{{\rm{c}}}{\left|\displaystyle \frac{\partial {\rho }_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{cond}}^{-1},\end{eqnarray} \tag{ 19 }$$

$$\begin{eqnarray}&&{\tau }_{\mathrm{coag}}={n}_{{\rm{c}}}{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coag}}^{-1},\end{eqnarray} \tag{ 20 }$$

$$\begin{eqnarray}&&{\tau }_{\mathrm{coal}}={n}_{{\rm{c}}}{\left|\displaystyle \frac{\partial {n}_{{\rm{c}}}}{\partial t}\right|}_{\mathrm{coal}}^{-1},\end{eqnarray} \tag{ 21 }$$

and the mixing timescale τ_mix is given by Equation (9). Generally, cloud particles grow if $\min ({\tau }_{\mathrm{cond}},{\tau }_{\mathrm{coag}},{\tau }_{\mathrm{coal}})\,\ll {\tau }_{\mathrm{mix}}$, and ascend without significant growth if $\min ({\tau }_{\mathrm{cond}},{\tau }_{\mathrm{coag}},{\tau }_{\mathrm{coal}})\,\gg \,{\tau }_{\mathrm{mix}}$. The mixing timescale τ_mix decreases with height as mentioned before, whereas the growth timescales increase with height because they are inversely proportional to the density. Hence the particle growth becomes relatively less effective as the particles ascend. In following subsections, we characterize the particle growth in three cases based on timescale argument.

Zoom In Zoom Out Reset image size

3.2.1. Depleted CCN Regime (τ _mix < τ _coag, τ_coal)

In the example shown in the left panel of Figure 4, τ_cond is much shorter than τ_mix and other growth timescales at the cloud base. The short τ_cond results in the quick growth of particles near the cloud base, as shown in Figure 3. At the same time, the rapid condensation also results in rapid depletion of condensing vapor. This depletion eventually suppresses the condensation growth, and hence the total cloud mass at the cloud base is limited by the total amount of condensing vapor there (i.e., ${\rho }_{{\rm{c}}}({z}_{{\rm{b}}})\approx {\rho }_{{\rm{v}}}({z}_{{\rm{b}}})={\rho }_{{\rm{s}}}({z}_{{\rm{b}}})$).

If n_CCN is so small that $\min ({\tau }_{\mathrm{coal}},{\tau }_{\mathrm{coag}})\gt {\tau }_{\mathrm{mix}}$ at the cloud base, the particles start to ascend as soon as the condensation growth is completed (the left column in Figure 4). In this case, the final particle size r_cond is determined by the deposition of available vapor onto CCNs—that is,

$$\begin{eqnarray}&&\displaystyle \frac{4}{3}\pi {r}_{\mathrm{cond}}^{3}{\rho }_{{\rm{p}}}{n}_{\mathrm{CCN}}\approx {\rho }_{{\rm{s}}}({z}_{{\rm{b}}}),\end{eqnarray} \tag{ 22 }$$

and thus

$$\begin{eqnarray}&&{r}_{\mathrm{cond}}\approx {\left[\displaystyle \frac{3{\rho }_{{\rm{s}}}({z}_{{\rm{b}}})}{4\pi {\rho }_{{\rm{p}}}{n}_{\mathrm{CCN}}}\right]}^{1/3},\end{eqnarray} \tag{ 23 }$$

where we have assumed the initial CCN mass density is much smaller than ρ_s(z_b). Figure 5 shows the final particle size and r_cond for each metallicity case. As shown in Figure 5, the final particle size approaches r_cond for lower CCN number density. Hence Equation (23) explains why the final particle size decreases with the increasing of CCN number density. Equation (23) also explains the results of Gao et al. (2018), who found that the efficient homogeneous nucleation (high particle number density) results in small particle size.

Zoom In Zoom Out Reset image size

3.2.2. Enriched CCN Regime (τ _coag < τ _mix < τ _coal)

Coagulation leads the further growth of cloud particles in addition to condensation if n_CCN is so high that τ_coag < τ_mix at the cloud base (see the middle panel of Figure 4). When coagulation is effective, the final particle size becomes larger than r_cond and eventually reaches the minimum value in the limit of high CCN number density, as seen in Figure 3. Hence we can expect the particle size must be larger than the minimum value determined by coagulation, even if the CCN number density is uncertain.

The minimum particle size can be analytically estimated in the following way. Because the final particle size ranges as r > 0.07 μm in Figure 3, the coagulation growth falls into diffusive regime, and the τ_coag is written by

$$\begin{eqnarray}&&{\tau }_{\mathrm{coag}}=\displaystyle \frac{3\eta }{4{k}_{{\rm{B}}}T\beta {n}_{{\rm{c}}}}.\end{eqnarray} \tag{ 24 }$$

Also, the slip factor can be approximated as $\beta \approx {\beta }_{\infty }{\mathrm{Kn}}_{{\rm{g}}}$, where ${\beta }_{\infty }=1.657$, because the mean free path near the cloud base (l ∼ 10 μm) is larger than the particle radius (i.e., Kn_g ≫ 1). Using the relation $4\pi {r}_{{\rm{c}}}^{3}{\rho }_{{\rm{p}}}{n}_{{\rm{c}}}/3={\rho }_{{\rm{c}}}$ and η = ρ_gv_thl/3 (Woitke & Helling 2003), where ${v}_{\mathrm{th}}\,=\sqrt{8{k}_{{\rm{B}}}T/\pi {m}_{{\rm{g}}}}$ is the mean thermal velocity, the coagulation timescale can be rewritten as

$$\begin{eqnarray}{\tau }_{\mathrm{coag}} & = & \displaystyle \frac{{\rho }_{{\rm{g}}}{v}_{\mathrm{th}}{r}_{{\rm{c}}}}{4{k}_{{\rm{B}}}T{\beta }_{\infty }{n}_{{\rm{c}}}}\\ & = & \displaystyle \frac{\pi {\rho }_{{\rm{p}}}{v}_{\mathrm{th}}}{3{k}_{{\rm{B}}}T{\beta }_{\infty }{q}_{{\rm{c}}}}{r}_{{\rm{c}}}^{4},\end{eqnarray} \tag{ 25 }$$

where q_c ≡ ρ_c/ρ_g is the cloud mass mixing ratio. Because ρ_c is determined by saturation vapor density at the cloud base (see Section 3.2.1), ${q}_{{\rm{c}}}={\rho }_{{\rm{s}}}({z}_{{\rm{b}}})/{\rho }_{{\rm{g}}}({z}_{{\rm{b}}})={m}_{\mathrm{KCl}}{q}_{\mathrm{KCl}}/{m}_{{\rm{g}}}$, where m_KCl is the mass of a KCl molecule. Coagulation growth completes when the vertical mixing becomes more efficient, and hence the final size is determined from the condition τ_coag = τ_mix. Equating Equations (25) and (9), the final particle size determined by coagulation r_coag is predicted as

$$\begin{eqnarray}&&{r}_{\mathrm{coag}}={\left(\displaystyle \frac{3{\beta }_{\infty }}{\sqrt{8\pi }}\displaystyle \frac{{m}_{\mathrm{KCl}}{q}_{\mathrm{KCl}}}{{\rho }_{{\rm{p}}}{K}_{{\rm{z}}}({z}_{{\rm{b}}})}{g}^{1/2}{H}^{5/2}\right)}^{1/4}{\left(\displaystyle \frac{{P}_{* }}{{P}_{{\rm{b}}}}\right)}^{1/10},\end{eqnarray} \tag{ 26 }$$

where P_b is the pressure of the cloud base and P_* is the pressure in which the coagulation growth is completed. Equation (26) implies the final particle size in this regime is almost independent of n_CCN, because P_* is insensitive to the choice of n_CCN, as seen in Figure 4. This explains why the final particle size is almost independent of n_CCN for high CCN number density in Figure 3. Particularly, we find that Equation (26) is in good agreement with the minimum final size derived from the numerical results if we assume P_* = 0.1P_b. In this case, Equation (26) can be rewritten as the following useful formula

$$\begin{eqnarray}{r}_{\mathrm{coag}} & = & 1.25\,\mu {\rm{m}}\\ & & \times {\left(\displaystyle \frac{g}{10{\rm{m}}{{\rm{s}}}^{-2}}\right)}^{1/8}{\left(\displaystyle \frac{H}{{10}^{2}\mathrm{km}}\right)}^{5/8}\\ & & \times {\left(\displaystyle \frac{{K}_{{\rm{z}}}({z}_{{\rm{b}}})}{{10}^{3}{{\rm{m}}}^{2}{{\rm{s}}}^{-1}}\right)}^{-1/4}{\left(\displaystyle \frac{{q}_{\mathrm{KCl}}}{{10}^{-5}}\right)}^{1/4}.\ \end{eqnarray} \tag{ 27 }$$

Figure 5 shows that the final particle size asymptotically reaches that predicted from Equation (27), except in the case of steam atmosphere. The deviation for steam atmosphere is caused by coalescence, as explained in next subsection.

3.2.3. Enriched Vapor Regime (τ _coal < τ _mix)

Coalescence is dominant only if condensing vapor is very abundant as in pure steam atmospheres, as shown in the right panel of Figure 4. When coalescence is dominant, the final particle size becomes larger than the lower limit set by coagulation r_coag (see the bottom panel of Figure 5). Because larger particles have larger settling velocity, coalescence suppresses the cloud-top height in the steam atmosphere.

Here we predict the threshold abundance of condensing vapor that induces the significant growth through coalescence. Because the particle size is smaller than the gas mean free path near the cloud base ∼10 μm in most of our calculations, the terminal velocity is expressed as the Epstein's law, approximated as

$$\begin{eqnarray}&&{v}_{{\rm{t}}}({r}_{{\rm{c}}})\approx \displaystyle \frac{2{\beta }_{\infty }g{\rho }_{{\rm{p}}}}{3{\rho }_{{\rm{g}}}{v}_{\mathrm{th}}}{r}_{{\rm{c}}}.\end{eqnarray} \tag{ 28 }$$

Therefore the coalescence timescale can be rewritten as

$$\begin{eqnarray}{\tau }_{\mathrm{coal}} & = & \displaystyle \frac{1}{2\pi {r}^{2}{\rm{\Delta }}{vn}}\\ & = & \displaystyle \frac{2{v}_{\mathrm{th}}}{{\beta }_{\infty }{{gq}}_{{\rm{c}}}}.\end{eqnarray} \tag{ 29 }$$

Because q_c = m_KClq_KCl/m_g (see Section 3.2.2), the coalescence timescale just above the cloud base is independent of n_CCN, and only depends on the mixing ratio of the condensing vapor. If τ_coal ≪ τ_mix, the cloud particles grow via coalescence in addition to condensation and coagulation. Comparing Equation (29) with τ_mix(z_b), we find that coalescence occurs near the cloud base if the condensate mixing ratio is much higher than

$$\begin{eqnarray}{q}_{* } & \approx & \displaystyle \frac{2{v}_{\mathrm{th}}{K}_{{\rm{z}}}{m}_{{\rm{g}}}}{{\beta }_{\infty }{{gH}}^{2}{m}_{\mathrm{KCl}}}\\ & \approx & 5\times {10}^{-7}{\left(\displaystyle \frac{g}{10{\rm{m}}{{\rm{s}}}^{-2}}\right)}^{-1/2}{\left(\displaystyle \frac{H}{{10}^{2}\mathrm{km}}\right)}^{-3/2}\\ & & \times \left(\displaystyle \frac{{K}_{{\rm{z}}}({z}_{{\rm{b}}})}{{10}^{3}\,{{\rm{m}}}^{2}\,{{\rm{s}}}^{-1}}\right)\left(\displaystyle \frac{{m}_{{\rm{g}}}}{2\,\mathrm{amu}}\right).\end{eqnarray} \tag{ 30 }$$

Substituting the parameters for the steam atmosphere of GJ1214 b, the mixing ratio of condensing vapor q_KCl = 2.61 × 10⁻⁴ exceeds the q_* ∼ 2 × 10⁻⁵ by an order of magnitude, and hence coalescence dominates the particle growth.

The prediction of the cloud-top height, defined as the height at which the atmosphere becomes opaque due to the cloud, is important because it determines the shape of observed spectra (e.g., Brown 2001). In order to predict it, we calculate the slant optical depth τ_s, defined as the optical depth for the path length of the transmitted starlight, using Equation (6) of Fortney (2005), given by

$$\begin{eqnarray}&&{\tau }_{{\rm{s}}}={\tau }_{{\rm{v}}}\sqrt{\displaystyle \frac{2\pi {R}_{{\rm{p}}}}{{H}_{{\rm{c}}}}},\end{eqnarray} \tag{ 31 }$$

where τ_v is the vertical optical depth of the cloud, R_p is the planetary radius, and ${H}_{{\rm{c}}}=| d\mathrm{ln}{n}_{{\rm{c}}}/{dz}{| }^{-1}$ is the cloud scale height. The cloud-top height can be estimated as the height at which τ_s exceeds unity. The vertical optical depth τ_v is given by

$$\begin{eqnarray}&&{\tau }_{{\rm{v}}}(z)={\int }_{z}^{\infty }{Q}_{\mathrm{ext}}({r}_{{\rm{c}}})\pi {r}_{{\rm{c}}}^{2}{n}_{{\rm{c}}}{dz},\end{eqnarray} \tag{ 32 }$$

where Q_ext is the extinction coefficient of the particles. To calculate Q_ext, we perform rigorous Mie calculations using BHMIE code (Bohren & Huffman 1983). We use the refractive index of KCl from Querry (1987) and assume an isotropic scattering for the calculations of scattering opacities. We assume the wavelength of λ = 1.4 μm, at which a prominent water feature is located. For GJ1214 b and GJ436 b, this feature is absent in the actual spectra (Knutson et al. 2014a; Kreidberg et al. 2014), and therefore the cloud-top height defined at 1.4 μm must be sufficiently high so that the clouds fully obscure the feature.

To estimate H_c, we use the fact that at high altitudes particle growth is negligible, and hence the vertical profiles are given by Equation (18). Since $H=| d\mathrm{ln}{n}_{{\rm{g}}}/{dz}{| }^{-1}$ and ${H}_{{\rm{c}}}\,=| d\mathrm{ln}{n}_{{\rm{c}}}/{dz}{| }^{-1}$, Equation (18) can be rewritten as

$$\begin{eqnarray}&&{H}_{{\rm{c}}}=H{\left(1+\displaystyle \frac{{v}_{{\rm{t}}}H}{{K}_{{\rm{z}}}}\right)}^{-1}.\end{eqnarray} \tag{ 33 }$$

Because the term v_tH/K_z is the ratio of τ_mix to τ_fall, H_c ≈ H for τ_mix ≪ τ_fall, and H_c ≈ K_z/v_t for τ_mix ≫ τ_fall.

Figure 6 shows the cloud-top pressure for GJ1214 b predicted by our calculations for different values of the metallicity and n_CCN. We find that the cloud-top height increases with n_CCN, but plateaus in the limit of high n_CCN. This means that one can predict the maximum height of the cloud top for given atmospheric metallicity. The presence of the maximum height results from the presence of the minimum particle size mentioned in Section 3.2.2. The impact of size distribution on the predicted maximum height is small, as we will discuss in Section 5.1.

Zoom In Zoom Out Reset image size

We also find that metal-enriched atmospheres are more likely to yield vertically extended clouds for the abundance of condensing vapor below the threshold (Equation (30)). Figure 6 shows that the cloud-top is placed at P ≳ 3 × 10⁻³ bar for 1× solar metallicity, P ≳ 3 × 10⁻⁴ bar for 10× solar metallicity, P ≳ 1 × 10⁻⁴ bar for 100× solar metallicity, and P ≳ 6 ×10⁻⁴ bar for the steam atmosphere. This trend arises because a higher metallicity atmosphere yields a higher total cloud mass and more efficient vertical mixing as mentioned before. For vapor abundance above the threshold, the case of the steam atmosphere, cloud-top height no longer increases with metallicities because coalescence leads significant growth for cloud particles, as mentioned in Section 3.2.3.

In the left panel of Figure 7, we plot the height (in pressure) of the cloud base and the maximum height of the cloud top predicted for GJ1214 b as a function of the metal mass fraction of the atmosphere.³ The blue-shaded area in Figure 7 thus indicates the heights where the top of the KCl cloud can exist for some CCN number density. For comparison, we also indicate by the gray shading the heights of the cloud top suggested by Kreidberg et al. (2014), based on Bayesian analysis on the observed transmission spectrum.

Zoom In Zoom Out Reset image size

We find that the maximum cloud-top height is too low to explain the flat spectrum for all plausible values of the atmospheric metallicity. In principle, a higher atmospheric metallicity provides a higher cloud-top height, as already mentioned in Section 3.3. However, even if we assume the steam atmosphere, the maximum cloud-top height (P_min = 6 × 10⁻⁴ bar in pressure) is still an order of magnitude higher than inferred by Kreidberg et al. (2014; cloud-top pressure ≤3 × 10⁻⁵ bar at 3σ confidence). This is because in the steam atmosphere, coalescence causes the significant growth of cloud particles of r_c ≥ 5 μm.

The above comparison clearly shows that a simple condensate cloud cannot explain the flat transmission spectrum of GJ1214 b. This fact might support the idea that the flat spectrum of GJ1214 b is caused by photochemical haze (Miller-Ricci Kempton et al. 2012; Morley et al. 2013, 2015; Kawashima & Ikoma 2018), rather than by mineral clouds. Alternatively, our cloud model might still be missing important physics of particle growth. For example, it is suggested both theoretically and experimentally (e.g., Dominik & Tielens 1997; Blum & Wurm 2000) that, unlike water cloud droplets, solid particles grow into highly porous particles through mutual sticking. This porosity evolution is neglected in Figure 7, but could help particles ascend to very high altitudes because porous particles have a lower settling velocity than compact particles of the same mass. We address this possibility in Section 5.3.

For GJ436 b, we find that the maximum cloud-top height for KCl clouds is high enough to be consistent with the transmission observations. The Bayesian analysis by Knutson et al. (2014a) indicates that the cloud top is present at atmospheric pressures of ≲10⁻² bar, except for metal-rich atmospheres of metal mass fraction ≳0.8, for which the location of the cloud top is not well constrained (see the gray-shaded area in the right panel of Figure 7). As shown in the right panel of Figure 7, the minimum cloud-top pressure P_min predicted from our model is 2 × 10⁻³, 3 × 10⁻⁴, 8 × 10⁻⁵, and 5 × 10⁻⁶ bar for the metallicities of 1×, 10×, 100×, and 1000× solar, respectively. Since we adopted the high K_z for the metallicity of 1000× solar compared with that for the steam atmosphere on GJ1214 b, cloud particles avoid the significant growth due to coalescence. This is a reason why the cloud-top height for 1000× solar is much higher than that for the steam atmospheres on GJ1214 b. Combining the Bayesian analysis results and our model prediction, we suggest that the flat spectrum of GJ436 b is likely caused by a KCl cloud with its top at ∼ 10⁻³–10⁻² bar for hydrogen-rich atmospheres (metal mass fraction ≲0.8) and at ∼10⁻²–10⁻⁵ bar for metal-rich atmospheres (metal mass fraction ≳0.8). However, because we here adopted a metallicity-independent eddy diffusion coefficient (see Section 2.3), we cannot conclude whether the atmosphere of GJ436 b is likely to be hydrogen-rich or metal-rich. Future three-dimensional modeling of GJ436 b's atmospheric circulation, like the one done by Charnay et al. (2015a) for GJ1214 b, would allow us to determine the atmosphere's metallicity.

Because the total particle cross section tends to be dominated by small particles rather than by the particles dominating the total cloud mass, the cloud-top height might be influenced by the size distribution, which is, however, not captured by our calculations. We here evaluate the impact of particle size distribution on the predicted cloud-top height by adding to our model a distribution of small particles. We assume that particles smaller than r_c obey the Hansen size distribution (Hansen 1971) given by

$$\begin{eqnarray}&&f(r)\equiv \displaystyle \frac{{dn}}{{dr}}={{Cr}}^{(1-3b)/b}\exp \left(-\displaystyle \frac{r}{{ab}}\right),\end{eqnarray} \tag{ 34 }$$

where f(r)dr is the number density of particles with radii between r and r + dr, C is a constant, a and b are the mean effective radius, and the effective variance defined by

$$\begin{eqnarray}&&a\equiv \displaystyle \frac{{\int }_{0}^{\infty }r\pi {r}^{2}f(r){dr}}{{\int }_{0}^{\infty }\pi {r}^{2}f(r){dr}},\end{eqnarray} \tag{ 35 }$$

$$\begin{eqnarray}&&b\equiv \displaystyle \frac{{\int }_{0}^{\infty }{(r-a)}^{2}\pi {r}^{2}f(r){dr}}{{a}^{2}{\int }_{0}^{\infty }\pi {r}^{2}f(r){dr}}.\end{eqnarray} \tag{ 36 }$$

The Hansen size distribution successfully reproduces the observed size distributions of terrestrial water clouds for b = 0.1–0.2 (Hansen 1971), and near-infrared spectral energy distributions of cloudy brown dwarfs for b > 0.5 (Hiranaka et al. 2016). The top panel of Figure 8 shows the Hansen size distributions for b = 0.2, 0.5, and 1.0. One can see that b < 0.5 yields log-normal-like size distributions, while b > 0.5 yields power-law-like size distributions. Therefore, the Hansen size distribution with various choices of b enables us to test the size distributions of various shapes.

Zoom In Zoom Out Reset image size

For each height, we determine the a and C so that the mass weighted size and the cloud mass density correspond to r_c and ρ_c calculated by our model, respectively. We calculate the a and C at each height using the following relations:

$$\begin{eqnarray}&&{r}_{{\rm{c}}}=\displaystyle \frac{{\int }_{0}^{\infty }{rm}(r)f(r){dr}}{{\int }_{0}^{\infty }m(r)f(r){dr}}=a(1+b),\end{eqnarray} \tag{ 37 }$$

$$\begin{eqnarray}&&{\rho }_{{\rm{c}}}=\displaystyle \frac{4\pi {\rho }_{{\rm{p}}}}{3}C{({ab})}^{(1+b)/b}{\rm{\Gamma }}\left(\displaystyle \frac{1+b}{b}\right),\end{eqnarray} \tag{ 38 }$$

where Γ(z) is the gamma function. However, r_c might be overestimated in our calculations at the height where τ_mix > τ_fall(r_c), because our model fails to trace the decreasing of mass weighted size due to the removal of large particles by gravitational settling. To avoid this issue, we use an analytical solution of the transport equation. In the upper atmosphere, the particle growth is negligible, as mentioned in Section 3, and hence, in a steady state, the particle number density obeys

$$\begin{eqnarray}&&-{n}_{{\rm{g}}}{K}_{{\rm{z}}}\displaystyle \frac{\partial }{\partial z}\left(\displaystyle \frac{{n}_{{\rm{c}}}}{{n}_{{\rm{g}}}}\right)-{v}_{{\rm{t}}}{n}_{{\rm{c}}}=0.\end{eqnarray} \tag{ 39 }$$

If we approximate $\beta =1+{\beta }_{\infty }{\mathrm{Kn}}_{{\rm{g}}}$,⁴ Equation (39) can be analytically solved as

$$\begin{eqnarray}{n}_{{\rm{c}}}(P) & = & {n}_{{\rm{c}}}({P}_{0}){ \mathcal P }\\ & & \times \exp \left[\displaystyle \frac{5\chi ({P}_{0})}{2(1+{\beta }_{\infty }{\mathrm{Kn}}_{0})}\right.\\ & & \times \left.\left(\left({{ \mathcal P }}^{2/5},-,1,),-,\displaystyle \frac{2{\beta }_{\infty }{\mathrm{Kn}}_{0}}{3},(,{{ \mathcal P }}^{-3/5},-,1\right)\right)\right],\end{eqnarray} \tag{ 40 }$$

where ${ \mathcal P }\equiv P/{P}_{0}$ and χ(P) ≡ τ_mix(P)/τ_fall(P) is the ratio of the mixing timescale to the falling timescale. We calculate the size distributions at the regions of τ_mix < τ_fall(r_c) using Equation (40) for each size bin. Figure 8 shows the constructed vertical size distributions for b = 0.2, 0.5, and 1.0. The reference pressure P₀ is set as a height where τ_mix = τ_fall(r_c), denoted as the white dotted lines in each panel. Figure 8 indicates that the larger b is, the more small particles are present at high altitude.

Figure 9 compares the cloud-top heights predicted by the characteristic size model with those by the model with the Hansen size distribution of b = 0.2, 0.5, and 1.0. With particle size distribution, the vertical optical depth is calculated as

$$\begin{eqnarray}&&{\tau }_{{\rm{v}}}(z)={\int }_{z}^{\infty }{\int }_{0}^{\infty }{Q}_{\mathrm{ext}}(r)\pi {r}^{2}f(r,z){drdz}.\end{eqnarray} \tag{ 41 }$$

We find that size distribution has little effect on the cloud top height, except for the case of the steam atmosphere. This is because KCl is a purely scattering material in near-infrared (i.e., the extinction is equivalent to the scattering). For purely scattering particles smaller than the wavelength, the extinction efficiencies steeply decrease with decreasing the particle size as Q_ext ∝ r⁴. Therefore, the contribution of such small particles (r ≪ λ) to the total cloud opacity is negligibly small. The most efficient extinction occurs at r ∼ λ/2π, which is 0.2 μm for λ = 1.4 μm. By contrast, the final characteristic sizes for the metallicities of 1×, 10×, and 100× solar are r_c ≈ 1 μm (see Figure 5), already close to 0.2 μm. Therefore, the addition of particles smaller than r_c has little effect on the optical depth and hence on the cloud-top height.

Zoom In Zoom Out Reset image size

The difference arising from size distributions becomes obvious only when the final characteristic size is orders of magnitude larger than r = λ/2π. This is the case for the steam atmospheres, in which the final characteristic size is r ≈ 5 μm ≫ 0.2 μm. In this case, varying the size distribution can decrease the cloud-top pressure by a factor of 3 from the prediction of the characteristic size method. However, we find that the cloud-top height for GJ1214 b with a steam atmosphere is still an order of magnitude lower than anticipated from the observation of Kreidberg et al. (2014). Thus we conclude that one cannot explain the flat spectrum of GJ1214 b solely by considering particle size distribution.

Our P–T structure neglects heat transport by convection, which is the process called "convective adjustment," and in reality becomes important when the temperature steeply declines with decreasing pressure (Manabe & Strickler 1964; Marley & Robinson 2015). We performed test calculations using the P–T structure provided by Charnay et al. (2015a), which includes the effect of the convective adjustment. As listed in Table 3, we confirmed that the maximum heights of cloud top are quantitatively similar to the results from radiative P–T profiles of Guillot (2010). The largest influence is only a factor of 2, which occurs for the case of the steam atmospheres. The reason why the cloud-top height is nearly the same for both P–T structures is that the convective adjustment only changes the cloud-base height slightly. Since the minimum particle size is not sensitive to the cloud-base height (r_coag ∝ P_b^1/10 from Equation (27)), the convective adjustment has little effect on the predicted cloud-top height.

In the calculations presented in Sections 3 and 4, we assumed that cloud particles are compact and their internal density is constant. This assumption would be valid for liquid droplets, but breaks down if solid KCl cloud particles grow into porous aggregates. As pointed out by Marley et al. (2013), porous aggregates are easily lofted to high altitude because they have large cross sections as compared with compact particles of the same mass. Therefore, the predicted cloud-top height could be influenced by particle porosity.

Here we quantify the impacts of particle porosity on the predicted cloud-top height. We introduce the volume filling factor ϕ defined by

$$\begin{eqnarray}&&\phi \equiv \displaystyle \frac{{\rho }_{\mathrm{int}}}{{\rho }_{{\rm{p}}}}.\end{eqnarray} \tag{ 42 }$$

The volume filling factor takes ϕ = 1 for compact particles and ϕ < 1 for porous particles. Snowflakes in the Earth are known to have ϕ = 0.5–0.005 (Magono & Nakamura 1965), while grains in protoplanetary disks could have an extremely low filling factor of ϕ ∼ 10⁻⁴ according to recent theoretical studies (Okuzumi et al. 2012; Kataoka et al. 2013). We repeat the calculations presented in Section 3.3 by replacing ρ_p as ρ_int and varying ϕ from ϕ = 1 to ϕ = 0.01.

To evaluate Q_ext for porous aggregates, we calculate the effective refractive index using the effective medium theory (EMT) with the Maxwell-Garnett mixing rule (Bohren & Huffman 1983). The EMT provides reasonable estimates for aggregate's absorption and scattering opacities when the particles that constitute the aggregates are smaller than the incident wavelength (Voshchinnikov et al. 2007; Shen et al. 2008).

In Figure 10, we show the final characteristic size and cloud-top height for various values of ϕ and atmospheric metallicities. We find that the final characteristic size increases with decreasing ϕ. This is because, in the limit of high N_CCN, the final characteristic size is proportional to ϕ^−1/4 as indicated by Equation (26). Since v_t ∝ ρ_pr_c ∝ ϕ^3/4 in upper atmospheres (see Equation (28)), porous aggregates are indeed easily lofted to high altitude as compared with compact particles.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

However, the cloud-top height does not appreciably increase with decreasing ϕ except for the steam atmosphere (see each bottom panel of Figure 10). When the particle porosity is taken into account, the maximum cloud-top height is ∼2 × 10⁻³ bar for the metallicity of 1×  solar, ∼2 × 10⁻⁴ bar for 10×  solar, and ∼3 × 10⁻⁵ bar for 100×  solar, respectively, which are only higher than those for compact particles by a factor of 2–3. The reason why the cloud-top height is insensitive to ϕ comes from the optical properties of porous aggregates, shown in Figure 11. The scattering mass opacity of a porous aggregate is proportional to ϕ as long as ϕ < λ/r, in which the aggregate itself becomes optically thin (Kataoka et al. 2014). The two effects of reducing opacities and increasing the cloud amounts at high altitude with decreasing ϕ largely cancel out, explaining why the impacts of particle porosities are not drastic for hydrogen-rich atmospheres.

By contrast, for the steam atmosphere, the cloud-top height for porous aggregates can be much higher than that for compact particles. The maximum cloud-top height is ∼1 × 10⁻⁴ bar for ϕ = 0.1 and ∼2 × 10⁻⁵ bar for ϕ = 0.01. The distinct increase in cloud-top height is caused by the efficient growth via coalescence. Because τ_coal only depends on the cloud mass mixing ratio (see Equation (29)), the efficient growth via coalescence occurs even for porous aggregates. Coalescence produces particles large enough to have a high scattering opacity, and hence the cancellation due to the effect of reducing opacities with decreasing ϕ does not occur appreciably. This is a reason why the particle porosity drastically increases the cloud-top height for the steam atmosphere.

Figure 12 shows the maximum extent of KCl clouds for GJ1214 b from the models that take into account size distribution (Section 5.1) and particle porosity. We find that the cloud-top height is still too low to be consistent with the flat spectrum of Kreidberg et al. (2014) for hydrogen-rich atmospheres (1×, 10×, and 100× solar metallicities). On the other hand, for the steam atmosphere with ϕ = 0.01, we find that KCl clouds can reach ≈1.5 × 10⁻⁵ bar, which is equivalent to within the 3σ Bayesian credible regions of the cloud-top height (P ≈ 3 × 10⁻⁵ bar) reported by Kreidberg et al. (2014). Since particle porosity naturally increases through coalescence, high-altitude cloud formation in the steam atmospheres, where coalescence is effective, might be a plausible explanation for the flat spectrum of GJ1214 b.

Zoom In Zoom Out Reset image size

We note that the estimates for the cloud-top height given above are based on the assumption of isotropic scattering. The forward scattering of cloud particles potentially reduces the effective cloud opacity (de Kok & Stam 2012; Robinson 2017) and hence produces a lower cloud top. This effect cannot be captured here correctly because EMT tends to overestimate the degree of forward scattering of porous aggregates (Shen et al. 2009; Tazaki et al. 2016; Tazaki & Tanaka 2018). The angular dependent properties of scattered light depend on the microstructure of an aggregate. Further understanding about the microstructure and optical property of aggregates is required to verify the possibility of high-altitude cloud formation by porous aggregates.