A&A 445, 747-758 (2006)
H. Klahr1,2 - W. Kley1
1 - Universität Tübingen, Institut für Astronomie und Astrophysik, Abt. Computational Physics, Auf der Morgenstelle 10, 72076 Tübingen, Germany
2 - Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany
Received 13 April 2005 / Accepted 16 August 2005
We study the evolution of an embedded protoplanet in a circumstellar disk using the 3D-Radiation Hydro code TRAMP, and treat the thermodynamics of the gas properly in three dimensions. The primary interest of this work lies in the demonstration and testing of the numerical method. We show how far numerical parameters can influence the simulations of gap opening. We study a standard reference model under various numerical approximations. Then we compare the commonly used locally isothermal approximation to the radiation hydro simulation using an equation for the internal energy. Models with different treatments of the mass accretion process are compared. Often mass accumulates in the Roche lobe of the planet creating a hydrostatic atmosphere around the planet. The gravitational torques induced by the spiral pattern of the disk onto the planet are not strongly affected in the average magnitude, but the short time scale fluctuations are stronger in the radiation hydro models.
An interesting result of this work lies in the analysis of the temperature structure around the planet. The most striking effect of treating the thermodynamics properly is the formation of a hot pressure-supported bubble around the planet with a pressure scale height of rather than a thin Keplerian circumplanetary accretion disk.
Key words: accretion, accretion disks - hydrodynamics - solar system: formation - radiative transfer - planets and satellites: formation
These hydrodynamical simulations have been extended recently to full 3D-MHD calculations though neglecting vertical gravity (Winters et al. 2003; Papaloizou & Nelson 2003; Nelson & Papaloizou 2003). Here, the turbulent state of the disk is modeled directly, and the migration rate is found to fluctuate strongly around a mean value.
All the above mentioned simulations treat the disk as isothermal locally, although polytropic disks have been considered in Kley (1999). Recently, D'Angelo et al. (2003a) added thermodynamical effects in a two-dimensional model. They assume that the dissipation produced in the disk (by viscosity) is radiated locally, while the equation of state includes dissociation and ionization effects. Using this method in addition to nested grid refinement around the planet it also has been possible for the first time to obtain an estimate of the two-dimensional temperature structure in the disk's midplane in the vicinity of the planet. However, these simulations still suffer from their two-dimensionality.
In the present paper we extend the radiation hydrodynamical planet-disk calculations to three dimensions, by applying a flux-limited diffusion approximation for the radiation. This allows us to study directly the dynamical influence of the planetary accretion luminosity and to determine the three-dimensional temperature structure in the vicinity the planet.
Here, we focus on numerical issues, and study the effects of the radiation transport, of resolution and mass accretion. By introducing radiation and opacities the model is no longer scale-free, as it is in the isothermal case. Hence, there is a difference, whether the planet is located at 5 AU or 0.5 AU, say. In addition to influencing the gravitational torques acting on the planet, the disk mass now has an effect on the temperature structure as well, because it provides the optical depth.
We find that for the test parameters (1 Jupiter mass at 5 AU) and the low numerical resolution that we apply, there is no significant acceleration or breaking of the mean radial drift in comparison to the isothermal approximation. However, the fluctuations in the radial drift appear an order of magnitude larger in the radiative models with realistic accretion than in the isothermal models.
We find the circumplanetary "accretion disk'' to be rather thick, i.e. the pressure scale height is in our models, which is even thicker than the results from 2D simulations (D'Angelo et al. 2003a). The deviation from a Keplerian profile in this "disk'' is about , which also indicates the strong support from the thermal pressure. Under such conditions the thin disk approximation breaks down and one has to deal with a circumplanetary cloud filling a sizable fraction of the Roche lobe rather than a disk.
The detailed modeling of the observational appearance of a planet interacting with a disk and a detailed parameter study are left to follow up studies.
In the next section we describe our model setup, the boundary and initial conditions. In Sect. 3 we present an analysis of our models, and we conclude in Sect. 4.
Due to the assumed local vertical isothermal disk as initial setup,
the vertical density structure is given to first order
by a Gaussian:
In the models we do not use any explicit physical disk viscosity and model the ideal hydrodynamical equations. However, we use an artificial bulk viscosity in the form of von Neumann & Richtmyer (1950), see Stone & Norman (1992). Applying an explicit viscosity, e.g. -viscosity would have introduced one further free parameter, which we seek to avoid. Although MHD turbulence is a good candidate to produce an -like viscosity (Shakura & Sunyaev 1973), it is not obvious that would be a constant in time and space (see e.g. Wünsch et al. 2005). It is also unclear in how far the magneto-rotational instability would operate in the gap region and in the circumplanetary bubble.
A second reason for not using an explicit viscosity is that it would have led to additional heating on top of the adiabatic contraction of the circumplanetary bubble, thus making it more difficult to identify the importance of the heating by the accretion onto the planet.
To calculate the potential due to a planet of mass
in the midplane of the accretion disk,
of the Hill sphere around the planet, which is a function of
the reduced mass
of the system
In our models we test the influence the parameter has on the mass of gas that accumulates in the Roche lobe of the planet.
Table 1: Parameters chosen in the different simulations and their results. These parameters are the dimensioning of the grid ( ), the kind of spacing (lc = logarithmically centered around the planet), the thermodynamics (iso = locally isothermal, rad = ideal gas plus radiation transport, rad+ = ideal gas, radiation transport and conservation of total energy for gas accreted onto planet), the gravitational smoothing radius ( = fraction of the Hill radius over that smoothing is applied), the accretion radius ( = fraction of the Roche lobe (respectively Hill radius) that mass is taken out of), the type of the initial model (G: starting with an existing gap), and the total number of orbits calculated ( ). The results we give are the accretion rate onto the Planet () in units of Jupiter masses per orbit, the mass contained in the Roche lobe ( ), the maximum density and temperature in the Roche lobe ( ) and the migration time-scale that we determined ( ).
The radiation flux
is calculated in flux-limited diffusion approximation
The code has been extensively tested in the context of accretion disks (Klahr et al. 1999). The new part that is implemented to the code is the inclusion and modification of the gravitational potential of the planet as described above.
In Table 1 we give an overview of all the models we compare in the following. The second column (grid) states the numerical resolution ( ) of the grid. There are low resolution (10020200) models with a logarithmic radial spacing and equidistant in the azimuthal and vertical direction. We also calculated some high resolution models, where the grid is in all dimensions logarithmically concentrated (indicated by lc) towards the location of the planet, which results in a resolution in the Roche lobe up to four times higher than in the low resolution case. In these models there are about 2400 grid points within the Roche lobe. We tested the influence of the logarithmically distorted high resolution grid on a disk without a perturbing planet. Despite the grid refinement around the planet and quite distorted grid cells towards the edges of the grid, there are no unstable growing modes.
Some models are calculated in the standard locally isothermal approximation (EOS = iso), some treat the internal energy equation plus flux-limited radiation transport (EOS = rad). We also state the parameters for the smoothing of the gravitational potential and for the mass accretion onto the planet . Some models start with no gap initially, others use an axisymmetric initial model with a gap (ini = G).
After briefly mentioning a test model for the non-equidistant high resolution grid, we will go through our models in the following order: first isothermal models, then radiative models, followed by non-accreting models, and finally accreting and energy conserving models.
This model uses only one grid-cell in the azimuthal direction but it includes the full potential of the planet, and mass is taken out of the grid from the center of the Roche lobe. It is used to run a model quickly into an initial state with a gap, which forms here not due to the gravitational torques by the planet but rather by the perturbation in gravity as it would result from a massive ring. However, as we shall see the obtained gap size and shape are in reasonable agreement with the realistic 3D case.
In Fig. 1 we introduce the kind of plot that we will subsequently use for some of the following models. It shows quantities measured in the midplane at the azimuthal position of the planet. Thus, it is one-dimensional data. We plot 10 snapshots, which are taken at equal time intervals, e.g. in the case of Fig. 1 every 26 orbits from t=60 to t=320 orbits. First we plot the density . One clearly sees how the gap forms and how stable it is. Most of the disk is still at the initial values, which are given by the dashed line. The next plot shows the azimuthally averaged surface density . This value is slightly higher than at 5 AU, because here we integrate over the total height of the computational domain and not only over . The temperature T in the third plot does not fluctuate as this run is isothermal.
|Figure 1: Model CI: snapshots taken every 26 orbits from t=60 to t=320 orbits, which is after the convergence of model CI. The radial slices of the density, temperature and azimuthal velocity are taken in the midplane at the location of the planet. Model CI defines the initial values for all the models that start off with an axisymmetric gap.|
|Open with DEXTER|
The deviation of the local azimuthal velocity ( ) from the Keplerian value is significant only in the vicinity of the "planet''. The shape of the rotational profile is determined by the gradient in the surface density and thus has no physical relevance for the real gap case where the gravitational torques are determining the shape of the gap. We also have to mention that this artificial gap structure is not unique, but there is a family of solutions, where either the surface density profile or the azimuthal velocity profile is a free parameter. Still this method generates a gap structure that is not too unrealistic to serve as an initial model for our simulations.
The stationary state of model CI is given by the balance of centrifugal forces, radial pressure gradients and the attraction of the combined stellar and ring potential, is used as the initial setup in all of the subsequent models that include a C in their names (see Table 1). Models without a C are started from the same initial condition as model CI is started, which is a disk without a gap.
|Figure 2: Model CDI: snapshots taken every 2 orbits from t=164 to t=184 orbits. The radial slices of the density are taken in the midplane at the location of the planet. The surface density is azimuthally averaged.|
|Open with DEXTER|
The second reason for the mass depletion is the accretion process onto the planet, which occurs from both sides, from the inner as well as from the outer part of the disk. But the mass reservoir in the outer part of the disk is much larger, thus the decrease in mass is not visible over the span of the simulation. In comparison to the simulations by Bryden et al. (1999), Kley (1999), and Lubow et al. (1999) the gap in the surface density is shallower. Those simulations used explicit tensor viscosity and started already with a wider gap than we used. So in our case it will take many hundreds of orbits before the gap will reach a comparable depth. In terms of computational time such a run would take several months of computation time on our fastest available workstations. The code in its present version is not parallelized, due to the implicit solver of the radiation part. In the future we will consider runs with a much longer evolution time.
|Figure 3: Model CDI: ( upper row) accreted mass and accretion rate onto the planet over time; ( lower row) evolution of the mass contained in the Roche lobe and mass change in the Roche lobe in units of Jupiter masses and orbital periods of the planet.|
|Open with DEXTER|
Plotting the mass accreted onto the planet and the mass accumulated in the Roche lobe (see Fig. 3) also shows that we are only slowly approaching a steady state of the solution. There is a mass of Jupiter masses in the Roche lobe and the accretion rate, i.e. the mass taken off the grid, is Jupiter masses per orbit. The time that mass resides in the Roche lobe in this model is only about a tenth of an orbital period. The non-zero accretion rate is interesting because we include no viscosity in our model. One might expect that planets opening a gap in an inviscid disk could not accrete any more mass. There is the possibility that the planet-disk torques together with the induced Reynolds stresses provide enough angular momentum transport to feed the planet via the circumplanetary accretion disk. Nevertheless, to clarify this issue of the accretion rate onto the planet the simulations will have to be continued for many more orbits at high resolution. It may still be possible that accretion will stop after a thousand orbits, which is the case in non-viscous 2D simulations of gap opening.
We also calculate the gravitational torques acting on the planet and translate them into the typical radial migration time (Kley et al. 2001). We measure the torques in three ways: excluding the mass contained in the inner half Hill radius of the Roche lobe (in the following the dashed curves), excluding all the Roche lobe (solid curves) and excluding 1.5 times the Hill radius (dotted curves). Bate et al. (2003) have used this method to give error bars to the measurements of the torques from the numerical simulations. This procedure is also very useful to study the relative importance of the torques generated from gas close to the planet. In Fig. 4 we see that the torques vary over time but the mean drift rate corresponds to the 105 yrs observed in other simulations and which fit the results from linear theory (Ward 1997; Tanaka et al. 2002; Goldreich & Tremaine 1980). It is remarkable how close the three measurements of the torques are at each time. This is a consequence of the fact that we take the mass out of a significant part of the Roche lobe producing a very clean gap. There is not enough matter left, to exert a torque on the planet. This changes in the models that do not treat mass accretion onto the planet ( ) or remove gas from a smaller part of the Roche lobe ( ).
|Figure 4: Model CDI: evolution of planet migration rate in terms of migration time. The solid curve is based on the torques excluding the Roche lobe and the solid horizontal line is the time average. In addition we plot the torques based on an exclusion of 1.5 and 0.5 Hill radii (dotted and dashed lines respectively). In this particular model the gap is very clean, such that all three torques mostly coincide.|
|Open with DEXTER|
|Figure 5: Model CD4I: snapshots taken each orbit from t=71 to t=81 orbits. The radial slices of the density are taken in the midplane at the location of the planet. The surface density is azimuthally averaged. The dotted line indicates the values for model CDI for comparison.|
|Open with DEXTER|
Nevertheless, the overall structure of the gap has not changed substantially. This model takes about 5 times more computational time than the low resolution case. This is mainly due to the fine resolution in the azimuthal direction, i.e. the Courant condition for the azimuthal velocity component. Table 1 shows that there is a little more mass accumulating in the Roche lobe (which was empty initially) than in the low resolution case. The torques are now strongly fluctuating while the average values are very similar to those of model CDI (see Table 1). The fluctuations are so strong that they even change their sign. Again the time averaged torques are insensitive to the radius around the planet, which has been excluded from the calculation of the torque. We suspect that the oscillations might be an effect of the high density concentration near the planet (see Fig. 5), because small density asymmetries induce large torque variation. Similar findings have been obtained for nested grid simulations by D'Angelo et al. (2002).
In future work we will return to mesh refinement around the planet to study these effects in more detail.
|Figure 6: Model CD4: snapshot at t=5 orbits after the refinement of the grid. Radial slices of density and temperature are taken in the midplane at the location of the planet, while surface density is the azimuthal average. (solid line = Model CD4; dotted line = Model CD; dashed line = Model CDS).|
|Open with DEXTER|
The isothermal approach used in the previous models is a strong simplification of the thermodynamics of the disk near the planet and may not be well suited for studying the gap opening problem. For the radiative model CD we keep all parameters identical to model CDI, but use the equation for internal energy including flux-limited radiation transport, see Eq. (6). Again, mass is removed from within the Roche lobe, which means that energy is taken out without releasing its accretion energy (this will be treated differently in models DR and DR4). The temperature and density structure of this model are displayed together with other models in Fig. 6 with dotted lines. In contrast to model CDI the disk has cooled dramatically because there is no internal heating in our inviscid model, while the region around the planet has become slightly hotter than in the isothermal gap case (Fig. 6). The rest of the gap is much cooler. There is also three times more mass piling up in the Roche lobe than in Model CDI (see also Table 1). This is a direct consequence of the radiative cooling. If the matter within the Roche lobe cools down it will contract leading to a higher density. The mass accretion rate onto the planet is roughly three to four times larger than for the corresponding isothermal model. Thus, there is more mass stored in the Roche lobe, which also influences the torque (see Fig. 7). While the overall average migration time is only 10% less than in model CDI, there are now variations between the torques measured by excluding more or less mass from the Roche lobe. Excluding either the entire Roche lobe or even 1.5 times the Hill radius gives very similar results, while excluding only 0.5 of the Hill radius doubles the amplitude of the fluctuations.
|Figure 7: Model CD: evolution of the planet migration rate in terms of migration time. The solid curve is based on the torques excluding the Roche lobe and the solid horizontal line is the time average. In addition we plot the torques based on an exclusion of 0.5 Hill radii (dotted lines) and 1.5 Hill radii (dashed lines). The horizontal lines refer to the time averaged values.|
|Open with DEXTER|
This is the first hint that radiative cooling changes the structure of matter in the Roche lobe and consequently the resulting torques on the planet.
The accretion rate is also close to the one in model CD. The average torques are still the same as in the isothermal case, model CDI. It also shows that the migration times are the longest if one excludes only half the Hill radius. This indicates the importance of the proper treatment of mass in the Roche lobe and its potential capability of altering the radial drift rate. This model is also used for comparison to model CD4 and CD in Fig. 6.
Here, we increase the resolution and decrease the gravitational smoothing again to 20% of the Hill radius (cf. model CD4I). The radius out of which mass may be taken is decreased in the same fashion. Otherwise this model (see Fig. 6) is identical to models CD (dotted line) and CDS (dashed line).
The most striking effect of the higher resolution is the higher temperature at the planet's location. Temperatures three times higher than in model CD are reached.
This can be explained by two effects which are: 1. a deeper potential well leading to stronger compression and release of accretion energy; 2. the much higher densities near the planet leading to higher optical depth and lower cooling rates.
Again the accretion rate is similar to model CD and the average torques are also not very different (see Table 1).
However, in the models with radiation transport one introduces a large error if one simply takes mass off the grid: internal energy is removed from the simulation along with the mass, but the accretion energy released during final accretion phase onto the planet is neglected. One can now either construct a model for this process and re-adjust the energy content of the gas after some mass has been taken away (see Models DR and DR4) or simply refrain from accreting mass and study the effects of the mass piling up. The first approach has one more free parameter than the second, which is the relative accretion rate (i.e. the fraction of the mass removed from the grid per orbit).
We will first study more closely the effect of taking mass arbitrarily from the grid. Therefore, we go back to the isothermal model CDI but do not remove any mass this time.
In comparison to model CDI the lack of accretion leads to a slightly higher concentration of mass around the planet, similar to the models with a deeper potential (see Table 1). The mass in the Roche lobe is now three times greater than in model CDI.
The major effect of not taking mass off the grid and letting it pile up in the Roche lobe is the introduction of strong fluctuations in the migration rate. The torques calculated including half of the Hill radius especially showed large fluctuations. Nevertheless, the time averaged migration rate is still dominated by the structures outside the Roche lobe and thus are of the same magnitude as in all the previous models (see Table 1).
|Figure 8: Model DR4: early stage of gap opening. Temperature in the plane of the protoplanetary disk at the azimuthal location of the planet. This snapshot is taken after 12 orbits, thus the gap is not completely empty yet. In this early phase the heating by the accretion onto the planet drives a strong convectional flow in the still optically thick disk material above the planet. Brightness is logarithmic temperature (lighter = warmer) ranging from 38 K to 1060 K and vectors give the logarithmic mass flux. We also plot iso-density lines with labels giving values in g cm-3. Density ranges from 1015 g cm-3 up to g cm-3. See the online edition for a color version of this plot.|
|Open with DEXTER|
In these models we start with a gap, do not remove any mass and perform the full radiation step. Models starting with a gap eventually converge to models without an initial gap. This convergence can be demonstrated for the low resolution models (see Table 1), while the high resolution models are too expensive in computational effort to show this in a reasonable computation time. The general trend can be seen in Table 1. Density and temperature increase with the resolution of the Roche lobe while the torques are not much different in the high and low resolution cases.
|Figure 9: Model DR: snapshot taken after t=121 orbits (solid line). Model DR4: after t=141 orbits (dotted line). Model DN: after t=121 orbits (dashed line). Model DN4: snapshot taken after t=55 orbits (dash dotted line). Radial slices of density and temperature are taken in the midplane at the location of the planet, while surface density is the azimuthal average.|
|Open with DEXTER|
|Figure 10: Model DR: ( upper row) accreted mass and accretion rate onto the planet over time; ( lower row) evolution of the mass contained in the Roche lobe and mass change in the Roche lobe in units of Jupiter masses and orbital periods of the planet.|
|Open with DEXTER|
|Figure 11: Model DR: evolution of planet migration rate in terms of migration time. The solid curve is based on the torques excluding the Roche lobe and the solid horizontal line is the time average. In addition we plot the torques based on an exclusion of 0.5 Hill radii (dotted lines) and 1.5 Hill radii (dashed lines), plus the time averaged values.|
|Open with DEXTER|
First we show the flow situation in model DR4 after 12 orbits past the initial set up. The gap has not cleared completely and the heating by accretion onto the planet drives strong vertical convection in the optically thick disk material above the planet (see Fig. 8). This "fountain'' flow is subsonic and will not escape from the disk because the velocity is too low in comparison to the escape velocity. Nevertheless, the effect might be very important for observations because even embedded planets still too small to open a gap will probably produce this effect, which shows as an extended bubble of hot gas to the planet thus raising the local pressure scale height in the disk.
Model DR approaches a quasi steady state after 100 orbits (see Fig. 9), which can also be inferred from the almost constant Roche lobe mass of and the accretion rate onto the planet of per orbit (see Fig. 10). Nevertheless, there are also short time scale fluctuations in the Roche lobe mass, which again show up in the fluctuations of the torques (see Fig. 11). Based on non-viscous 2D runs we know that it will take up to a thousand orbits before the gap is so deep that no accretion onto the planet occurs any longer. In 3D with radiation transport we presently cannot afford such long runs. The resulting maximal temperature in the Roche lobe is 800 K for model DR and 1500 K for model DR4. A finer grid resolution might even lead to higher temperatures.
|Figure 12: Model DR: density in the midplane of the protoplanetary disk after 121 orbits. Dark means low, light means high density. The density ranges from to . See the online edition for a color version of this plot.|
|Open with DEXTER|
|Figure 13: Model DR: temperature in the midplane of the protoplanetary disk after 121 orbits. Brightness corresponds to the logarithm of temperature. Black indicates 10 K and white 800 K. See the online edition for a color version of this plot.|
|Open with DEXTER|
|Figure 14: Model DR4: pressure scale height in the midplane of the circumplanetary cloud after 141 orbits. Radial location is given in units of Hill radii.|
|Open with DEXTER|
For model DR we also display the density structure that results in the midplane (see Fig. 12) after 122 orbits. The gap has clearly formed but still mass is left in the horseshoe orbit. This mass will decrease gradually over time. The density structure does not change significantly from any isothermal simulation, but the temperature of course does (see Fig. 13). The planet is the hottest spot in the simulation. We adopted a grey approximation for the radiation transport, thus continuum (wavelength dependent) radiation transport simulations (using our temperature and density structure) will show whether this hot circumplanetary gas-dust mixture will be observable. The extent of the warm gas around the planet as well as the two mass feeding streams of warm gas hitting the circumplanetary cloud are clearly visible now. We do not resolve any circumplanetary accretion disk yet, for two reasons: first our resolution is too low and second the pressure dominates the hot gas around the planet strongly enough so that no Keplerian disk forms, and one is left with a pressure-supported planetary envelope filling most of the Roche lobe.
In Fig. 14 we display the pressure scale height of the gas in the midplane around the planet. In the circumplanetary "disk'' is the inverse Mach number of the flow around the planet. Due to the high temperatures the relative vertical thickness is about 0.5, which is larger than those values typically used in thin Keplerian accretion disks with H/R < 0.1. In the case of this structure is better described by a rotating envelope which additionally shows strong deviations from a Keplerian rotation law of about .
For model DR4 we display some details of the flow onto the planet (see Fig. 15). We plot the temperature in the plane of the protoplanetary disk at the azimuthal location of the planet. The temperature is hottest around the planet and the gas in the gap is also heated by the accretion luminosity of the planet. The flow indicates that accretion mainly occurs via the poles of the planet and that there is convection in the envelope around the planet. We can see no inflow along the equatorial plane. The mass falling onto the central region around the planet originates only in part from the gap material above the planet. The largest contribution in mass is via the convective flow in the upper part of the circumplanetary cloud.
This simulation also suggests that at least at this stage of gap opening and planet disk interaction there is no formation of moons from the circumplanetary material possible. The reason is simple: with and the resulting strongly sub-Keplerian rotation, all solid material will quickly rain down onto the planet. The maximum radial drift velocity is given by the difference between the actual azimuthal velocity and and the Keplerian value (Weidenschilling 1977). Thus before solid material has orbited the planet once, it will fall onto the planet. The fastest drifting solids are given by unity in the parameter , with the Keplerian frequency and the friction time of the solid body. For the physical conditions in our model DR4 (e.g. density and temperature) the fastest drifting objects would be 1 m sized boulders, which drift into Jupiter in about a month if released at 0.1 Roche radii. Thus the drift problem for the formation of moons from circumplanetary material is even more crucial than in the case of planet formation. We conclude that at least for the time when Jupiter is accreting its mass there can be no formation of a satellite system.
Mean torques and migration rates are not strongly affected by the treatment of the thermodynamics in the case of Jupiter mass planets. This might change for lower mass planets, which are more embedded in the disk. The fluctuations of the torques on the other hand are much stronger, in particular in higher resolution cases. Similar effects have been observed in high resolution nested grid simulations (D'Angelo et al. 2003b,2002) and also MHD simulations of planet-disk interactions (Winters et al. 2003; Papaloizou & Nelson 2003; Nelson & Papaloizou 2003). In some cases the torques even change their sign for a short period.
We find that planets are most likely to form a circumplanetary pressure-supported envelope rather than an accretion disk around them, with strong convective vertical flows. The relative pressure scale height in the circumplanetary material is at least 0.5, in which case the approximations for a thin Keplerian accretion disk are no longer valid. What results is a cloud which rotates at only 50 percent of the Keplerian value.
This slow rotation of the circumplanetary gas makes it impossible to form the moons of Jupiter from circumplanetary material as the solid content will rain down the planet in less than an orbit around Jupiter.
|Figure 15: Model DR4: Temperature in the plane of the protoplanetary disk at the azimuthal location of the planet after 141 orbits. Brightness is logarithmic temperature (lighter = warmer) between 30 K and 1500 K, contours are equi-density lines (in g cm-3) and vectors give the logarithmic mass flux. See the online edition for a color version of this plot.|
|Open with DEXTER|
How does this thick and slowly rotating envelope behave with respect to magneto-centrifugal wind as discussed for a circumplanetary disk by Fendt (2003)? The high thermal pressure could be beneficial for launching the magneto-centrifugal wind, but on the other hand we observe a strong mass infall onto the poles of the planet and the surface of the circumplanetary cloud, thus a wind might be suppressed. Detailed investigations on that issue are needed.
The deviation from symmetry with respect to the midplane is very small in our simulations. Some convection cells show eventual overshooting across the midplane but in general, treating only one hemisphere would have led to the same results.
We observe strong vertical flows in the early gap opening phase. Considering the entire mass accretion phase of the young planet starting from a few earth masses (Pollack et al. 1996), one finds that these vertical fountains may last for hundred-thousands of years before the planet opens its gap, making the effect clearly relevant for observations. As a result there would be locally a stronger flaring of the disk, and more radiation could interact with the small dust grains in the surface layer above the planet. Future instruments like SIM should be able to observe this asymmetry in the scattered light from disks in which planets are forming (Geoffrey Bryden, personal communication). On the other hand, the effect might be weaker if we allow for a gradual buildup of the planet. This has to be investigated in future simulations.
We do not observe the formation of vortices at the edges of the gap, which can often be found in the equivalent 2D simulations (Kley 1999; Koller et al. 2003). Further investigations will have to show whether this is related to the 3D nature of our simulations or a consequence of the comparable low resolution of the grid.
Finally, we speculate that the extended warm gas and dust which fills most of the Roche lobe should be detectable in the far infrared. Investigations such as by Wolf & D'Angelo (2005) indicate this already for 2D flat simulations that include radiative cooling. Similar studies using 3D density and temperature distributions will follow in the near future.
The authors wish to thank Geoff Bryden for fruitful discussions. This work was sponsored by the European Community's Research Training Networks Programme The Origin of Planetary Systems under contract HPRN-CT-2002-00308.