The Astrophysical Journal Letters, 780:L8 (5pp), 2014 January 1 C 2014.
doi:10.1088/2041-8205/780/1/L8
The American Astronomical Society. All rights reserved. Printed in the U.S.A.
ON THE VIGOR OF MANTLE CONVECTION IN SUPER-EARTHS 1
Takehiro Miyagoshi1 , Chihiro Tachinami2 , Masanori Kameyama3 , and Masaki Ogawa4
Institute for Research on Earth Evolution, Japan Agency for Marine-Earth Science and Technology, 3173-25 Showa-machi, Kanazawa-ku, Yokohama 236-0001, Japan;
[email protected] 2 Department of Earth and Planetary Sciences, Tokyo Institute of Technology, Meguro, Tokyo 152-8551, Japan;
[email protected] 3 Geodynamics Research Center, Ehime University, 2-5 Bunkyo-cho, Matsuyama, Ehime 790-8577, Japan;
[email protected] 4 Department of Earth Sciences and Astronomy, University of Tokyo at Komaba, 3-8-1 Komaba, Meguro, Tokyo 153-8902, Japan;
[email protected] Received 2013 September 2; accepted 2013 November 20; published 2013 December 10
ABSTRACT Numerical models are presented to clarify how adiabatic compression affects thermal convection in the mantle of super-Earths ten times the Earth’s mass. The viscosity strongly depends on temperature, and the Rayleigh number is much higher than that of the Earth’s mantle. The strong effect of adiabatic compression reduces the activity of mantle convection; hot plumes ascending from the bottom of the mantle lose their thermal buoyancy in the middle of the mantle owing to adiabatic decompression, and do not reach the surface. A thick lithosphere, as thick as 0.1 times the depth of the mantle, develops along the surface boundary, and the efficiency of convective heat transport measured by the Nusselt number is reduced by a factor of about four compared with the Nusselt number for thermal convection of incompressible fluid. The strong effect of adiabatic decompression is likely to inhibit hot spot volcanism on the surface and is also likely to affect the thermal history of the mantle, and hence, the generation of magnetic field in super-Earths. Key words: convection – Earth – methods: numerical – planets and satellites: interiors – planets and satellites: terrestrial planets adiabatic compression (e.g., Jarvis & McKenzie 1980), phaserelationship (e.g., Tackley et al. 2013), and other properties of mantle materials (e.g., Stamenkovi´c et al. 2011) as well as plate tectonics. This means that a simple extrapolation may make it difficult to separate these effects from each other when all of them come together in one numerical model. A step-by-step approach is crucial to develop reliable images on the style of mantle convection in super-Earths. Because of the high hydrostatic pressure (∼TPa; Valencia & O’Connell 2009), mantle convection is most likely to be strongly affected by adiabatic compression of mantle materials in super-Earths. Jarvis & McKenzie (1980) carried out linear and nonlinear analyses of a thermal convection of moderately compressible fluids with constant viscosity in two-dimensional space and found that density stratification by the adiabatic compression enhances the convective stability when the depth of the fluid layer far exceeds the thermal scale height. Their findings are confirmed later for convection in a three-dimensional spherical shell (Bercovici et al. 1992) and as well as in the presence of much stronger compression (Liu & Zhong 2013). However, this stabilizing effect of adiabatic compression is not noticeable in a more sophisticated model of convection with variable viscosity (Tackley et al. 2013). It is still an open issue to clarify how important the stabilizing effect is in the mantle of super-Earths. In this Letter, we present numerical models of thermal convection of highly compressible fluid with strongly temperaturedependent viscosity and depth-dependent thermal expansivity at a high Rayleigh number expected in super-Earths. We concentrate on clarifying how adiabatic compression affects thermal convection beneath the lithosphere, which develops along the surface boundary owing to the strong temperature-dependence of viscosity, leaving other complications for future works.
1. INTRODUCTION Super-Earths are extra-solar planets which have high mean density (larger than 5000 kg m−3 ; Howard et al. 2010) and small masses (up to ten times the Earth’s). Since the first discovery of the extra-solar planet (Mayor & Queloz 1995), a large number of super-Earths have been detected (Borucki et al. 2011). The high mean density suggests that these planets contain a silicate mantle and iron core like the Earth (Valencia et al. 2010). Inspired by these discoveries, many studies have been done to clarify the internal structure (e.g., Seager et al. 2007; Valencia et al. 2007b; Sotin et al. 2007; Wagner et al. 2011) and the thermal evolution (e.g., Papuc & Davies 2008; Kite et al. 2009; Gaidos et al. 2010; Tachinami et al. 2011; Stamenkovi´c et al. 2012) of these planets. Here, we present numerical models of thermal convection of compressible fluid with strongly temperaturedependent viscosity to discuss how mantle convection operates in the mantle of super-Earths. Many earlier numerical studies on mantle dynamics in superEarths focus on the possible plate tectonics on those planets. Earlier studies, either by simple scaling or by numerical simulations (Valencia et al. 2007a; Valencia & O’connell 2009; Korenaga 2010; van Heck & Tackley 2011), conclude that plate tectonics is likely to occur in super-Earths. Valencia et al. (2007a) suggests that plates become thinner as super-Earths become larger, and hence, that plate tectonics occurs more easily. Korenaga (2010) and van Heck & Tackley (2011) suggest that the presence of water which works as a lubricant at plate boundaries is more important than the planetary size to inducing plate tectonics. Implicit in these pioneering studies is an assumption that mantle convection takes place in super-Earths. Indeed, these studies are based on a simple extrapolation to the relevant conditions from the models of Earth’s mantle convection. However, given the complications involved in the existing numerical models, it is not straightforward but rather dangerous to simply extrapolate these models to super-Earths. The large size of super-Earths can influence mantle convection in many ways through its effects on
2. MODEL Mantle convection of a super-Earth with ten times the Earth’s mass is modeled by a thermal convection of compressible 1
The Astrophysical Journal Letters, 780:L8 (5pp), 2014 January 1
Miyagoshi et al.
when M = 10 M⊕ . The parameter values adopted in our numerical experiments are, Di = 5, Ra = 109 or 1010 , Ts = 0.1, and r = 104 . Notice that Di is defined with α0 . The average value, D¯ i , defined by 1 Di = Di α/α0 dz, (2)
z/d 1 0.8
0
0.6
is 0.79 and is close to the value assumed in Tackley et al. (2013). The basic equations are discretized by a finite difference method and are solved by the ACuTEMAN numerical code (Kameyama et al. 2005; Kameyama 2005; King et al. 2010). We continued the calculations until the convection reached its statistically steady state. The employed mesh is uniform with a grid number of 1024 (horizontal) and 256 (vertical).
0.4 0.2 0 0
0.2
0.4
0.6
0.8
1
α/α0
3. RESULTS
Figure 1. Thermal expansivity plotted against height z/d, where d is the depth of the mantle; z/d = 0 is the bottom of the mantle; and z/d = 1 is the surface of the mantle.
Figure 2 shows our result calculated at Ra = 1010 and r = 104 . (The results are presented in non-dimensional form hereinafter.) Figure 2(a) shows the distribution of potential temperature, Tp (x, z), derived from the temperature distribution T(x, z) shown in Figure 2(b) as 1 Tp (x, z) = T (x, z) exp −Di α(ξ )/α0 dξ . (3)
fluid with an infinite Prandtl number and strongly temperaturedependent viscosity in a two-dimensional rectangular box. The basic equations are numerically solved under the truncated anelastic liquid approximation (Jarvis & McKenzie 1980). There is no internal heat source. The temperature is fixed at Tb and Ts on the bottom and the surface boundaries, respectively. The sidewalls are adiabatic. All boundaries are shear stress free and impermeable. The employed basic equations are the same as those of Jarvis & McKenzie (1980) and Tachinami et al. (2013) except that the viscosity is temperature-dependent. The reference density depends on depth as in Tachinami et al. (2013); the density linearly increases with depth by a factor of about three across the mantle (Valencia et al. 2006). The thermal expansivity, α, depends on depth as shown in Figure 1; α decreases with depth to about 1/10 of its surface value α0 = 3 × 10−5 (K−1 ) on the bottom of the mantle. This profile of α is estimated for MgO from an earlier ab-initio calculation (Tsuchiya & Tsuchiya 2011) and is consistent with earlier laboratory measurements of α under high pressure (Chopelas & Boehler 1992) as well as the profile assumed in earlier models of mantle convection (Tackley et al. 2013). The reduction of α in deep mantle is crucial for convective instability to take place in super-Earths (Tachinami et al. 2013). In contrast, we assume that thermal diffusivity, κ, is constant for simplicity; the assumption of constant κ employed here may yield a rapid cooling of the mantle and core compared to the cases with variable κ (van den Berg et al. 2005). The viscosity, η, depends on temperature, T, as η = η0 exp[E(Tb − T )], (1)
z
The overall structure of the convective flow pattern is dominated by the lithosphere that develops as a thermal boundary layer (TBL) of cold and highly viscous fluid along the surface boundary. The lithosphere laterally moves along the surface and then sinks to the bottom boundary as plumes. The potential temperature contrast across the lithosphere is 0.37, about 93% of that across the mantle; the potential temperature of the isothermal core is higher than that expected for convection with constant viscosity. The potential temperature in the upper part of the isothermal core is slightly higher than that of the lower part, and the isothermal core is convectively stable. In addition to the lithosphere, another TBL develops along the hot bottom boundary. The temperature contrast across the bottom TBL is smaller than that across the lithosphere because of the high temperature in the isothermal core (Figure 2(c)). Small hot plumes ascend from the bottom TBL, as we will discuss further below. The profile of the vertical gradient of horizontally averaged temperature, −∂T /∂z, plotted in Figure 2(d) is almost equal to that of the adiabatic temperature gradient, −∂T /∂zad = Di · α(z) · T , except in the TBLs along the surface and the bottom boundaries. Note, however, that −∂T /∂z is slightly less than −∂T /∂zad just above the bottom TBL (z 0.1), reflecting the thermal structure of the isothermal core. The most prominent feature of the convection shown in Figure 2 is that hot plumes from the bottom boundary stop ascending in the middle of the mantle, and do not penetrate further upward into the shallower part of the mantle. The velocity field in Figure 2(b) clearly shows this feature; the velocity vector is almost horizontal at depths less than about 0.6 except around the cold sinking plumes on the right hand side of the figure. (Notice that the hot region around the left side boundary in Figure 2(a) is induced by viscous dissipation due to the cold descending plume along the boundary, and that it does not induce upwelling flow as shown in Figure 2(b)). Hot plumes stop ascending, because they lose their buoyancy as they ascend owing to the temperature decrease due to adiabatic decompression. The potential temperature Tp in the plume indicated in Figure 2(a), for example, approaches the ambient
where η0 and E are constants and η0 is the viscosity at T = Tb . When converted into its non-dimensional form, four free parameters arise in the present formulation: the surface temperature Ts normalized by the temperature difference across the mantle, ΔT = Tb − Ts ; the dissipation number, Di = α0 gd/Cp ; the Rayleigh number, Ra = ρ0 α0 gΔT d 3 /η0 κ; and the viscosity contrast across the mantle, r = exp(EΔT ). Here, g and ρ 0 are the gravity and the density on the surface boundary, respectively, d is the depth of the mantle, and Cp = 103 (J kg−1 K−1 ) is the specific heat. According to Valencia et al. (2006), g and d depend on the planetary mass, M, as g = g⊕ (M/M⊕ )0.5 and d = d⊕ (M/M⊕ )0.28 in an Earth-like super-Earth, where the subscript ⊕ stands for the Earth. With these relations, g⊕ = 9.8 (m s−2 ), and d⊕ = 2900 (km), Di becomes about 5 2
The Astrophysical Journal Letters, 780:L8 (5pp), 2014 January 1
Miyagoshi et al.
z
(c) (a)
1
Tp(x,z)
0.8 0.6 0.4 0.2 0 0.1
(d)
0.2
0.3
0.4
Tp
z 1
0.8 0.450
0.485
0.6 0.4
(b)
T(x,z) and velocity
0.2 0 0
0.5
1.0
1.5
2.0
− ∂ T ∂z , − ∂T ∂z ad (e)
z
0.1
0.3
1.1
0
0.47
0.49
0.51
Tp
Figure 2. Numerical model calculated at Ra = 1010 and r = 104 . (a) A snapshot of potential temperature distribution. (b) The corresponding snapshot of temperature (colors) and velocity fields (arrows). The length of arrows is proportional to the magnitude of the velocity. (c) The profile of the horizontally averaged potential temperature Tp as a function of height z. (d) The profiles of horizontally averaged temperature gradient −∂T /∂z (red) and the adiabatic temperature gradient −∂T /∂zad = Di · α(z) · T (blue) as a function of height z. (e) Crosses show the potential temperature at the diamond points in the inset, which is a zoom-in figure of the plume indicated by the dashed circle in panel (a), plotted against z. Tp is also shown by the dashed black line.
is reduced at depth and makes the buoyancy of the hot plumes smaller than expected from the nominal value of Ra mentioned above by an order of magnitude; the effective Rayleigh number defined with the α at the bottom boundary is ∼109 . At this Rayleigh number, however, hot plumes still ascend readily to the surface boundary under the Boussinesq approximation where the effect of adiabatic cooling is neglected (e.g., Yuen et al. 1993). The global circulation induced by the sluggish movement of the lithosphere may also affect the growth and ascent of hot plumes to some extent. The global circulation is however, significant only around the cold plumes (see Figure 2(b)), and hardly affects the growth of the hot plumes away from the cold plumes; those remote hot plumes still stop ascending on the way to the surface boundary (see, for example, the plumes around the one indicated by the dashed circle in Figure 2(a)). A similar convective flow pattern is obtained at a lower Rayleigh number of 109 (Figure 3). The lower Ra , however, does affect the detail of the convective flow pattern. The ascending plume around the left sidewall reaches a shallower level of z ∼ 0.8 (Figures 3(a) and (b)), much higher than hot plumes do in Figure 2. Furthermore, the lithosphere (marked by a super-adiabatic vertical temperature profile) is thicker, and the stratosphere-like layer is much thinner than that observed in Figure 2, as can be seen from Figure 3(d). There are almost no regions in the isothermal core where the −∂T /∂z is equal to the −∂T /∂zad ; the thickness of the stratosphere-like layer depends on the Rayleigh number. Figure 4 shows that adiabatic compression strongly affects not only the convective flow pattern but also the efficiency of convective heat transport measured by the Nusselt number, N u = −∂T /∂z|z=1 ; Nu is the ratio of the actual surface heat flux to the heat flux expected when heat is transported across the mantle all by thermal conduction. Under the Boussinesq
potential temperature as the plume ascends (Figure 2(e)). The present result is in a stark contrast to that in our earlier models with constant viscosity (Tachinami et al. 2013), where hot plumes readily ascend to the surface boundary. The difference arises because the lithosphere observed in Figure 2 makes the ambient temperature higher than that for the isoviscous convection. The higher ambient temperature makes the excess temperature of hot plumes smaller, and hence, the effect of adiabatic cooling of hot ascending plumes more prominent; the lithosphere also plays an important role in stopping hot ascending plumes. Because of the absence of hot plumes, the shallower part of the isothermal core looks like a stratosphere in spite of the cold plume that penetrates downward through the mantle. In contrast to hot ascending plumes, the cold plume sinks through almost the entire mantle (Figures 2(a) and (b)). The velocity field shows that the downward flow of the cold plumes is much stronger than the upward flow of hot plumes and does not significantly decay as the plume goes down to the bottom. One reason for this difference is that the effect of adiabatic compression depends on the temperature of plumes. Adiabatic heating in descending cold plumes is weaker than adiabatic cooling in hot ascending plumes, because of the lower temperature of cold plumes. Thus, cold plumes do not lose so much negative buoyancy as hot plumes lose their positive buoyancy. The larger temperature contrast between the cold plume and the surrounding isothermal core and the higher viscosity of cold plumes due to the temperature-dependence of viscosity are also important reasons for the difference. A comment is necessary here on the reason why hot plumes stop ascending on the way to the surface boundary in Figure 2. Besides the adiabatic cooling mentioned above, the depthdependent α might impede the ascent of hot plumes. The α 3
The Astrophysical Journal Letters, 780:L8 (5pp), 2014 January 1
Miyagoshi et al.
(a) Tp(x,z)
z
(c)
1
0.8 0.6 0.4 0.2 0 0.1 0.450
0.2
(b) T(x,z) and velocity
0.3
0.4
0.5
Tp
0.485
(d)
z
1 0.8 0.6 0.4 0.2
0 -0.5 0.3
1.1
0
0.5
1.0
1.5
− ∂ T ∂z , − ∂T ∂z ad
Figure 3. Same as Figures 2(a)–(d), but for the numerical model calculated at Ra = 109 and r = 104 .
on the way to the surface boundary, though cold plumes sink to the bottom boundary through the whole mantle. A stratospherelike layer, where the convective flow is induced only by the drag force from moving lithosphere, develops in the shallower part of the mantle. The efficiency of convective heat transport is lowered below the value predicted from models of mantle convection under the Boussinesq approximation by a factor of around four, and the lithosphere becomes thicker accordingly. Before discussing implications of these results for super-Earths, however, it is important to compare our model with earlier models for better understanding of the effects of adiabatic compression on thermal convection and the limitation of our model. The overall features of the convective flow patterns obtained here are similar to that of the “sluggish lid” convection identified for the incompressible (Boussinesq) fluids, where the lithosphere develops along the surface boundary and moves by its own negative buoyancy (Kameyama & Ogawa 2000; Ogawa 2008). However, a closer look at Figures 2 and 3 reveals that the flows lack secondary cold plumes originating from the base of the lithosphere, which are commonly observed in incompressible cases; adiabatic compression may modify the detail of the lithosphere on the “sluggish lid” regime from that observed in the incompressible cases (Moresi & Solomatov 1995; Kameyama & Ogawa 2000). Further systematic numerical simulations are necessary to elucidate the entire picture of the effect of adiabatic compression on thermal convection of variable viscosity fluid. The convective flow patterns shown in Figures 2 and 3 are quite different from the ones observed in earlier more sophisticated models of mantle convection with variable viscosity in super-Earths of various size (Tackley et al. 2013). Their cold plumes sink only to the middle of the mantle rather than to the core–mantle boundary, as we observe. Besides, there is no stratosphere-like layer like the one observed in Figure 2 in their models. However, it is not clear how the complicated rheology employed in Tackley et al. (2013) affects the convective flow pattern, and further systematic research, where the factors included in the rheology of Tackley et al. (2013) are taken into account one-by-one, is necessary to predict the dynamics of the mantle in super-Earths with certainty. Keeping this caveat in
Nu 100
10
1
107
10 8
109
1010 Ra
Figure 4. Nusselt numbers calculated for the cases shown in Figures 2 and 3 are plotted against the Rayleigh number (crosses). The dashed line shows the relation of N u ∝ Ra0.31 (Christensen 1984) while the star shows the result of Kameyama & Ogawa (2000), both obtained under the Boussinesq approximation.
approximation, Kameyama & Ogawa (2000) found Nu = 4.5 at r = 104 and Ra = 6 × 106 . When the scaling law of Nu ∝ Ra0.31 (Christensen 1984) is employed, the Nu at higher Ra is expected to be as shown by the dashed line in the figure. The Nusselt number Nu we calculated for the convection shown in Figures 2 and 3 is, however, lower than the expected Boussinesq value by a factor of ∼3.7 at Ra = 1010 . The efficiency of heat transport by thermal convection is significantly reduced by the strong effect of adiabatic compression. The lower Nusselt number implies that the lithosphere calculated in Figures 2 and 3 is much thicker than the lithosphere that develops under the Boussinesq approximation. 4. SUMMARY AND DISCUSSION The numerical models presented here suggest that adiabatic compression strongly affects mantle convection with strongly temperature-dependent viscosity at high Rayleigh number expected in super-Earths. Hot plumes induced by basal heating lose their thermal buoyancy by adiabatic cooling and stop ascending 4
The Astrophysical Journal Letters, 780:L8 (5pp), 2014 January 1
Miyagoshi et al.
mind, we suggest the possible implications of the effect of adiabatic compression on the mantle convection we found for the dynamics of the mantle in super-Earths, below. A straightforward implication of Figures 2 and 3 for sufficiently large super-Earths is that the strong effect of adiabatic decompression is likely to suppress hot spot volcanism caused by hot plumes ascending from deep mantle. Though hot plumes do develop in deep mantle (Figure 2(a)), they are unlikely to ascend to the surface. Besides, these plumes are so weak that they are almost invisible in the snapshots of temperature-distribution (Figures 2(b) and 3(b)). The sluggish movement of the lithosphere observed in Figures 2 and 3 may be important in understanding the surface environment of super-Earths. The moving lithosphere, if it occurs, may drive carbonate–silicate cycle and the circulation of water in super-Earths, and hence may stabilize the surface environment as moving plates do in the Earth (e.g., Kasting & Catling 2003). The average velocity of the lithosphere is about 1 cm yr−1 when scaled to a super-Earth with a mantle depth of 6000 km and is comparable to the Earth’s plate velocity. For this lithospheric speed to be attained, the viscosity contrast between the surface and the asthenosphere must be kept moderate owing to, say, the effects of water or elevated surface temperature by the blanket effect of the atmosphere. In future work, we will carry out further numerical experiments to find the range of viscosity contrast for the sluggish lid regime in super-Earths. The low efficiency of heat transport by thermal convection observed in Figure 4 also contains important implications. This figure suggests that the lithosphere does not become so thin as suggested from parameterized convection models in superEarths (e.g., Valencia et al. 2007a). The thick lithosphere would make it difficult for plates to bend and subduct at trenches. The inferred thick lithosphere is hence likely to significantly affect the dynamics of the lithosphere on super-Earths. The low efficiency of convective heat transport in the mantle may make it difficult to cool the core and to maintain magnetic field in the planet. It is important to carry out further numerical experiments with more realistic mantle rheology in order to more quantitatively estimate the thermal history of super-Earths under the influence of strong adiabatic compression. The thick lithosphere observed in Figures 2 and 3 may also indirectly affect the evolution of super-Earths. The thermal diffusion time of the lithosphere l 2 /κ (l is the thickness of the lithosphere) is on the order of 10 billion years for the convection of Figure 2, when the depth of the mantle is 6000 km as expected for a super-Earth of ten times the Earth’s mass. This time gives a measure of the time for mantle convection to go through with the initial transient stage affected by the process of planetary
formation. Indeed, we found that the initial transient stage is longer than 10 billion years, which is comparable to the age of the universe, in the calculation of Figure 2. This suggests that mantle convection on the initial transition stage, rather than the statistically steady stage, is more relevant to the evolutionary history of sufficiently large super-Earths. We thank the anonymous referee for helpful comments. Numerical simulations were performed on the Earth Simulator at Japan Agency for Marine-Earth Science and Technology. This work was supported by JSPS KAKENHI grant number 25287110. REFERENCES Bercovici, D., Schubert, G., & Glatzmaier, G. A. 1992, JFM, 239, 683 Borucki, W. J., Koch, D. G., Basri, G., et al. 2011, ApJ, 736, 19 Chopelas, A., & Boehler, R. 1992, GeoRL, 19, 1983 Christensen, U. R. 1984, PEPI, 35, 264 Gaidos, E., Conrad, C. P., Manga, M., & Hernlund, J. 2010, ApJ, 718, 596 Howard, A. W., Marcy, G. W., Johnson, J. A., et al. 2010, Sci, 330, 653 Jarvis, G. T., & McKenzie, D. P. 1980, JFM, 96, 515 Kameyama, M. 2005, J. Earth Simulator, 4, 2 Kameyama, M., Kageyama, A., & Sato, T. 2005, JCoPh, 206, 162 Kameyama, M., & Ogawa, M. 2000, E&PSL, 180, 355 Kasting, J. F., & Catling, D. 2003, ARA&A, 41, 429 King, S. D., Lee, C., van Keken, P. E., et al. 2010, GeoJI, 180, 73 Kite, E. S., Manga, M., & Gaidos, E. 2009, ApJ, 700, 1732 Korenaga, J. 2010, ApJL, 725, L43 Liu, X., & Zhong, S. 2013, GeoJI, 194, 125 Mayor, M., & Queloz, D. 1995, Natur, 378, 355 Moresi, L.-N., & Solomatov, V. S. 1995, PhFl, 7, 2154 Ogawa, M. 2008, FlDyR, 40, 379 Papuc, A. M., & Davies, G. F. 2008, Icar, 195, 447 Seager, S., Kuchner, M., Hier-Majumder, C. A., & Militzer, B. 2007, ApJ, 669, 1279 Sotin, C., Grasset, O., & Mocquet, A. 2007, Icar, 191, 337 Stamenkovi´c, V., Breuer, D., & Spohn, T. 2011, Icar, 216, 572 Stamenkovi´c, V., Noack, L., Breuer, D., & Spohn, T. 2012, ApJ, 748, 41 Tachinami, C., Ogawa, M., & Kameyama, M. 2013, Icar, submitted Tachinami, C., Senshu, H., & Ida, S. 2011, ApJ, 726, 70 Tackley, P. J., Ammann, M., Brodholt, J. P., Dobson, D. P., & Valencia, D. 2013, Icar, 225, 50 Tsuchiya, T., & Tsuchiya, J. 2011, PNAS, 108, 1252 Valencia, D., Ikoma, M., Guillot, T., & Nettelmann, N. 2010, A&A, 516, A20 Valencia, D., & O’Connell, R. J. 2009, E&PSL, 286, 492 Valencia, D., O’Connell, R. J., & Sasselov, D. D. 2006, Icar, 181, 545 Valencia, D., O’Connell, R. J., & Sasselov, D. D. 2007a, ApJL, 670, L45 Valencia, D., Sasselov, D. D., & O’Connell, R. J. 2007b, ApJ, 656, 545 van den Berg, A. P., Rainey, E. S. G., & Yuen, D.A. 2005, PEPI, 149, 259 van Heck, H. J., & Tackley, P. J. 2011, E&PSL, 310, 252 Wagner, F. W., Sohl, F., Hussmann, H., Grott, M., & Rauer, H. 2011, Icar, 214, 366 Yuen, D. A., Hansen, U., Zhao, W., Vincent, A. P., & Malevsky, A. V. 1993, JGR, 98, 5355
5
