Free Access
Issue
A&A
Volume 560, December 2013
Article Number A40
Number of page(s) 10
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/201322503
Published online 02 December 2013

© ESO, 2013

1. Introduction

The early stages of planet formation take place in protoplanetary disks around young stars and the growing and evolving planets shape the structure of the disks. One important observational goal is the detection of signatures that give direct hints of the presence of planets in such disks. One class of systems that have been linked to planets are the so-called transitional disks that show a deficiency in the IR excess in their spectral energy distribution (SED) and/or show an extended inner hole in the disk in direct images. In recent years many of these transitional disks have been discovered, for example by the Spitzer Space Telescope (Muzerolle et al. 2010; Kim et al. 2013). The occurrence rate of transitional disks rises monotonically with age from about 1% at an age of one million years to about 20% at 8 million years (Muzerolle et al. 2010) and so they are potentially very interesting indicators of the late phase in the planet formation process. Disks with gaps or completely cleared inner cavities usually fshow mass accretion onto the star, but with accretion rates typically an order of magnitude smaller than that of a continuous disk (Kim et al. 2013). As there is a clear decrease in from continuous disks over disks with gaps to disks with cavities, Kim et al. (2013) suggest that the presence of sub-stellar companions (planets) may be the most likely explanation of the properties of transitional disks.

Growing planets can open gaps within the disk after accumulating a sufficient amount of mass, and so they constitute a natural possible explanation (Lin & Papaloizou 1993). Indeed, single planet simulations of embedded planets could show the formation of large inner cavities with a reduced mass flow through the inner hole (Quillen et al. 2004; Varnière et al. 2006). In these simulations, a massive Jupiter-like planet opened a gap in the disk and the inner cavity was cleared out rapidly through mass flow through the inner boundary onto the star. While it was possible to produce a deep cavity, it has been pointed out that this may be due to the special outflow boundary condition used at the inner edge (Crida et al. 2007). Additional clearing of the inner disk in the case of one embedded planet region may also be mediated by disk photo-evaporation (Alexander & Armitage 2009), which has been investigated in simulations by Rosotti et al. (2013). Here, the timescale of the evaporation process has to be chosen such that it is compatible with the observed long life times of accretion in T Tau disks (Zhu et al. 2011). Because it is known that multiple planets can clear out much deeper and wider gaps than single planets (Kley 2000), it has been suggested that the large gaps and inner cavities may be caused by the presence of multiple (up to four) planets (Dodson-Robinson & Salyk 2011; Zhu et al. 2011). In this scenario, however, mass transfer has to occur across the planets and the system has to be dynamically stable over long timescales.

The recent observations presented by Mayama et al. (2012) that show a kind of spiral arm connecting the outer disk through the cleared area with the star strengthened these hints of planets within transitional disks. In either case, multiple planets or a single planet, the important question is, how strong is the mass flow from the outer disk across the planet into the inner cavity onto the star?

Even though the option of multiple planets has been suggested as a possible cause for the observed properties of transitional disks, we will nevertheless focus in this paper of the effect of single embedded planets in the disk. As pointed out above, even in the presence of multiple planets it is necessary for mass to flow from the outer disk into the inner gap region to provide for the observed accretion signatures onto the star. Here, we focus on the dynamics of single massive planets that are surrounded by an outer accretion disk and no inner counterpart. For smaller mass planets (<Mjup) the mass accretion rate across a gap has been analysed by Lubow & D’Angelo (2006). This single planet could be considered the outermost in a multi-planet system. Our simulations will give important information about the mass flow from the outer disk into this planetary system.

We are particularly interested in the amount of mass that can flow across massive planets into the inner cavity. For this purpose we have performed a detailed analysis of this process considering planets of different masses. Additionally, we have improved on the thermodynamics of the disk and consider isothermal as well as radiative disks with and without irradiation from the central star.

In Sect. 2 we describe our physical and numerical modelling of the process. The standard model is described in Sect. 3. In Sect. 4 we present the results of the isothermal runs for different planet masses together with some numerical tests. This is followed in Sect. 5 by the radiative results that are compared to even more extended models that include irradiation in Sect. 6.

2. Model set-up

2.1. Physics and equations

We assumed an infinitesimally thin disk around the star and therefore solved the vertically integrated versions of the hydrodynamical equations. For the coordinate system we chose cylindrical coordinates (r,ϕ,z), centred on the star where the disk lies in the equatorial z = 0 plane.

The vertically integrated versions of the continuity equation and the equations of motions in the r − ϕ plane can be found in Müller & Kley (2012, Sect. 2.1). Here, we included an energy equation to allow for a more realistic thermodynamic treatment of the disk. The vertically integrated energy equation reads

(1)where e is the internal energy density, Q+ the heating source term, and Q the cooling source term.

The cooling term Q is given by the vertically integrated divergence of the radiative flux,

(2)

where σR is the Stefan-Boltzmann constant, ρ the density, κ the Rosseland mean opacity and T the temperature. It can be written as (3)Here Qrad describes the radiative losses from the disk surfaces which is calculated using suitably averaged opacities (Müller & Kley 2012, Sect. 2.1), and the second term corresponds to the radiative diffusion in plane of the disk; H = csK is the pressure scale height given by the midplane sound speed, cs, and the Keplerian rotational angular velocity, ΩK, of the disk. The coefficient of the radiative diffusion in the rϕ plane is given by .

The heating term Q+ consists of the viscous dissipation Qvisc and the irradiation from the star Qirrad. The stellar irradiation to one surface of the disk can be approximated by (4)(Günther et al. 2004), where Tstar is the effective photospheric temperature of the star and Rstar the radius of the star. The factor β accounts for the non-perpendicular impact of the irradiation from the star onto the disk.

2.2. Numerical considerations

We used the adiabatic version of the FARGO code (Masset 2000; Baruteau 2008). For the radiative cooling Qrad from Eq. (3) we used our implementation from Müller & Kley (2012). The radiative diffusion source term is solved separately after the rest of the source-terms as a flux-limited diffusion equation for the temperature (Kley 1989; Kley & Crida 2008). We solve it implicitly using a successive over-relaxation (SOR) method. For the stellar irradiation we implemented the approximation by Günther et al. (2004).

The planet in the disk has a Plummer type gravitational potential to account for the vertical extent of the disk and to avoid numerical problems of a point-mass potential. We used a smoothing value of ϵ = 0.6H as this describes the vertically averaged forces very well (Müller et al. 2012).

To avoid numerical problems the density cannot fall below a minimum value of Σfloor = 10-7·Σ0, where Σ0 is the reference density at r = rjup = 5.2 au, and the temperature is always at least Tfloor = 3   K, which is about the temperature of the cosmic background radiation.

Because the disk has an abrupt end at the inner edge of the computational domain, the stellar irradiation produces here artificial high temperatures because of the large amount of energy deposited. To account for this we introduced an estimated horizontal optical depth by rescaling the vertical optical depth τ (see Müller & Kley 2012, Eq. (12)) by the radial extent of the corresponding grid cell (5)where h = H/r is the aspect ratio of the disk. If we set Qirrad = 0. This simple procedure keeps the temperatures within the gap region regular and allows for standard heating in the outer parts of the disk.

2.3. Boundary conditions

To maintain a given disk structure in the outer parts of the disk we implemented the damping mechanism by de Val-Borro et al. (2006), where specified quantities like radial velocity and angular velocity are damped towards their initial values. This damping is described by (6)where ξ ∈  { Σ,vr,vϕ,e }, τ is the damping timescale, and R(r) is a linear ramp-function rising from 0 to 1 from the damping star radius to the outer radius of the computational domain. To be able to measure the mass loss rate in our models we only damped the velocity components within 0.9   rmax to rmax using a time factor of τ = 3 × 10-2·2π   ΩK(rmax)-1.

In addition to the damping we used a reflecting boundary condition at the outer edge so that no mass could escape through the outer boundary. For the inner edge we used a zero-gradient outflow boundary condition and measured the amount of mass lost through this boundary to calculate the mass accretion rate onto the star. The details of the implementation of the boundary conditions can be found in Müller & Kley (2012, Sect. 2.2).

3. The standard model

Our model consists of a planet with a varying mass of 1 to 16 Jupiter masses (Mjup) orbiting a Sun-like star with a semi-major axis of 1   rjup = 5.2   au. The planet is fixed in its circular orbit and therefore cannot change its orbital parameters during the simulation. Beyond rjup the star is surrounded by a gaseous disk that extends up to rmax = 5   rjup, the outer edge of the computational domain. The inner radius depends on the model and ranges from 0.7 to 0.9   rjup, where we chose for the main part of our models rmin = 0.7. The gas disk is set up initially as if there were no planet with a simple surface density power-law profile of r-0.5 and a temperature power-law profile of r-1. To prevent strong shock waves in the initial phase of the simulations the planet’s mass is ramped up slowly over 50 orbits. All the models were run for a total simulation time of 5000 planetary orbits. For the viscosity we used a constant kinematic viscosity ν = 1015  cm2   s-1 that corresponds to an α-value of about 4 × 10-3 at the reference radius. Table 1 summarizes the parameters of the standard model.

Table 1

Parameters of the standard disk model.

4. Isothermal simulations

For the first set of simulations we used a locally isothermal approach keeping the initial temperature stratification throughout the whole simulation fixed. We assumed a constant aspect ratio H/r = 0.05, that corresponds to a temperature profile of r-1 with 137   K at the inner edge and 24   K at the outer edge of the disk. Similar simulations have been performed by Kley & Dirksen (2006) for a constant surface density profile; they damped the surface density at the outer edge of the disk.

4.1. Inner radius dependency

In this first set of simulations, a fixed planet mass of Mp = 4  Mjup was assumed, while the location of the inner radius was varied to check its influence. The L1 point lies at about 0.896 au and so all radii were chosen in such a way that the Roche-Lobe of the planet is fully included in the computational domain. Kley & Dirksen (2006) showed, that for planets with Mp ≥ 3  Mjup the disks become eccentric after a very long time (>1000 orbits). Figures 1 and 2 show the azimuthally averaged surface density and disk eccentricity profiles after 5000 orbits for Mp = 4  Mjup and five different inner radii. For the eccentricity we calculated first a mass weighted average of each grid cell and then averaged over the azimuthal direction. We note that the models are in different phases of their eccentricity oscillation as shown in Fig. 5. The gap opened by the planet at r = 1   rjup is clearly visible and the gas in the vicinity has eccentric orbits with an eccentricity of up to 0.22. The eccentricities for all models are very similar. The model with the smallest inner radius shows a slightly smaller eccentricity just outside of the planet inside the gap. This does not influence the mass flow across the planet, however. The eccentricity can also be noticed in the surface density plot show in Fig. 3. In this r-ϕ plot, circular motions are represented by horizontal lines. However, the gas at the outer edge of the gap has wavelike perturbations which visualize the eccentricity of its orbit. This is similar to Fig. 1 of Kley & Dirksen (2006) where the outer disk also becomes eccentric.

thumbnail Fig. 1

Azimuthally averaged radial profile of surface density of the isothermal set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER

thumbnail Fig. 2

Azimuthally averaged radial profile of eccentricity of the isothermal set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER

thumbnail Fig. 3

Surface density of the isothermal set-up after 5000 orbits with a planet of Mp = 4   Mjup and an inner radius of rmin = 0.7   jup. The planet is located at the angle of π in this r − φ representation.

Open with DEXTER

These quasi-stationary states of the disks need a rather long time to establish. Figure 4 shows the disk mass evolution and Fig. 5 the disk eccentricity evolution of the disks over the total simulation time of 5000 orbits. The planet’s mass is ramped up during the first 50 orbits which is indicated by the dotted vertical line. For about the first 1000 orbits (depending on the inner radius), the disk evolution is very smooth with a nearly constant mass content and small eccentricity. After this time the disk becomes eccentric and thereby loses more mass through the inner boundary. The simulations with smaller computational domain (larger rmin) need some more time to develop the disk eccentricity, but after that display a similar evolution.

Figure 6 shows the normalized mass accretion rate through the inner boundary of the computational domain using a moving average over 50 orbits. The normalized mass flow − (t)/Mdisk(t) is shown where the actual mass of the disk is used, as the surface density Σ and therefore the disk mass Mdisk cancel out of the equations in locally isothermal simulations. During the first ~1000 orbits the mass accretion rate is almost constant at about . However, as soon as the disk becomes eccentric (see Fig. 5) the mass accretion rate starts to increase as more mass is crossing the planetary orbit and flows through the inner boundary at rmin.

The initial mass loss through the inner boundary is smaller than the standard mass accretion rate of stationary, viscous accretion disks as given by (7)After the disk becomes eccentric, it seems to converge against the value, which is indicated by the horizontal dashed line.

We note that all the figures with normalized mass accretion rates are displayed in units of , whereas accretion rates are usually given in M   a-1. The 0.7   rjup models have an initial disk mass of 0.12   M, and because the disks lose less than 20% of their mass during the simulation time of 5000 orbits (see Fig. 4), one can get a rough estimate of the mass accretion rate in M   a-1 by dividing the values in the figures by a factor of (8)The shapes of the time-dependent curves (Figs. 46) appear to be shifted in time, because the disks with the smaller inner computational radii need more time to loose the mass in the outer disk. This leads to a surface density enhancement at the outer edge of the gap and subsequently to stronger eccentricity excitations.

thumbnail Fig. 4

Total disk mass over time for the isothermal set-up with a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 5

Global disk eccentricity over time for the isothermal set-up with a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 6

Normalized mass accretion rate over time for the isothermal set-up with a planet mass of Mp = 4   Mjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER

4.2. Planet mass dependency

The final state of the simulations are not strongly dependent on the minimum radius of the computational domain. So we chose rmin = 0.7   rjup for the following simulations to compare different planet masses. As shown by Kley & Dirksen (2006), the outer disk turns eccentric for planetary masses larger than about 3   Mjup. So we varied the planet mass from 1   Mjup to 16   Mjup by doubling it each step. Figures 7 and 8 show again the radial profiles of the surface density and the eccentricity at the end of simulation after 5000 orbits. We note that the models are in different phases of their global eccentricity oscillation as shown in Fig. 10. For example, in Fig. 8 the 16   Mjup model shows a maximum eccentricity of about 0.35, whereas the maximum eccentricity varies between 0.32 and 0.4 during one period of the oscillation. In agreement with the previous results, the simulations with the small mass planets, 1 and 2 Mjup, do not show any significant disk eccentricity. Beyond a planet mass of about 3−4 Mjup the disk becomes eccentric. As expected, the more massive planets form a much wider gap in the disk, and the eccentricity of the disk is much higher for the 4   Mjup and more massive planets.

thumbnail Fig. 7

Azimuthally averaged radial profile of surface density of the isothermal set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER

thumbnail Fig. 8

Azimuthally averaged radial profile of eccentricity of the isothermal set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER

In Figs. 9 and 10 the disk mass and disk eccentricity over time are shown. The average disk eccentricity increases clearly with the planetary mass. It is easy to see that the higher the planet mass is, the faster the disk eccentricity starts to grow and that the mass loss starts earlier. Whereas the previously examined example of the 4   Mjup planet needs about 1000 orbits until the disk eccentricity grows, the lower mass planet models do not develop a noteworthy eccentricity at all. The models with more massive planets start to develop much large eccentricities almost instantly. The normalized mass accretion rate is displayed in Fig. 11. There is no clear trend that with increasing planet mass the mass accretion rate grows or shrinks because there are two competing effects. The growing disk eccentricity favours the mass transport across the planet into the inner cavity of the disk whereas the larger gap due to the larger planet mass hinders mass transport through this gap. All the mass accretion rates are of the order of the typical disk accretion rates from Eq. (7).

thumbnail Fig. 9

Total disk mass over time for the isothermal set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 10

Global disk eccentricity over time for the isothermal set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 11

Normalized mass accretion rate over time for the isothermal set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER

5. Radiative simulations

In Sect. 4 we performed locally isothermal simulations and maintained the initial temperature stratification throughout the whole simulation. In this section we extend the set-up and also solve for the energy equation (Eq. (1)), but neglect the radiative diffusion (2H   ∇(KT) = 0), because in standard accretion disk models this term is much smaller than the vertical cooling. Additionally, we neglect stellar irradiation (Qirrad = 0). We redid the same set of calculations (except for the rmin = 0.4   rjup case) as in the isothermal case for a direct comparison of the results.

5.1. Inner radius dependency

Figures 12 and 13 show the surface density and eccentricity profiles at the end of the simulation. Whereas the shape of surface density profile does not look too different from the isothermal models shown in Fig. 1, the absolute values within the gap are about a magnitude smaller. All models display a much lower total eccentricity. The planets also open a gap that can be seen in Fig. 14, but compared with to the isothermal case (Fig. 3) the outer edge of the gap is much less distorted because of the much smaller eccentricity. Only within the gap does a significant eccentricity develop and it is a factor of three to four smaller than in the isothermal simulations (see Fig. 3).

thumbnail Fig. 12

Azimuthally averaged radial profile of surface density of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER

thumbnail Fig. 13

Azimuthally averaged radial profile of eccentricity of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER

thumbnail Fig. 14

Surface density of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup and an inner radius of rmin = 0.7   jup. The planet is located at the angle of π.

Open with DEXTER

The global disk eccentricity is more or less constant over time and it reaches its final, very small value shortly after the planet has been ramped up. The value of 0.013 at the end of the simulation is about the same as the value for the first stage in the isothermal evolution (Fig. 5). Because of the small constant value of the disk eccentricity there is a small mass loss through the inner edge of the simulation region, which can be seen in the disk mass over time simulation shown in Fig. 15. As in the isothermal models, the mass loss at the start of simulation is dependent on the exact position of the inner boundary, but for the rest of the simulation the results are only shifted.

thumbnail Fig. 15

Total disk mass over time for the radiative set-up with a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

The normalized mass accretion rate shown in Fig. 16 is therefore almost constant except for some jitter. For the 0.7   rjup model we obtain an average value of after 1000 orbits. Assuming a constant disk mass (see Eq. (8)), this corresponds to a mass accretion rate of (6.680 ± 0.008) × 10-9M   a-1, which fits very well to observational data (Kim et al. 2013). Interestingly, it is now much smaller than the value given by Eq. (7) in contrast to the isothermal simulations, as the disk is much less eccentric and therefore pushes less mass through the gap. In addition, the disk can radiate its energy away and is therefore much cooler which also hinders mass transport through the gap. In the next section we present a comparison of the temperatures in the isothermal and radiative disks.

thumbnail Fig. 16

Normalized mass accretion rate over time for the radiative set-up with a planet mass of Mp = 4   Mjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER

Figure 17 shows a zoomed version of Fig. 14 with a different scaling of the surface density. Within the gap almost no gas is left (<0.001   g   cm-2), but the spiral arm induced by the planet is clearly visible.

thumbnail Fig. 17

Surface density of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup and an inner radius of rmin = 0.7   jup. The planet is located at the angle of π.

Open with DEXTER

5.2. Planet mass dependency

For the comparison of different planet masses we again chose rmin = 0.7   rjup for the inner boundary. Figures 18 and 19 display the surface density and eccentricity profiles at the end of the simulation. As in the isothermal, case the more massive planets open a deeper gap in the disk. In the isothermal simulations, the disk eccentricity profiles were much larger for masses Mp ≥ 4   Mjup, whereas in the radiative models this is only the case for masses Mp ≥ 8   Mjup. In addition, all eccentricities are smaller than those seen in the isothermal case. We note again, that the models are in different phases of their eccentricity oscillation as shown in Fig. 22. Because the temperature profile can also change during the simulation, Fig. 20 shows the radial temperature profile at the end of the simulation. The temperature in the gap is rather low and therefore we can see the deeper gaps for the more massive planets even in the temperature profile. The maximum temperature is always at the outer edge of the gap and is dependent on the planet’s mass. The edge temperature is the highest for the 4   Mjup and lower for higher and smaller planet masses.

thumbnail Fig. 18

Azimuthally averaged radial profile of surface density of the radiative set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER

thumbnail Fig. 19

Azimuthally averaged radial profile of eccentricity of the radiative set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER

thumbnail Fig. 20

Azimuthally averaged radial profile of temperature of the radiative set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER

The global disk mass evolution is shown in Fig. 21. The mass loss through the inner edge now depends on two factors. On the one hand, the mass flow through the gap will be increased by the disk’s eccentricity as in the isothermal case. The global disk eccentricity evolution in time as shown in Fig. 22 is similar to that of the isothermal case. In the radiative case the planet mass required to reach a certain disk eccentricity is higher than in the isothermal case. On the other hand thermal pressure changes the depth of the gap which may lead to an enhanced mass flow across the gap. This is also mass dependent as already seen in Fig. 20. Figure 23 now shows the normalized mass accretion rate for the radiative models. As both effects interact we can see no clear trend for the mass accretion rate. For the two lowest mass planets (1 and 2 Mjup), the mass accretion rates are highest because they have the smallest gap size (see Fig. 18). The 4   Mjup case has a wider and deeper gap with no disk eccentricity and shows a much smaller accretion rate. For the 8   Mjup case the disk becomes eccentric and the mass flow increases again, while for the 16   Mjup case the gap widens substantially such that the mass flow rate becomes very small.

thumbnail Fig. 21

Total disk mass over time for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 22

Global disk eccentricity over time for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 23

Normalized mass accretion rate over time for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER

5.3. Disk mass dependency

In the isothermal models the normalized mass accretion rate was independent of the disk mass Mdisk because as the surface density Σ cancels out of the equations. Therefore, the mass accretion rate scales with the same factor as the surface density, for example if we decrease the surface density by a factor of 10, the mass accretion rate decreases by a factor of 10.

In the radiative case, however, this is no longer true. Figure 24 shows the normalized mass accretion rates and the mass accretion rates for different disk masses. The normalized mass accretion rates (upper panel of Fig. 24) seem to scale in the wrong direction at first sight. For example, the Mdisk = 0.12   M model has a smaller normalized mass accretion rate than the Mdisk = 0.012   M model, but if we calculate the mass accretion rates (lower panel of Fig. 24) with Eq. (8) the trend flips as expected and the Mdisk = 0.12   M model has a larger mass accretion rate than the Mdisk = 0.012   M model, but it does not scale with the same factor as the disk mass. For the rest of the simulations we used the Mdisk = 0.12   M model.

thumbnail Fig. 24

Normalized mass accretion rate over time (upper panel) and mass accretion rate over time (lower panel) for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

thumbnail Fig. 25

Total disk mass over time for different set-ups with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

6. Comparison of different models

In this section we extend our models again for radiative diffusion and stellar irradiation and compare the following four different models: the isothermal set-ups from Sect. 4 with a constant temperature stratification throughout the simulation, the radiative set-ups from Sect. 5 where we solve the energy equation with viscous dissipation and radiative cooling, the radiative with diffusion set-ups where we also solve the radiative diffusion in the rϕ plane, and the radiative with irradiation set-ups where we include radiative diffusion and stellar irradiation. As radiative diffusion and stellar irradiation are very costly in terms of computation time we ran those models only with an inner radius of rmin = 0.7   rjup and planet masses of Mp = 4   Mjup and Mp = 8   Mjup.

Figures 25 and 26 show the disk mass and disk eccentricity evolution over the whole simulation for the Mp = 4   Mjup models. The isothermal model shows a huge difference in both plots. The disk loses much more mass and becomes significantly more eccentric. Between the radiative model and the radiative with diffusion model there is barely a visible difference, so that the neglect of the radiative diffusion in accretion disk simulations is a valid assumption. The addition of irradiation increases the mass loss again because the disk is much hotter, but it does not play an important role for the global disk eccentricity.

thumbnail Fig. 26

Global disk eccentricity over time for different set-ups with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER

The normalized mass accretion rate of all models in shown in Fig. 27. As already seen in Fig. 25, the mass loss of the isothermal model is much larger than all the other models. The radiative with diffusion model has a slightly larger accretion rate than the radiative model, but both are smaller than the standard mass accretion rate of stationary accretion disks (Eq. (7)). The accretion rate of the radiative with irradiation model needs much more time to converge and is significantly larger than those of the other two radiative models.

thumbnail Fig. 27

Normalized mass accretion rate over time for different set-ups with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER

Therefore, the mass accretion in the irradiative model must be driven by thermal pressure through the gap, as the eccentricity is the same. This can be seen in Fig. 28, which shows the temperature profiles of all models after 5000 orbits. The isothermal model still has the unchanged initial condition, whereas the radiative models show a temperature drop in the gap. The radiative with diffusion set-up shows a slightly higher temperature in the gap as heat is diffused into the gap from the outer edge. The temperature of radiative with irradiation is much higher overall, especially in the outer regions, but the gradient on the outer edge of the gap is also much steeper and therefore has pushed more material through the gap.

thumbnail Fig. 28

Azimuthally averaged radial profile of temperature of the isothermal, radiative, radiative with diffusion and radiative with irradiation set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup.

Open with DEXTER

thumbnail Fig. 29

Azimuthally averaged radial profile of the surface density of the isothermal, radiative, radiative with diffusion and radiative with irradiation set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup.

Open with DEXTER

In Fig. 29, we display the surface density at the final state after 5000 orbits for all our models. In the outer part of the disk, all models show comparable surface densities. In the inner gap region, the densities are typically 3 to 4 orders of magnitude smaller than in the outer disk. As expected by the mass accretion rates, the isothermal and irradiated models show an enhanced density in the gap region. In contrast, the radiative models have a surface density about an order of magnitude smaller inside the gap. For the most realistic irradiated models, the density contrast is about 4 orders in magnitude. For models with the 8 Mjup planet the principal surface density structure is very similar for the radiative models.

7. Summary and conclusions

We have performed hydrodynamic simulations of transitional disks where the inner edge is formed by a very massive planet and studied the mass flow rate of the gas from the outer disk into the inner cavity. For comparison, we have studied locally isothermal disks as well as viscously heated radiative disks including irradiation from the central star.

We find that the mass flow past the planet is largest for the isothermal models. These have a much higher temperature in the inner regions of the disk and the corresponding pressure gradient is partly responsible for the high mass flow. The second contribution comes from the eccentricity of the disk for planet masses beyond roughly 3   Mjup. The disk has its largest eccentricity near its inner edge, and reaches about 0.2 for the 4 and 8 Mjup planets and up to 0.35 for the 16   Mjup planet in the isothermal disk. For larger radii the disk eccentricity slowly declines. For eccentric disks the planet passes periodically through the disk which leads to the enhanced mass accretion into the inner gap. The asymmetries in recent observations (e.g. Mayama et al. 2012) can be explained by these eccentricities or by a vortex within the disk (Ataiee et al. 2013).

In the radiative cases, there are two reasons why the mass mass accretion rate is generally much smaller. First, the temperature is much smaller in the gap region which leads to a smaller pressure effect, and second the radiative models have a smaller disk eccentricity than the isothermal models. The maximum disk eccentricity is around 0.2 for the most massive 16 Mjup planet. The reduction in eccentricity comes about because the emission of radiation from the disk surfaces leads to an energy leakage of the eccentric modes of the disk. The model including stellar irradiation shows again a larger mass flow which is still smaller than the purely isothermal disk.

The surface density contrast between the outer and inner disk is interesting. In the eccentric disk case we find a contrast of about four orders of magnitude which is even larger in the non-eccentric disk cases. Despite this large contrast in surface density, the mass accretion rate past the planet into the gap changes by much less. For the isothermal models it is always comparable to the typical disk accretion rate (see Fig. 11) in contrast to the observations that show a reduction of about an order of magnitude compared to stationary disks (Kim et al. 2013). However, for the more realistic radiative cases, the mass accretion rate can be reduced by over an order of magnitude, if the disk is not eccentric, i.e. if the planet is smaller than about 5 Mjup. This mass flow rate is in good agreement with the observations (Kim et al. 2013). For larger mass planets and eccentric disks the mass flow rate is about 40% of the standard equilibrium disk value (see Fig. 23).

Our model contains several restrictions that could be improved on in subsequent works. The models are only 2D because 3D runs over this extended simulation time of several thousand orbits are not feasible at the present time. We estimate that because of the large mass of the planets and the deep gaps that are opened, the differences between 2D and 3D may be not too large. For Jupiter-mass planets the gaps are nearly indistinguishable (Kley et al. 2001) in 2D and 3D. Furthermore, the planets cannot accrete material from their surroundings. Because of the strong decline of the density in the gap, this effect may be more important for low mass planets as well. It may play a role, however, in the eccentric disk case where it could lead to higher planet masses than otherwise reachable (Kley & Dirksen 2006).

In our models we kept the planet fixed in a circular orbit, while it is known that the planets will tend to migrate because of the torques of the outer disk. For circular planets in non-eccentric disks, the migration will be directed inwards. However, the eccentric disk for more massive planets will make the orbit the slightly eccentric and the migration rate slows considerably (D’Angelo et al. 2006). For e ≥ 0.2, the migration may

be even outwards, but these high planetary eccentricities are only reached for very massive planets beyond 10   Mjup (Papaloizou et al. 2001). However, all these results are only obtained using locally isothermal disks. The question of planetary migration and possible eccentricity growth will have to be addressed in future simulations. In our case, we might consider the simulated single planet being the outermost planet in a multiple system. Our results should then give a good approximation of the expected mass flow of the gas into the gap.

We also did not include any dust particles within the models. Rice et al. (2006) showed that the outer edge of the gap acts as a filter for larger particles (≳10 ?m). However, the simulations by Rice et al. are only locally isothermal. The question whether this filter mechanism also works within more realistic, radiative simulations will also have to be addressed in future simulations.

Acknowledgments

Tobias Müller received financial support from the Carl-Zeiss-Stiftung. Wilhelm Kley acknowledges the support of the German Research Foundation (DFG) through grant KL 650/8-2 within the Collaborative Research Group FOR 759: The Formation of Planets: The Critical First Growth Phase. We thank Andres Carmona for very fruitful discussions. We thank the referee for providing constructive comments which helped to improve this paper. Most of the simulations were performed on the bwGRiD cluster in Tübingen, which is funded by the Ministry for Education and Research of Germany and the Ministry for Science, Research and Arts of the state Baden-Württemberg, and the cluster of the Collaborative Research Group FOR 759.

References

All Tables

Table 1

Parameters of the standard disk model.

All Figures

thumbnail Fig. 1

Azimuthally averaged radial profile of surface density of the isothermal set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER
In the text
thumbnail Fig. 2

Azimuthally averaged radial profile of eccentricity of the isothermal set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER
In the text
thumbnail Fig. 3

Surface density of the isothermal set-up after 5000 orbits with a planet of Mp = 4   Mjup and an inner radius of rmin = 0.7   jup. The planet is located at the angle of π in this r − φ representation.

Open with DEXTER
In the text
thumbnail Fig. 4

Total disk mass over time for the isothermal set-up with a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 5

Global disk eccentricity over time for the isothermal set-up with a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 6

Normalized mass accretion rate over time for the isothermal set-up with a planet mass of Mp = 4   Mjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER
In the text
thumbnail Fig. 7

Azimuthally averaged radial profile of surface density of the isothermal set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER
In the text
thumbnail Fig. 8

Azimuthally averaged radial profile of eccentricity of the isothermal set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER
In the text
thumbnail Fig. 9

Total disk mass over time for the isothermal set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 10

Global disk eccentricity over time for the isothermal set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 11

Normalized mass accretion rate over time for the isothermal set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER
In the text
thumbnail Fig. 12

Azimuthally averaged radial profile of surface density of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER
In the text
thumbnail Fig. 13

Azimuthally averaged radial profile of eccentricity of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup.

Open with DEXTER
In the text
thumbnail Fig. 14

Surface density of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup and an inner radius of rmin = 0.7   jup. The planet is located at the angle of π.

Open with DEXTER
In the text
thumbnail Fig. 15

Total disk mass over time for the radiative set-up with a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 16

Normalized mass accretion rate over time for the radiative set-up with a planet mass of Mp = 4   Mjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER
In the text
thumbnail Fig. 17

Surface density of the radiative set-up after 5000 orbits with a planet of Mp = 4   Mjup and an inner radius of rmin = 0.7   jup. The planet is located at the angle of π.

Open with DEXTER
In the text
thumbnail Fig. 18

Azimuthally averaged radial profile of surface density of the radiative set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER
In the text
thumbnail Fig. 19

Azimuthally averaged radial profile of eccentricity of the radiative set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER
In the text
thumbnail Fig. 20

Azimuthally averaged radial profile of temperature of the radiative set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup.

Open with DEXTER
In the text
thumbnail Fig. 21

Total disk mass over time for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 22

Global disk eccentricity over time for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 23

Normalized mass accretion rate over time for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER
In the text
thumbnail Fig. 24

Normalized mass accretion rate over time (upper panel) and mass accretion rate over time (lower panel) for the radiative set-up with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The values are plotted as a moving average over 50 orbits to remove jitter. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 25

Total disk mass over time for different set-ups with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 26

Global disk eccentricity over time for different set-ups with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass.

Open with DEXTER
In the text
thumbnail Fig. 27

Normalized mass accretion rate over time for different set-ups with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup. The vertical dotted line at t = 50   Pjup marks the time when the planet has reached its full mass and the horizontal dashed line is the gas accretion rate given by Eq. (7).

Open with DEXTER
In the text
thumbnail Fig. 28

Azimuthally averaged radial profile of temperature of the isothermal, radiative, radiative with diffusion and radiative with irradiation set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup.

Open with DEXTER
In the text
thumbnail Fig. 29

Azimuthally averaged radial profile of the surface density of the isothermal, radiative, radiative with diffusion and radiative with irradiation set-up after 5000 orbits with an inner radius of the computational domain of rmin = 0.7   rjup and a planet mass of Mp = 4   Mjup.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.