Free Access
Issue
A&A
Volume 579, July 2015
Article Number A105
Number of page(s) 19
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201526117
Published online 08 July 2015

© ESO, 2015

1. Introduction and motivation

It is commonly assumed that planets form in circumstellar disks around young pre-main-sequence (PMS) stars (e.g., Mordasini et al. 2010). During this formation process planets perturb the disk because of their gravitational potential (Goldreich & Tremaine 1980). The results of this interaction are characteristic structures such as spiral arms (Ward 1997; Tanaka et al. 2002) and gaps (Papaloizou & Lin 1984). Therefore, observations of these large-scale structures provide indirect constraints on the existence of protoplanets and on the relevant physical processes involved in the protoplanet-disk interaction.

Wolf & D’Angelo (2005) and Ruge et al. (2013) showed the possibility of tracing gaps in continuum re-emission light for various disk configurations and ALMA configurations. In addition, Ruge et al. (2014) investigated the feasibility of tracing gaps in continuum re-emission maps and scattered light images simultaneously. In contrast to this, only a few studies exist investigating the feasibility of tracing planet-induced structures through molecular line observations. Semenov et al. (2008) showed the feasibility of drawing conclusions about disk properties such as large-scale temperature gradients and chemical stratification from simulated HCO+ (43) observations. Dutrey et al. (2008) used spatially unresolved but spectroscopically resolved line observations to trace inner cavities due to the radial dependence of the (Keplerian) velocity field. Cleeves et al. (2011) demonstrated the feasibility of detecting a huge inner cavity (R ~ 45 AU) using molecular lines (CO, H2CO2). Using ALMA, Casassus et al. (2013) were able to image multiple molecular gas species within the dust free inner hole of HD 142527 and to utilize the varying optical depths of the different species to detect a flow of gas crossing the gap opened by a planet.

However, the main component of protoplanetary disks, cold molecular hydrogen (H2), is very difficult to observe directly in these disks because of the missing electrical dipole moment. Thus, other molecules with lower excitation temperatures must be observed to derive the dynamics, temperature, density, and chemical structure in these disks. The most frequently observed molecule is CO, mainly because of its high abundance and low excitation energy. Observations of various isotopologues of this molecule provided the first constraints on the physical properties such as the (vertical) gas temperature structure of young disks (e.g., Piétu et al. 2007). Through comprehensive surveys of protoplanetary disks around solar-type stars other molecules such as HCO+, H2CO, CS, CN, and HCN were detected (Kastner et al. 1997; Öberg et al. 2010; Thi et al. 2004). Recently, Chapillon et al. (2012a) provided evidence of HC3N in the disks around GO Tau and MWC 480 using the IRAM 30 m telescope. However, because of insufficient spatial resolution and the lack of sensitivity these observations were restricted to the simplest and most abundant molecules. The Atacama Large (Sub)-Millimeter Interferometer (ALMA) will potentially allow the direct study of the molecular composition and distribution in protoplanetary disks in the submillimeter wavelength regime both at global scale as well as in the potential planet-forming region. The first ALMA observations of HD 163296 show that CO radiation originates from a thin layer with an opening angle of about 15° with respect to the disk midplane (Rosenfeld et al. 2013; de Gregorio-Monsalvo et al. 2013). More recently with the unprecedented sensitivity of ALMA, cyclic C3H2 (which has been proposed as a useful probe of radiation penetration in disks) has been detected in the disk surrounding HD 163296 (Qi et al. 2013).

The motivation to search for gaps in protoplanetary disks using molecular line observations is given not only by the aim to indirectly detect planets, but even more to derive complementary constraints on the disk physics from the planet-disk interaction. For example, the resulting gap structure, given by its shape and depth in the gas and dust distribution, depends not only on the mass of the planet and grain size, but also on the interaction of the planet with the gas and dust phase (Paardekooper & Mellema 2004). As a result, larger dust grains (150 μm) are expected to decouple from the gas, pile up at the edges of the gap and significantly alter the gas-to-dust mass ratio (Paardekooper & Mellema 2006). Moreover, we emphasize that a gap in the density distribution of a cirumstellar disk is not necessarily an irrevocable sign of an embedded planet. Flock et al. (2015) showed that MRI turbulence is sufficient to generate a gap in the gas phase without the presence of a planet. Recently, Bruderer et al. (2014) were able to detect CO emission inside the dust cavity of IRS 48. They concluded that a massive inner companion is the main driver for the clearing of the cavity and excluded photoevaporation and grain-growth only. Therefore, complementary observations of the gas and the dust phase are crucial in order to improve and verify our understanding of the predominant physics; e.g., the interaction between the gas and dust component in the planet forming regions of protoplanetary disks. As a first step into this potentially new field of studying disk physics through planet-disk interaction, we investigate the feasibility of tracing planet-induced gaps through spatially and spectroscopically resolved molecular line observations. For this purpose we developed the N-LTE line and dust continuum radiative transfer code Mol3D, which is based on the radiative transfer code MC3D (Wolf et al. 1999; Wolf 2003).

Our studies are based on spatial density and velocity distributions resulting from 3D hydrodynamic simulations of planet-disk interaction using the PLUTO code (Mignone et al. 2007). The temperature structure and corresponding molecular line maps are calculated with Mol3D. These velocity-channel maps are used as an input for the CASA (ver. 4.2) simulator (Petry & CASA Development Team 2012). Finally, we investigate the feasibility of tracing planet-induced gaps in the simulated ALMA maps. Based on these results, our goal is to prepare future ALMA observations of similar protoplanetary disks. In contrast to earlier studies, which were often focused on exemplary case studies (e.g., Cleeves et al. 2011), we explore a large parameter space covering various disk properties such as their total mass, size, and abundant molecules of typical PMS T Tauri and Herbig Ae stars.

This study is structured as follows. In Sect. 2 we briefly present the radiative transfer code Mol3D, discussing its underlying physical assumptions and test cases. In Sect. 3 we describe the study of tracing gaps using molecular lines with ALMA. Finally, in Sects. 4 and 5 we discuss our results and provide conclusions about the feasibility of tracing planet-induced gaps though molecular line observations with ALMA.

2. Molecular line radiative transfer code

In this section, we briefly describe our line and dust continuum radiative transfer code Mol3D.

Mol3D features

The N-LTE 3D parallelized line radiative transfer code Mol3D features full 3D spherical, cylindrical, and Cartesian grids with several types of boundary spacing (e.g., logarithmic, linear). Mol3D typically uses three steps to produce velocity-channel maps, spectra, dust continuum maps, and spectral energy distributions (SED). In the first step, the dust temperature is calculated based on the assumption of local thermal equilibrium. To achieve a high efficiency, we use a combination of the continuous absorption method proposed by Lucy (1999) and immediate temperature correction by Bjorkman & Wood (2001).

In the second step the level populations are calculated. As discussed in Appendix B.2.1, this is the main problem of line radiative transfer. In the case of circumstellar disks, the problem can be significantly simplified by using adequate approximation methods. Mol3D allows the local thermodynamic equilibrium (LTE), the free-escape probability (FEP), or the large velocity gradient (LVG) method to be applied. It has been shown that LTE is a reliable assumption for lower transitions of abundant molecules, but it is recommended that LVG be used, especially for high transitions of complex molecules with low abundances (Pavlyuchenkov et al. 2007).

The last step is to generate velocity spectra or velocity-channel maps using a new efficient ray-tracing algorithm. Here, we have to consider that the source function is not constant inside a single cell, which is due to the underlying velocity field.

As our code solves the radiative transfer equation including both gas and dust, we have to take the optical properties of both into account. For the intensity integration, we implemented an embedded Runge-Kutta-Fehlberg solver of order 4(5) with automatic step control. The main advantage is that with providing a relative (10-8) and an absolute (10-20) error limit for the intensity integration, the algorithm calculates the corresponding step-width automatically. Conservation of energy is ensured by a recursive pixel refinement algorithm based on the underlying model grid.

Mol3D is tested and benchmarked against other line and continuum radiative transfer codes. The results are presented in Appendix B along with a detailed derivation of the underlying concepts and physical assumptions. The code will soon be publicly available1 under an open-source license.

3. Model setup

In this section we investigate the feasibility of observing planet-induced structures using molecular lines with ALMA. For this purpose we perform 3D hydrodynamic (HD) simulations, follow-up line-radiative-transport (LRT) calculations, and simulate ALMA observations using the CASA simulator.

3.1. Three-dimensional hydrodynamic simulations

thumbnail Fig. 1

Example of HD simulations with embedded Jupiter-mass planet (planetmass/starmassratio = 10-3). Color-coded is the H2 number density in log scale. Left: midplane density. The planet has cleared out its orbit, producing a prominent gap. Owing to the perturbations, complex structures like spiral arms are also visible. Right: cut perpendicular to the disk midplane. The gap is clearly visible in the gas distribution. We note that the upper layers of the disk have nearly the same density as in the undisturbed case.

Open with DEXTER

In preparation for this study we started with an analytic density distribution including a radially symmetric sharp gap completely depleted of dust and gas. We find that this approach is insufficient because it is based on an unrealistically high density contrast between the gap and its surroundings, resulting in a gap detection rate of almost 100% for the planet-disk configurations considered in our study (see Sect. 3.2). Instead, justified assumptions about the ratio between the density inside and outside the gap, the gap width, and the depth of the gap are needed.

To reach this goal, we make use of full 3D hydrodynamic simulations obtained with the well-known finite-volume fluid dynamics code PLUTO (Mignone et al. 2007). The code uses the HLLC2 solver for the hydrodynamics, a second-order Runge-Kutta time integration and second-order variation diminishing (TVD) scheme. The setup of the simulations was originally presented in Uribe et al. (2011) and the models used in this work are shown in Ruge et al. (2013). From their entire set of planet-disk configurations we consider two scenarios. The first model includes an embedded Jupiter-mass planet (fixed mass ratio Mplanet/Mstar = 0.001) on a fixed circular orbit. The second disk with the same initial parameter setting does not contain a planet, and thus serves as an undisturbed smooth reference disk for comparison.

We use a spherical grid where θ is the angle from the vertical axis and φ is the azimuthal angle. To save computation time, the vertical extent of our PLUTO simulations is chosen to cover only θ [π/20.3, π/2 +0.3], covering four pressure scale heights above the disk midplane.

As initial conditions, it is assumed that the disk gas has a sub-Keplerian rotation profile, and that the density profile ρ0ra and the speed of sound cs = c0(rsinθ)b follow a power-law with exponents a = 3/2 and b = 0.5, respectively. The initial density distribution is given by (1)The disk is described by a locally isothermal equation of state . The ratio of the pressure scale height h to the radial coordinate of the disk is taken to be a constant such that h = H/ (rsinθ) = 0.07. A softened gravitational potential is used for the planet, (2)where ϵ is the softening parameter which is defined as a fraction of the Hill radius ϵ = l·rp(Mp/ 3)1/3 with l = 0.3.

The original disk configuration has an inner and outer edge of 2 AU and 9 AU, respectively. The planet is located at rp = 5 AU. For the hydrodynamics simulations, buffer zones are included to avoid boundary effects, so the complete radial domain is from 1 − 10 AU. For the follow-up radiative transfer, we need to convert the resulting density distribution. We average every two cells in each direction to save memory. Thus, we are able to reduce the required number of grid cells to . For the subsequent radiative transfer (RT) simulations with Mol3D we add two cells (one above and one below the disk) to fill the entire θ [0, π] domain. The grid parameters for both codes are summarized in Table 1.

Table 1

Comparison of grid parameter space of both codes.

3.2. Disk model

The hydrodynamic simulations discussed in Sect. 3.1 rely on arbitrary units for length and mass. This allows us to introduce linear scaling parameters to scale the disk mass and size to our needs for protoplanetary and/or transitional disks.

Disk size

For the inner and outer disk radius we chose the linear scaling parameter k ∈ [7,..,25]. Hence, the inner radius of the disk is in the range of Rin = k·2 AU (i.e., Rin [14,...,50 AU]) and the outer radius in the range of Rout = k·9 AU (i.e., [63,...,225 AU]. As a result, our disk models mimic different evolutionary states from early protoplanetary disks to more evolved transitional disks or circumbinary disks with huge inner cavities. The disk density distribution resulting from the simulation of a Jupiter-mass planet located at 80 AU in a disk ranging from Rin = 32 AU to Rout = 144 AU (corresponding to a scaling parameter of k = 16) is presented in Fig. 1. Several specific features in the disk density distribution that result from planet-disk interaction can be seen. The most prominent one is the distinct gap between 68 and 92 AU. In the cut through the density distribution perpendicular to the disk midplane, the gap appears as a tunnel in the protoplanetary disk, as the planet has hardly any influence on the disk surface (see Fig. 1, right). Furthermore, the midplane density distribution (Fig. 1, left) shows pronounced spiral arms; these features are barely visible in the surface of the disk in this case as well.

Mass range

We consider a wide range of disk masses (gas+dust) from 2.67 × 10-1 M to 2.67 × 10-5 M which approximately covers the range of protoplanetary and transitional disk masses derived from submm/mm observations (e.g., Andrews et al. 2010, 2011).

Dust properties and stellar parameters

For the dust we use a mixture of 62.5% astronomical silicate and 37.5% graphite (Weingartner & Draine 2001). Furthermore, we assume the dust grains to be spherical following a grain size distribution given in Eq. (3) with an exponent of d = 3.5 (Dohnanyi 1969): (3)We consider dust grains with a minimum radius of amin = 5 nm and a maximum radius of amax = 250 nm as found in the interstellar medium (Mathis et al. 1977). The radiation source is the central star. For this study we consider a T Tauri and a Herbig Ae star with parameters given in Table 2).

Table 2

Considered pre-main-sequence stars and their properties.

Distribution of the different molecular species

For the gas-phase, we assume a gas-to-dust-mass ratio of 100:1 and perfect thermal coupling between gas and dust.

For the spatial distribution of different molecular species relative to each other only weak constraints are available from observational and numerical studies. For example, observations have shown carbon monoxide emission from the cold outer parts of protoplanetary disks where CO should be frozen out (<~20 K; Dartois et al. 2003; Goto et al. 2012). A possible explanation has been proposed by Aikawa & Nomura (2006) and Semenov et al. (2006). They argue, that the CO abundance could be significantly higher in the cold (midplane) regions of protoplanetary disks owing to vertical and radial mixing. Motivated by this idea, Hersant et al. (2009) computed the chemical evolution of circumstellar disks. They included grain surface reactions with and without turbulent mixing and CO photodesorption. However, the resulting CO column densities are still not consistent with observations. As indicated by Cleeves et al. (2011), the disk midplane in the gap region can be irradiated by the stellar radiation and X-rays. Thus, CO would not be frozen out onto grains, which would result in better gap detection rate thanks to higher contrast.

Because the debate on relative molecular abundances in circumstellar disks is still open, we decided to use fixed abundances relative to molecular hydrogen. In particular, we follow Pavlyuchenkov et al. (2007) and assume their reported abundances for all molecules considered in this study (see Table 3). These values are obtained by using a gas-grain chemical model with surface reactions (Semenov et al. 2004, 2005), which is based on the UMIST 95 database of gas-phase reactions (Millar et al. 1997). This approach keeps the model simple and hence reduces the complexity of our analysis.

Table 3

Properties of the different molecules considered in this study.

Line radiative transfer setup

The level populations are calculated by using the LVG method (see Appendix B.2.1) which has proven its reliability for protoplanetary disks (Pavlyuchenkov et al. 2007). Using these results, we produce synthetic velocity-channel maps and spectra for every model configuration. We assume a pure Keplerian rotation profile. In this study we fix the turbulent velocity vturb = 100 m/s for all of our models, following the results of earlier molecule observations (e.g., Piétu et al. 2007; Chapillon et al. 2012b). We choose a fixed bandwidth of 50 m/s and calculate 141 velocity-channel maps for each of the models with 10° inclination and 71 for every face-on orientated model in our sample. An overview of the considered parameter space of our model is compiled in Table 4.

4. Results

4.1. Ideal line maps

In this section we present ideal velocity-channel maps, analyze the appearance of planet-induced disk perturbations, and investigate the feasibility of tracing these structures with ALMA (Sect. 4.2). In the following we use the term “depth of the gap” as the flux density contrast in the velocity-channel maps caused by the gap.

Optical depth effects

We begin our discussion based on ideal 12C16O (32) velocity-channel maps of our T Tauri disk model. To start with, we assume a total disk mass of 2.67 × 10-5 M, which is at the very low end of typical disk masses of these

Table 4

Parameter space of the disk model.

objects (e.g., HD 166191, Coka Tau/4; Kennedy et al. 2014; Forrest et al. 2004; D’Alessio et al. 2005). The inner rim of the chosen disk has a radius of 32 AU while the outer radius amounts to 144 AU. The planet is located at ~80 AU and the planet-induced gap has a width of ~32 AU. The disk is seen face-on. In this case, the line shows a symmetric gap throughout the different velocity channels, because the line width is given only by thermal and microturbulent velocities (see Fig. 2; see also Appendix B.2). We find that the contrast between the gap and its edges is not constant throughout the individual velocity channels, but varies with the optical thickness of the chosen line. In the optically thin case, the gap depth increases as it approaches the line center owing to the increasing line flux. In our case, the 12C16O (32) transition is optically tick at its center even for this lowest disk mass of our sample. Thus, the τ = 1 surface seen in the line wings belongs to parts of the disk closer to the disk midplane than the disk layers traced at the line center (i.e., at 0 m/s). Therefore, for this specific configuration we observe the gap best at velocities ±300 m/s, while the main contribution at the line center comes from the upper layers of the disk. As these upper layers are less disturbed by the planet-disk interaction, the gap is less pronounced in the center of an optically thin line, i.e., the nearly undisturbed upper disk layers mask the gap. In analogy, further signatures of the planet-disk interaction like the pronounced spiral arm in the outer disk can be identified best in the line wings (here, in the ±400 m/s channel maps).

thumbnail Fig. 2

Ideal 12C16O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet with an orbital distance of 80 AU. The disk mass amounts to 2.67 × 10-5 M, which is the lowest disk mass in our sample, and the outer radius is 144 AU. The flux density in Jy/mas2 is color coded for each individual channel. The disk is seen face-on (0°). The gap is clearly visible, particularly in the line wings. Even more complex features like the huge spiral arm (e.g., at 300 m/s) can be recognized. The total width of the line is only given by thermal and microturbulent velocities (Sect. B.2).

Open with DEXTER

If we increase the optical depth by increasing the disk mass by two orders of magnitude (i.e., 2.67 × 10-3 M), which represents a typical disk mass of protoplanetary disks, we find that the gap vanishes nearly completely in the velocity-channel maps (Fig. 3). As only the envelope above the gap can be observed, the disk appears to be undisturbed without any clear signs of planet-disk interaction. This is the main problem when observing 12C16O. For disk masses corresponding to typical protoplanetary disks masses we find that the 12C16O molecule is not an adequate tracer for gaps. In this case, other less abundant molecules/isotopologues that allow deeper disk layers to be observed are better suited. We illustrate this on the basis of the 12C18O (32) transition (Fig. 4). Because of the low abundance of this carbon monoxide isotopologue, this line is less optically thick and thus traces deeper disk layers. However, its excitation temperatures of the lower energy levels are very similar to those of the 12C16O molecule. In contrast to the optically thick 12C16O line of the same model, the gap is clearly detectable.

thumbnail Fig. 3

Ideal 12C16O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet (at 80 AU). The disk mass amounts to 2.67 × 10-3 M and the outer radius is 144 AU. The flux density in Jy/as2 is color coded for each individual channel. The disk is seen face-on (0°). The gap is barely visible. The line is optically thick owing to the material above the gap. However, the flux density is comparable to the model with 100 times lower mass (Fig. 2).

Open with DEXTER

thumbnail Fig. 4

Ideal 12C18O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet (at 80 AU). The disk mass amounts to 2.67 × 10-3 M and the outer radius is 144 AU. The flux density in Jy/as2 is color coded for each individual channel. The disk is seen face-on. Mostly because of the low abundance of this CO isotopologue, this line traces deeper parts of the disk. Thus, the gap is clearly detectable in contrast to the optically thick 12C16O line of this model (cf. Fig. 2).

Open with DEXTER

thumbnail Fig. 5

Ideal 12C12O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet. The flux density in Jy/as2 is color coded for each individual channel. The inclination of the disk is 10° and the disk has a total mass of 2.67 × 10-5 M. The gap is visible but does not appear as a symmetric ring like in the continuum case because the parts of the disk contributing in the individual channel maps change due to Doppler shift.

Open with DEXTER

The most massive disks (2.67 × 10-1 M) in our sample do not show any signs of a gap in any of the considered transitions (see Table 3). For the reasons outlined above, for a specific molecule observation with the aim to reveal gaps, we recommend estimating the optical thickness of the disk (for example from its mass derived through continuum observations) in order to choose the appropriate molecule and (optically thin) transition.

Disk inclination

In contrast to the ideal face-on models, the velocity-channel maps of a disk with a small inclination of 10° is shown in Fig. 5. Negative velocities (top left corner) belong to parts of the disk moving toward the observer and areas with positive (bottom right) velocities are moving away from the observer. In this figure the main problem for detecting gaps in inclined disks, and thus the majority of potential targets, becomes obvious. Only fractions of the gap can be seen, but it does not appear as a symmetric ring as the gas covers a broad range of velocities (vturb) along the gap throughout the individual velocity channels. This effect is illustrated in Fig. 6. The underlying disk model (Fig. 6, left) has an outer radius of 117 AU and the planet is located at 65 AU. Owing to the disk inclination of 10°, the line is significantly Doppler shifted. The velocity-channel map shows a cavity in the orbital region of the planet. However, this not an irrevocable sign of a gap. Owing to the projected velocity in the observer’s direction, the undisturbed model (see Fig. 6, right) also reveals a prominent local minimum in the same disk region. Thus, the velocity patterns effectively hide the gap and hinder unambiguous identification, especially for optically thick lines. A possible solution is to use velocity integrated intensity maps for gap detection instead (see Fig. 7). For low disk inclinations, these images do not deviate much from face-on inclined disks, i.e., for an inclination of 10° as considered in this study the flux deviation compared to the disk model seen face-on is less than 5%. The velocity integration has the advantage that the patterns disappear and the gap can be identified easily. However, the velocity dependency of the gap depth (see previous section) is lost. Consequently, the velocity-integrated maps are dominated by the line flux at the line center, which traces only the upper disk layers in the case of optically thick lines.

Central star

We now discuss the impact of using the Herbig Ae star, which is ~45 times more luminous, as the central star.

As the disk temperature is higher than the T Tauri case, the molecules are exited to higher energy states. The line radiation of the outer parts of the disk increases because a larger fraction of molecules is sufficiently exited to the relevant energy levels due to the warmer disk. In contrast to the situation in the outer disk regions, the higher temperature in the innermost region results in a depopulation of the lower energy states. Consequently, the low line transitions we consider in this study produce less flux in this region.

Second, the mass of the Herbig Ae star is about 3.5 times higher, resulting in a higher Keplerian velocity at a given radial distance from the central star. Thus, for the inclined disks the linewidth of the Doppler broadened line is much higher than in the T Tauri star case. Consequently, at a certain velocity channel (bandwidths of 50 m/s in this study) the parts of the disk with the appropriate velocities to contribute to this channel appear precise and well defined (see Fig. 8). Thus, assuming the same disks for both PMS stars, individual disk regions can be better studied in the case of the more massive central star, which increases the chances of gap detection. This effect is a direct result of the dependence of the linewidth on the kinetic temperature (cf. Eq. (B.8)).

Detectability of gaps in the ideal case

We describe our approach to identifying a gap and to deriving its depth in the velocity-channel maps.

We follow the procedure outlined by Ruge et al. (2013). We derive the contrast between the inner/outer edge of the gap and the gap itself by using radial cuts throughout the simulated maps (see Fig. 9). To reduce the probability of false detection (e.g., local minimum in the brightness distribution due to the inclination effect discussed in Sect. 4.1), 64 individual radial cuts are considered in the analysis. Only if a gap is indicated along more than 50% of successive radial cuts is it considered as detected. We successfully tested our algorithm with smooth undisturbed reference hydrodynamic models to verify that no substructure is misleadingly considered as a gap. However, because of the very general approach due to the huge parameter space, we emphasize that for specific observations better results may be achieved by adjusting the algorithm to the specific observation (e.g., a well-adapted compensation of background noise/sparse uv-coverage, which leads to asymmetric structures).

Overview charts

In this section, the feasibility of detecting gaps through observations of selected molecular lines for all considered disk configurations (see Table 4) is summarized. For this purpose, we derive the gap depth from the ideal velocity-channel maps for all five molecules and 32 assumed transitions using the algorithm described above. For every transition we derive the expected ALMA background noise of the reconstructed velocity-channel maps (σ) using the ALMA sensitivity calculator (see Sect. 4.2.1, Table 3). Owing to the atmospheric transmission, σ mostly rises with frequency resulting in a worse signal-to-noise ratio (S/N). To take care of this effect, we derive the gap depth and weight it by the expected sensitivity of ALMA at every line transition considered in this study. Thus, we obtain a criterion of how well individual transitions allow gap detection with ALMA. These results allow those transitions to be chosen that are best suited for the simulated ALMA observations in Sect. 4.2.

Figure 10 shows the results for an exemplary disk with an outer radius of 144 AU and a disk mass of 2.67 × 10-3 M around the Herbig Ae star. The gap depth measured in every velocity channel for the 12C16O (32), 12C18O (32), and CS (76) transition is shown. All transitions are optically thick at the line center. While the 12C16O (32) and CS (76) transitions fail to detect the gap at the line center, the optically thin 12C18O (32) transition allows for gap detection even at the central frequency. Nevertheless, the highest gap depth is achieved by the CS (76) transition at the line wings (at ±350 m/s).

To summarize our findings for all molecules, we search for the transition that provides the highest gap depth for given disk mass and disk size and present the result in one overview map. The resulting charts for both considered PMSs are shown in Fig. 11. The full overview charts showing the individual gap depth of all 32 transitions considered in this study can be found in Appendix A.

As already indicated in Sect. 4.1, the most promising molecule for identifying gaps in low-mass disks (~10-4 − 10-5 M) is 12C16O. The gap can be detected either at the line center for the lowest disk masses or at the line wings for the intermediate-mass disks. In most configurations the flux density and consequently the gap depth is higher even in the line wings compared to the gap depth which can be achieved with other molecules.

For the disks with intermediate masses (~10-2 − 10-3 M) of our sample gap detection becomes very challenging. In these cases 12C16O is optically thick for almost all considered disk configurations even at the line wings. Thus, transitions of less abundant molecules like HCO+, HCN, or CS are better suited for gap detection. However, because of the low abundance of these molecules the net flux is very low. Consequently, the gap depth that relates the flux density to the background noise is too low for gap detection in many cases. The most massive (~10-1 M) disks in our sample do not show any reliable sign of a gap for any molecule transition considered. We find that the gaps of these disks can be traced best outside the line in the dust continuum, even with the low bandwidth considered in this study.

When comparing the molecular line emission from the disk around both PMSs, we find very similar results in terms of the probability of gap detection. The effects due to their different stellar mass and luminosity (see Sect. 4.1 for a discussion) lead to slightly different molecules, which offer the best probability of detecting a gap (e.g., because of different excitation states 12C16O allows for gap detection even for the disks with a mass of ~10-4 M).

Therefore, our results demonstrate the feasibility of gap detection using molecular lines, but it is only efficient for low-mass protoplanetary disks or disks with huge inner cavities (e.g., transitional disks or circumbinary disks).

thumbnail Fig. 6

HCO+ (43) velocity-channel map at 150 m/s. The disk has a slight inclination of 10°. Disturbed (left) and smooth disks reveal similar results. Both models show indications of a gap. Thus, unambiguous gap detection requires detailed analysis of the origin of those structures.

Open with DEXTER

thumbnail Fig. 7

12C16O (32) velocity integrated intensity map of the channel map shown in Fig. 5. Owing to the integration, the patterns due to the underlying velocity field have vanished.

Open with DEXTER

thumbnail Fig. 8

Ideal 12C16O (32) velocity-channel map of a Herbig Ae disk model with an embedded Jupiter-mass planet. The disk mass amounts to 2.67 × 10-5 M and the outer radius is 144 AU. The flux density in Jy/as2 is color coded for each individual channel. The inclination of the disk amounts to 10°.

Open with DEXTER

thumbnail Fig. 9

Illustration of the method of detecting gaps. The blue line indicates a radial cut through the HCO+ (43) velocity-channel map (Fig. 14) at the central frequency (v = 0.0 m/s). For the simulated ALMA observations, the gap is considered as being detected if the depth of the gap is deeper than three times the noise level σ of the selected transition. In this particular case, the gap is clearly detectable, resulting in a gap contrast of ~7σ.

Open with DEXTER

thumbnail Fig. 10

Depth of a gap for three different molecules/transitions over velocity and weighted by the ALMA noise level σ. The disk has an outer radius of 144 AU and the disk mass amounts to 2.67 × 10-3 M and the disk heated by the Herbig Ae star. For this disk model, all considered transitions are optically thick at the line center and the detection of the gap is only feasible with 12C18O due to its low abundance. Nevertheless, at the line wings other molecules transitions allow for gap detection too, and in this particular case the gap depth has a maximum value when using the CS (73) transition.

Open with DEXTER

4.2. Simulated ALMA observations

4.2.1. ALMA simulation setup

thumbnail Fig. 11

Overview plot of the maximum gap depth in the ideal velocity-channel maps when comparing all transitions considered in this work. The central radiation source is a T Tauri (left) or a Herbig Ae star (right). For the most massive disks in our sample. All considered transitions are optically thick and only trace the surface of the disks. Thus, for the smallest disks the gap is detectable outside the line in the dust continuum (marked with Dust).

Open with DEXTER

thumbnail Fig. 12

12C16O (32) velocity-channel map of a model with an embedded Jupiter-mass planet as it could be observed using ALMA configuration 14 (1.6 km baseline). The inclination of the disk is 0° (face-on). The outer disk radius amounts to 144 AU and the disk mass to 2.67 × 10-5 M. The flux density in Jy/beam is color coded for each individual channel. In contrast to the ideal velocity-channel maps, complex structures are eliminated owing to beam convolution and noise. However, the gap can still be distinguished at the line slopes (±300 m/s).

Open with DEXTER

thumbnail Fig. 13

12C16O (32) velocity-channel map of the same model shown in Fig. 12, but with a disk inclination of 10°. The flux density in Jy/beam is color coded for each individual channel.

Open with DEXTER

We use the Common Astronomy Software Application package (CASA 4.2, Petry & CASA Development Team 2012) to simulate observations based on our ideal velocity-channel maps. We consider a total of six medium and extended configurations (baselines from 0.7 to 16 km) using 50 antennas, a total observing time of 3 h, a target distance of 140 pc, and the sky-coordinates of the “Butterfly Star” (α = 04h33m17sδ = + 22°53′20′′, Wolf et al. 2003; Gräfe et al. 2013) which is representative of young stellar objects in the Taurus-Auriga star-forming region. The precipitable water vapor (PWV) in the atmosphere has a major influence. Therefore, we also include thermal noise in our simulations and use the recommended column densities for the line transitions we consider in our sample (see Table 3). First, we employ the simobserve task to generate realistic visibility measurements. The final simulated observations are obtained using the simanalyze task which uses the CLEAN-algorithm for image reconstruction (Högbom 1974). The expected noise levels σ are estimated with the ALMA sensitivity calculator3 because the measured noise in the simulated images strongly depend on the quality of the deconvolution4. The properties of the considered molecules and the expected sensitivity are summarized in Table 3. We consider the ALMA polarization configuration in dual-mode. With respect to the results from Sect. 4.1 we select those molecules that offer the highest probability of gap detection, i.e., 12C16O(32), 12C18O(32), CS (76), HCO+(32), and HCN (32). In addition, their corresponding transitions around 345 GHz (ALMA band 7, see Table 3), offer a reliable comparability of the different molecules and yield a good compromise between spatial resolution (higher frequencies) and better atmospheric transmission (lower frequencies)5. The 12C16O (43) transition is added to provide a second 12C16O transition offering higher spatial resolution at the cost of sensitivity.

4.2.2. Simulated line observations

In Figs. 12 and 13 a face-on and an inclined disk model are shown, respectively, as they could be observed with ALMA using configuration 14 (~1.6 km baseline) within three hours (corresponding ideal maps are shown in Figs. 2 and 5, respectively). More compact baselines result in higher flux densities and consequently higher S/Ns, while larger baselines allow gaps of compact disks to be traced at the cost of decreasing S/N. We find for the most extended disk configurations (Rin ~ 45 AU, Rout ~ 200 AU) in our model setup, maximum baselines of about 1 km are sufficient to identify the gap in these disks, but fail for the smaller disks.

thumbnail Fig. 14

HCO+ (43) velocity-channel map of a model with an embedded Jupiter-mass planet as it could be observed using ALMA configuration 14 (1.6 km baseline). The flux density in Jy/beam is color coded for each individual channel. The inclination of the disk amounts to 0° (face-on). Compared to the same model but simulated 12C16O (32) observation, the line flux is about 35% lower owing to the distinct lower abundance of HCO+. However, because of the low abundance the gap is clearly detectable, even at the central wavelength.

Open with DEXTER

For the compact disks in our parameter space, higher spatial resolution and thus longer baselines are required. In our sample the smallest disks have an inner/outer radius of 14 AU/63 AU, respectively. We find that unambiguous gap detection is feasible using intermediate baselines of about 13 km. Longer baselines are not required and only result in sparser uv-coverage and thus a lower S/N of the reconstructed maps. Furthermore, many of the smaller planetary signatures (inner and outer spiral arms and disk inhomogenities) are eliminated/smoothed out because of beam convolution and worse uv-coverage. However, for this particular model, and for most of the expanded models in our sample, the gap can be identified in both inclination cases using ALMA configuration 14.

For the same disk, the resulting HCO+ (4-3) velocity-channel map for ALMA configuration 14 is shown in Fig. 14. The HCO+ molecule is less abundant than 12C16O and consequently the line is less optically thick. Therefore, regions closer to the disk midplane can be observed and the gap appears more pronounced. In contrast, the line flux itself is about 35% weaker than the simulated 12C16O (32) line of the same model.

4.2.3. Gap detection results

Our main intention is to investigate the feasibility of tracing gaps in molecular lines with ALMA. To allow this for a large range of disk models, we define a general criterion for a quantitative measure of gap detection in the simulated ALMA maps. We apply the same algorithm as in the case of ideal maps (see Sect. 4.1). The resulting gap depth is given in units of σ to clarify the meaning of the gap detection. Only if the depth of the gap is greater than three times the noise level σ do we consider the gap as detected.

thumbnail Fig. 15

Measured gap depth in units of σ over velocity. Gap detection with the optical thick 12C16O (32) transition (blue) is feasible, but only at the line edges (e.g., ± 700 m/s). In contrast, the less abundant 12C18O (32) transition (green) allows for gap detection at the line center. However, the CS (76) gives the most significant gap depth of about 17σ at ± 350 m/s.

Open with DEXTER

First, we demonstrate that the optical depth effect (see Sect. 4.1 for a discussion of the ideal case) has direct observational impact. As shown in the ideal case, the choice of the molecule for a given disk mass/size is most important for unambiguous gap detection. In Fig. 15 we present the gap depth for the same transitions as in Fig. 10, but now for a simulated ALMA observation utilizing configuration 14 (~1.6 km baseline). Both results are closely comparable as the gap depth in the ideal case shows the same trend as the simulated ALMA observations, i.e., 12C16O (32) cannot be used as a tracer for gaps, and 12C18O (32) allows for gap detection even at the line center, but the CS (76) transition traces the gap best and results in a gap depth of about 17σ at ±350 m/s.

Second, our results are summarized in the overview charts for the Herbig Ae (left) and T Tauri (right) star (Fig. 16). Compact ALMA configurations (9, 11) result in high S/Ns, but lack of spatial resolution. As a consequence, utilizing the most compact ALMA configuration nine gaps are only identifiable for the most extended models of our sample (e.g., disk models with Rout = 255 AU). If we consider longer baselines (higher spatial resolution) in the range of ~1.62.3 km we are able to identify gaps for more compact models with a σ> 3. Intermediate baselines (e.g., configuration 1114, ~1.01.6 km) potentially allow gaps to be traced for the wide range of models considered in this study.

thumbnail Fig. 16

Overview maps of the detectability of gaps around Herbig Ae (left) and T Tauri (right) stars. Each column represents one ALMA configuration and each row represents one molecule transition. The x-axis of the inner plot gives the total mass of the disk, and on the y-axis the scaling factor k is linearly connected to the spatial dimension of the disks. The depth of the gap is color coded in units of σ. We only consider a gap as detected when the depth is at least three times the background noise σ. For both stars, gap detection is feasible for a wide range of models in our setup within three hours of observation. The difference between the two stars is mostly due to the much higher luminosity. Similar to the results in the ideal case, gap detection is feasible but is restricted to low-mass disks (Mdisk ≤ 10-3). Gap detection using the most extent configuration (max. baseline ~16 km) is not feasible within three hours of observation, which is mostly due to the small bandwidth of 50 m/s and is not shown in this overview plot.

Open with DEXTER

The comparison between the molecules considered shows that CS, HCN, and HCO+ allow the detection of gaps even for disks with intermediate masses (~10-3 M), while 12C16O tends to fail. As already discussed in Sect. 3.2, owing to the optical thick 12C16O line only the surface of the disks can be observed, which is barely perturbed by the planet-disk interaction. In contrast, for lines of molecules with lower abundances the surface is still optically thin even for more massive disks. We find that especially CS and HCN are proper candidates to serve as a tracer for gap detection in the mass range from ~10-3 M to 10-4 M.

For more massive disks (Mdisk ≥ 10-2 M) the 12C16O isotopologue 12C18O is an adequate candidate. However, owing to the resulting low line flux gaps these massive disks should be observed in the continuum because in the submm wavelength range provided by ALMA the disk is still optically thin in the gap region (see Ruge et al. 2013).

Finally, we find that because of the assumed observation time of 3 h, long-baseline configurations (20, 28) are not suited to identifying gaps because of the low intrinsic flux and inadequate uv-coverage.

5. Summary

Our intention was to explore how planet-induced structures, in particular planetary gaps, appear in high-angular resolution molecular line maps, and subsequently if ALMA would allow the detection of these structures. For this reason, we developed the line radiative transfer code called Mol3D. It allows the calculation of the level-populations with several N-LTE approximate methods and the production of ideal velocity-channel maps. First, we tested the new code against other established (molecular line) RT codes.

Based on 3D hydrodynamic simulations carried out with the PLUTO code, we simulated synthetic velocity-channel maps for five molecules (12C16O, 12C18O, HCO+, HCN, and CS) for a total of 32 rotational transitions. We summarize our findings:

  • The patterns of the global (Keplerian) velocity field in thecorresponding channel maps of inclined disks maskplanet-induced structures. Gap identification can be ambiguousin these cases.

  • Molecules that tend to become optically thick (like 12C16O) only trace the upper layers of the disk, which are rarely affected by the planet-disk interaction. In consequence, gap detection for massive disks (10-2 M) is barely feasible with molecular lines.

  • For lower mass disks, depending on the molecule and its optical properties (abundances, excitation state) the gap can be observed for a broad range of considered disk configurations, not necessarily at the central wavelength of this transition but at the line wings. Because of the higher optical depth of the line at its center, the τ = 1 surface traces only the upper unperturbed disk layers, while the τ = 1 surface of the line wings traces deeper and thus more strongly perturbed layers.

  • With ALMA it is feasible to identify gaps using molecular lines for several disk configurations, if one respects the effects we outlined above (i.e., by choosing an appropriate optically thin molecule and corresponding transition). In particular, for low-mass disks (10-3 M) the feasibility of identifying gaps using molecular lines is very promising. For this mass range molecular line observations fill the discovery space for gaps where continuum observations are supposed to fail (see Fig. 8 in the complementary study of Ruge et al. 2013).

  • Our simulations are based on an assumed distance of 140 pc and a fixed sky coordinate (max. elevation ~45°) to give reliable results for the Taurus star-forming region. Thus, for many disks located closer or at more favorable positions on the sky the feasibility of detecting gaps is expected to be even higher.

We emphasize that the observation of molecules with different abundances allow for unambiguous identification of gaps and the physical and chemical composition of protoplanetary disks. Especially for low-mass disks where gap detection in the continuum is difficult owing to the low intrinsic flux, molecular line observations are a promising alternative with which to identify gaps.

Finally, we note that our study is based on HD simulations. In addition, Uribe et al. (2011) have remarked that magnetohydrodynamics (MHD) simulations result in more extended gaps owing to the presence of the weak magnetic field. Including intrinsic stellar radiation, the form of the gaps also changes because the gap walls can be directly radiated by the host star, resulting in deeper gaps. This effect is even higher for the outer wall (Jang-Condell & Turner 2012). Using smoothed particle hydrodynamic (SPH) simulations, Gonzalez et al. (2012) found that gaps are deeper for larger dust grains and bigger planets. In addition, we only considered a Jupiter-mass planet in this study. More massive planets would enlarge the gap width owing to the higher gravitational potential. If these effects are present it would therefore be even easier to observe gaps with ALMA than proposed within this study.


2

Harten-Lax-van Leer-Contact (see Toro 2009).

Acknowledgments

This presented work is funded through the DFG grant: WO 857/10-1. Special thanks go to Yaroslav Pavlyuchenkov for his comments and advice, which were very helpful in the early development phase of Mol3D. We also want to acknowledge Peter Scicluna and Tommy Breslein for language editing and Jan Philipp Ruge for the helpful advice, especially regarding the use of the CASA simulator.

References

Online material

Appendix A: Gap depth complete overview charts

The full overview charts for all molecules considered in our study are summarized in Fig. A.1 for the T Tauri case and in Fig. A.1 for the Herbig Ae case. The derived σ-weighted gap depth is color coded for each considered disk configuration.

thumbnail Fig. A.1

Overview plot of gap depth in the ideal velocity-channel maps for five different molecules. The central radiation source is a T Tauri star. The σ-weighted gap depth is color coded.

Open with DEXTER

thumbnail Fig. A.2

Overview plot of gap depth in the ideal velocity-channel maps for five different molecules. The central radiation source is a Herbig Ae star. The σ-weighted gap depth is color coded.

Open with DEXTER

Appendix B: Mol3D

In this section we present a detailed derivation of the underlying methods, physics, and assumptions used by our new N-LTE line and dust continuum radiative transfer code Mol3D.

Appendix B.1: Radiative transfer equations

We start with the transfer equation in the following form, (B.1)with the monochromatic intensity Iν, the emission factor jν(s), and the total absorption factor αν(s). The intensity from point p1 to p2 along its path s can be calculated using the integral form of Eq. (B.1), (B.2)where τν is the optical depth between two points, p1 and p2, along the ray of sight defined as (B.3)

Appendix B.2: The gas-phase

For the gas, e.g., for a selected transition from lower level i to upper level j, the emission (jν) and absorption (αν) coefficients have the following form: We follow the ansatz of Rybicki & Hummer (1991) for the non-overlapping multi-level line treatment. For a molecule with N levels and a selected transition between lower i and upper j level, we consider the spontaneous downward transition rates Aij, the Einstein coefficients Bij for stimulated transition, and the collision rate Cij. The collision rate can be calculated from the number density of the collision partner and the downward collision rate coefficient γij, which is the Maxwellian average of the collision cross section σ measured in laboratory: (B.6)All molecule data we use in our software package (e.g., frequencies, transition rates Aij, collision rate coefficients) are obtained from the Leiden Atomic and Molecular Database (LAMDA) (Schöier et al. 2005). As a consequence, Mol3D uses the same common input format for the molecule data as many other line radiative transfer programs available like RADEX (van der Tak et al. 2007) or URAN(IA) (Pavlyuchenkov & Shustov 2004). This approach offers the advantage of easily extending the code for several common molecules available in the LAMDA database (3 atomic and 33 molecular species at the time of writing). The dominant broadening effects are thermal and turbulent broadening. Hence, we assume a Gaussian line profile function, but in principle any other profile function can be included as well, (B.7)where νij is the central frequency of the transition and c the speed of light. The parameter vtot gives the total line broadening factor taking into account the assumed microturbulent velocity vturb and the kinetic velocity vkin due to the kinetic temperature Tkin: (B.8)In our code, the microturbulent velocity vturb can either be set to a common value, ~0.10.2 km s-1 for protoplanetary disks (Piétu et al. 2007; Hughes et al. 2011), or for example in the case of an outcome of (M)HD/MRI simulations, it can be provided as an input for every grid cell.

Finally, to describe the full line radiative transfer problem, a set of balance equations is needed: (B.9)In fact, Eq. (B.9) is not one single equation, but a whole set of linear equations, one for every energy level. In our code this linear matrix equation is solved with a simple Gaussian elimination algorithm. This set has to be solved locally (for every grid cell), but it has a global character due to its dependency on the line integrated mean intensity defined as (B.10)Equations (B.9) and (B.10) are directly coupled to Eq. (B.1).

To solve the line radiative transfer problem, one needs to start with estimated values for the level populations, solve the radiative transfer equation, calculate the mean intensity, and then solve Eq. (B.9) again to create a new set of level populations. Then one needs to iterate until the correct level populations have been found.

Appendix B.2.1: Level populations

The main problem in radiative line transfer is the calculation of the level populations. Equation (B.9) shows that the level populations are directly coupled to the mean intensity and therefore to the gas temperature, which is unknown for complex structures as in the case of protoplanetary disks.

We implemented dissimilar algorithms to calculate the level populations with three approximate methods, for example, local thermodynamical equilibrium LTE, FEP, and LVG. As Pavlyuchenkov et al. (2007) have shown for the case of protoplanetary disks, for the most common molecules LVG is a very accurate method of calculating the level populations and gives a very good ratio between computational efficiency and reliability of the results.

LTE level populations

In this assumption the level populations follow a Boltzmann distribution and the excitation temperature of the selected transition is assumed to be equal to the kinetic temperature Texc = Tkin: (B.11)Here gi and gj are the statistical weights of the ij line transition and νij the corresponding central frequency.

LVG level populations

In the case of protoplanetary disks, the radial velocity gradients are usually much larger than the local thermal and microturbulent velocities (Weiß et al. 2005; Castro-Carrizo et al. 2007). This means that photons emitted at a certain disk region can only interact and consequently get absorbed locally. In this approach, the mean intensity J is approximated from the local source function S in combination with an external radiation field Jext: The quantity β gives the probability of a photon beeing able to escaping the model. We note that in the optically thin case a photon can escape the model without interaction. For this case, β equals 1 and the LVG method would be equal to the full escape probability method. Therefore, the FEP method (also included in Mol3D) uses the same formalism as the LVG method, but with β always set to 1 (see also van der Tak et al. 2007). To calculate the local β, our program uses Eq. (B.14) introduced by Mihalas et al. (1978); de Jong et al. (1980) for accretion disks in general, (B.14)where τ is the effective optical depth of the observed line. We note that this formula is well suited in the case of proto-circumstellar disks. For other geometries and scenarios other approximations might give better results (see, e.g., de Jong et al. 1975 or Osterbrock & Ferland 2006). At this point an approximation of the optical line depth is made. As it is connected to the local velocity field V, microturbulence vturb, and gas density, we use the following expression (see also Pavlyuchenkov et al. 2007): (B.15)

Appendix B.3: Benchmarks

We demonstrate the applicability of the code Mol3D. For this purpose we compare its results to those obtained with the line radiative transfer code URAN(IA) (Pavlyuchenkov & Shustov 2004) and the dust continuum radiative transfer code MC3D (Wolf et al. 1999; Wolf 2003).

Appendix B.3.1: URAN(IA)

URAN(IA) is a 2D radiative transfer code for molecular lines. It features a Monte Carlo algorithm for calculating the mean intensity and an accelerated lambda iteration (ALI) method for self-consistent molecular excitation calculations. This code has been extensively tested in 1D and 2D, for example against the RATRAN code (Hogerheijde & van der Tak 2000). It successfully passed all tests formulated in the benchmark paper for N-LTE molecular radiative transfer by van Zadelhoff et al. (2002). It has been used in several applications, for example modeling the starless core L1544 (Pavlyuchenkov & Shustov 2004), and in dedicated parameter space studies of molecular line formation of molecular line formation of prestellar cores (Pavlyuchenkov et al. 2008).

In this code, several approximate methods of calculating the level populations are included. In a comprehensive study, Pavlyuchenkov et al. (2007) have shown the reliability of these methods in protoplanetary disk scenarios.

In Sect. 2 we introduced the different methods of calculating level populations implemented in Mol3D. In Fig. B.1 we present the spectrum of the HCO+ (43) transition of a typical protoplanetary disk (Mdisk = 0.07 M), comparing the different methods assuming a uniform N(HCO+)/N(H) ratio of 1 × 10-8. For simplicity, we choose a flared disk model around a T Tauri star utilizing a temperature distribution with a gradient in radial and vertical direction and a power-law model for the density distribution adapted from Shakura & Sunyaev (1973): (B.16)The disk is orientated face-on. We assume pure Keplerian rotation and neglect dust re-emission. All spectra are normalized to the maximum value of the LTE solution. To allow comparison, the resulting spectra obtained with the URAN(IA) code are also included.

URAN(IA) and Mol3D produce comparable HCO+ (43) spectra with differences of less than about 0.5%. The most realistic solution (black line) is obtained using the accelerated Monte Carlo (ART), to classify the approximation methods. As discussed in greater detail in Appendix B.2.1, the LTE method overestimates the net flux and the FEP method underestimates it. The fluxes in the line obtained with the LVG method are in between those obtained with the other methods and are closely comparable to the ART method. This result has also been found by Pavlyuchenkov et al. (2007) and we refer the interested reader to their study.

thumbnail Fig. B.1

HCO+ (43) transition of a protoplanetary disk, oriented face-on. Different colors are obtained with different methods. Solid lines are obtained using the URAN(IA) software package and dotted lines represent solutions obtained with Mol3D. Shown is the normalized flux over velocity. Both codes produce comparable results within acceptable differences (less than 0.5%). LVG is a good approximation in the case of circumstellar disks (Pavlyuchenkov et al. 2007).

Open with DEXTER

Appendix B.3.2: MC3D

MC3D is a 3D continuum radiative transfer code that uses a Monte Carlo method to calculate self-consistent dust temperature distributions, continuum re-emission/scattering maps, and SEDs. Among a wide variety of radiative transfer studies it has been applied extensively for the analysis of multiwavelength high-angular-resolution observations of circumstellar disks (e.g., Schegerer et al. 2008; Sauter & Wolf 2011; Madlener et al. 2012; Gräfe & Wolf 2013). MC3D was successfully tested against other continuum RT codes (e.g., Pascucci et al. 2004).

In this section, we compare Mol3D with MC3D for a typical protoplanetary disk scenario. The density distribution is described by Eq. (B.16), assuming a total disk mass of 0.04  M. The temperature structure is calculated self-consistently with

each code using their Monte Carlo scheme. Both codes use a similar approach based on the assumption of local thermal equilibrium and immediate temperature correction as proposed by Bjorkman & Wood (2001).

thumbnail Fig. B.2

Midplane temperature of a typical protoplanetary disk around a pre-main-sequence T Tauri star. This disk inner radius amounts to 2 AU and the outer radius to 200 AU, respectively. The total disk mass amounts to 0.04  M. Mol3D (red) and MC3D (blue) obtain comparable results. In this case, the differences depend significantly on the optical disk properties and subsequently on the number of photon packages simulated in the Monte Carlo process. Thus, the error is higher in the dense region at the inner rim of the disk (max. 10%) and lower at the outer disk regions (max. 2%).

Open with DEXTER

The resulting midplane temperature around a PMS T Tauri star is shown in Fig. B.2. The calculation of the temperature in the dense midplane is one of the most crucial RT problems, because these disk regions can hardly be reached by stellar photons directly. Thus, it is mostly heated indirectly by thermal re-emission radiation from the innermost disk regions and upper disk layers. Both codes produce very smooth and similar temperature distributions. The statistical nature of the Monte Carlo method results in deviations of about 10% for the innermost disk parts and less than 2% for the outer regions, which is mainly due to the probability that the photon packages can reach these regions.

We also compared the dust re-emission/scattered light maps and SEDs and find that both codes produce comparable results with maximum deviations of about 10%.

All Tables

Table 1

Comparison of grid parameter space of both codes.

Table 2

Considered pre-main-sequence stars and their properties.

Table 3

Properties of the different molecules considered in this study.

Table 4

Parameter space of the disk model.

All Figures

thumbnail Fig. 1

Example of HD simulations with embedded Jupiter-mass planet (planetmass/starmassratio = 10-3). Color-coded is the H2 number density in log scale. Left: midplane density. The planet has cleared out its orbit, producing a prominent gap. Owing to the perturbations, complex structures like spiral arms are also visible. Right: cut perpendicular to the disk midplane. The gap is clearly visible in the gas distribution. We note that the upper layers of the disk have nearly the same density as in the undisturbed case.

Open with DEXTER
In the text
thumbnail Fig. 2

Ideal 12C16O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet with an orbital distance of 80 AU. The disk mass amounts to 2.67 × 10-5 M, which is the lowest disk mass in our sample, and the outer radius is 144 AU. The flux density in Jy/mas2 is color coded for each individual channel. The disk is seen face-on (0°). The gap is clearly visible, particularly in the line wings. Even more complex features like the huge spiral arm (e.g., at 300 m/s) can be recognized. The total width of the line is only given by thermal and microturbulent velocities (Sect. B.2).

Open with DEXTER
In the text
thumbnail Fig. 3

Ideal 12C16O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet (at 80 AU). The disk mass amounts to 2.67 × 10-3 M and the outer radius is 144 AU. The flux density in Jy/as2 is color coded for each individual channel. The disk is seen face-on (0°). The gap is barely visible. The line is optically thick owing to the material above the gap. However, the flux density is comparable to the model with 100 times lower mass (Fig. 2).

Open with DEXTER
In the text
thumbnail Fig. 4

Ideal 12C18O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet (at 80 AU). The disk mass amounts to 2.67 × 10-3 M and the outer radius is 144 AU. The flux density in Jy/as2 is color coded for each individual channel. The disk is seen face-on. Mostly because of the low abundance of this CO isotopologue, this line traces deeper parts of the disk. Thus, the gap is clearly detectable in contrast to the optically thick 12C16O line of this model (cf. Fig. 2).

Open with DEXTER
In the text
thumbnail Fig. 5

Ideal 12C12O (32) velocity-channel map of a T Tauri disk model with an embedded Jupiter-mass planet. The flux density in Jy/as2 is color coded for each individual channel. The inclination of the disk is 10° and the disk has a total mass of 2.67 × 10-5 M. The gap is visible but does not appear as a symmetric ring like in the continuum case because the parts of the disk contributing in the individual channel maps change due to Doppler shift.

Open with DEXTER
In the text
thumbnail Fig. 6

HCO+ (43) velocity-channel map at 150 m/s. The disk has a slight inclination of 10°. Disturbed (left) and smooth disks reveal similar results. Both models show indications of a gap. Thus, unambiguous gap detection requires detailed analysis of the origin of those structures.

Open with DEXTER
In the text
thumbnail Fig. 7

12C16O (32) velocity integrated intensity map of the channel map shown in Fig. 5. Owing to the integration, the patterns due to the underlying velocity field have vanished.

Open with DEXTER
In the text
thumbnail Fig. 8

Ideal 12C16O (32) velocity-channel map of a Herbig Ae disk model with an embedded Jupiter-mass planet. The disk mass amounts to 2.67 × 10-5 M and the outer radius is 144 AU. The flux density in Jy/as2 is color coded for each individual channel. The inclination of the disk amounts to 10°.

Open with DEXTER
In the text
thumbnail Fig. 9

Illustration of the method of detecting gaps. The blue line indicates a radial cut through the HCO+ (43) velocity-channel map (Fig. 14) at the central frequency (v = 0.0 m/s). For the simulated ALMA observations, the gap is considered as being detected if the depth of the gap is deeper than three times the noise level σ of the selected transition. In this particular case, the gap is clearly detectable, resulting in a gap contrast of ~7σ.

Open with DEXTER
In the text
thumbnail Fig. 10

Depth of a gap for three different molecules/transitions over velocity and weighted by the ALMA noise level σ. The disk has an outer radius of 144 AU and the disk mass amounts to 2.67 × 10-3 M and the disk heated by the Herbig Ae star. For this disk model, all considered transitions are optically thick at the line center and the detection of the gap is only feasible with 12C18O due to its low abundance. Nevertheless, at the line wings other molecules transitions allow for gap detection too, and in this particular case the gap depth has a maximum value when using the CS (73) transition.

Open with DEXTER
In the text
thumbnail Fig. 11

Overview plot of the maximum gap depth in the ideal velocity-channel maps when comparing all transitions considered in this work. The central radiation source is a T Tauri (left) or a Herbig Ae star (right). For the most massive disks in our sample. All considered transitions are optically thick and only trace the surface of the disks. Thus, for the smallest disks the gap is detectable outside the line in the dust continuum (marked with Dust).

Open with DEXTER
In the text
thumbnail Fig. 12

12C16O (32) velocity-channel map of a model with an embedded Jupiter-mass planet as it could be observed using ALMA configuration 14 (1.6 km baseline). The inclination of the disk is 0° (face-on). The outer disk radius amounts to 144 AU and the disk mass to 2.67 × 10-5 M. The flux density in Jy/beam is color coded for each individual channel. In contrast to the ideal velocity-channel maps, complex structures are eliminated owing to beam convolution and noise. However, the gap can still be distinguished at the line slopes (±300 m/s).

Open with DEXTER
In the text
thumbnail Fig. 13

12C16O (32) velocity-channel map of the same model shown in Fig. 12, but with a disk inclination of 10°. The flux density in Jy/beam is color coded for each individual channel.

Open with DEXTER
In the text
thumbnail Fig. 14

HCO+ (43) velocity-channel map of a model with an embedded Jupiter-mass planet as it could be observed using ALMA configuration 14 (1.6 km baseline). The flux density in Jy/beam is color coded for each individual channel. The inclination of the disk amounts to 0° (face-on). Compared to the same model but simulated 12C16O (32) observation, the line flux is about 35% lower owing to the distinct lower abundance of HCO+. However, because of the low abundance the gap is clearly detectable, even at the central wavelength.

Open with DEXTER
In the text
thumbnail Fig. 15

Measured gap depth in units of σ over velocity. Gap detection with the optical thick 12C16O (32) transition (blue) is feasible, but only at the line edges (e.g., ± 700 m/s). In contrast, the less abundant 12C18O (32) transition (green) allows for gap detection at the line center. However, the CS (76) gives the most significant gap depth of about 17σ at ± 350 m/s.

Open with DEXTER
In the text
thumbnail Fig. 16

Overview maps of the detectability of gaps around Herbig Ae (left) and T Tauri (right) stars. Each column represents one ALMA configuration and each row represents one molecule transition. The x-axis of the inner plot gives the total mass of the disk, and on the y-axis the scaling factor k is linearly connected to the spatial dimension of the disks. The depth of the gap is color coded in units of σ. We only consider a gap as detected when the depth is at least three times the background noise σ. For both stars, gap detection is feasible for a wide range of models in our setup within three hours of observation. The difference between the two stars is mostly due to the much higher luminosity. Similar to the results in the ideal case, gap detection is feasible but is restricted to low-mass disks (Mdisk ≤ 10-3). Gap detection using the most extent configuration (max. baseline ~16 km) is not feasible within three hours of observation, which is mostly due to the small bandwidth of 50 m/s and is not shown in this overview plot.

Open with DEXTER
In the text
thumbnail Fig. A.1

Overview plot of gap depth in the ideal velocity-channel maps for five different molecules. The central radiation source is a T Tauri star. The σ-weighted gap depth is color coded.

Open with DEXTER
In the text
thumbnail Fig. A.2

Overview plot of gap depth in the ideal velocity-channel maps for five different molecules. The central radiation source is a Herbig Ae star. The σ-weighted gap depth is color coded.

Open with DEXTER
In the text
thumbnail Fig. B.1

HCO+ (43) transition of a protoplanetary disk, oriented face-on. Different colors are obtained with different methods. Solid lines are obtained using the URAN(IA) software package and dotted lines represent solutions obtained with Mol3D. Shown is the normalized flux over velocity. Both codes produce comparable results within acceptable differences (less than 0.5%). LVG is a good approximation in the case of circumstellar disks (Pavlyuchenkov et al. 2007).

Open with DEXTER
In the text
thumbnail Fig. B.2

Midplane temperature of a typical protoplanetary disk around a pre-main-sequence T Tauri star. This disk inner radius amounts to 2 AU and the outer radius to 200 AU, respectively. The total disk mass amounts to 0.04  M. Mol3D (red) and MC3D (blue) obtain comparable results. In this case, the differences depend significantly on the optical disk properties and subsequently on the number of photon packages simulated in the Monte Carlo process. Thus, the error is higher in the dense region at the inner rim of the disk (max. 10%) and lower at the outer disk regions (max. 2%).

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.