RESOLVED SUBMILLIMETER OBSERVATIONS OF THE HR 8799 AND HD 107146 DEBRIS DISKS

We detect the CO(3–2) emission reported by Williams & Andrews (2006) and demonstrated by Su et al. (2009) to be at an LSR velocity coincident with that of HR 8799. The line emission is quite narrow, appearing in only two of the 0.7 km s⁻¹ channels, indicating a linewidth of ≲1.4 km s⁻¹. The line does not image cleanly in the interferometric data, exhibiting striping from northeast to southwest across the image. There are two likely reasons for this: (1) the position angle of the stripes matches that of the CO filament in Figure 4 of Williams & Andrews (2006), and the separation between stripes is equal to that between the peak of the synthesized beam and its largest (30%) sidelobes, and (2) there is a knot of bright CO emission visible in the James Clerk Maxwell telescope (JCMT) map to the northwest of the SMA map center at a distance of approximately the width of the SMA primary beam, which coincides with a peak in the SMA map. Not only are interferometers inherently poorly suited to imaging extended emission like the CO filament observed with the JCMT, but the presence of a bright source outside the primary beam on its own would be expected to cause striping across the image. Nevertheless, there is no evidence of a velocity gradient across the line, which is characteristic of rotation and would be expected if the CO(3–2) emission were associated with the disk. Given the spatial resolution of the SMA, such a gradient should be easily detectable if it were present. Assuming a gas disk spatially coincident with the observed dust emission, the CO(3–2) line should be concentrated well within the 30'' SMA primary beam and should exhibit centroid shifts across the line by ∼5'' between the red and blue channels, which is not observed. We therefore conclude that the CO emission does not originate from the disk around HR 8799, although the coincidence in spatial location and velocity are suggestive that the disk may be associated with the CO filament.

The bottom right panel of Figure 1 shows the HD 107146 SED assembled from the literature (Landolt 1983; Williams et al. 2004; Carpenter et al. 2005; Moór et al. 2006; Hillenbrand et al. 2008; Corder et al. 2009, plus IRAS and Two Micron All Sky Survey (2MASS) fluxes). We model the SED with three components: (1) a Kurucz model stellar photosphere with surface gravity log g = 4.5 and effective temperature T_eff = 5859 K (Carpenter et al. 2008; Hillenbrand et al. 2008), (2) a dust belt with a temperature of 69 K which is required to reproduce the mid-IR fluxes but does not contribute substantially to the millimeter-wavelength flux, and (3) an outer debris belt with the properties described below. It is this latter component that accounts for effectively all of the 880 μm flux, and it must therefore be modeled as a spatially extended component to reproduce both the observed visibilities and the SED.

As in Williams et al. (2004), we assume a dust grain emission efficiency Q_λ = 1 − exp[ − (λ₀/λ)^β], where λ₀ is a critical wavelength; this has the desired asymptotic behavior that Q_λ = (λ₀/λ)^β for λ ≫ λ₀ and Q_λ = 1 for λ ≪ λ₀, while varying smoothly between the two extremes. However, unlike Williams et al. (2004) who fit the SED with a single dust temperature, we require the dust grains to be in radiative equilibrium with the star and allow the disk to be spatially extended. With these assumptions, using the textbook formulation from Tielens (2005), the dust temperature is given by

$$\begin{eqnarray} &&T_r = \left[ \frac{L_*}{\sigma r^2} \frac{\pi ^3}{240} \frac{1}{(\beta +3)! \zeta (\beta +4) Q_0} \left(\frac{\lambda _0 k}{hc}\right)^{-\beta }\right]^\frac{1}{4+\beta }, \end{eqnarray} \tag{ 1 }$$

where L_* is the stellar luminosity, σ is the Stefan–Boltzmann constant, r is the distance from the star, β is the dust grain emission efficiency power-law index, ζ is the Riemann zeta function, and Q₀ is the dust grain efficiency at the critical wavelength λ₀ = 2πa. We assume a single characteristic grain size a, since this is sufficient to model the data and more complexity in the grain size distribution is not required to fit currently available data. The flux as a function of wavelength is then

$$\begin{eqnarray} &&F_\lambda =\frac{\pi a^2 Q(\lambda)}{d^2} \int _{r_{\mathrm in}}^{r_{\mathrm out}} \! 2 \pi r \, B_\lambda (T_r) \, n(r)\, {d}r, \end{eqnarray} \tag{ 2 }$$

where d is the distance to the source, r_in and r_out are the inner and outer radii of the dust disk, B_λ is the Planck function, and n is the number density of dust grains. The number density of grains is related to the surface mass density of dust grains as Σ = nm_g, where m_g is the mass per dust grain, and we parameterize Σ as a function of radius in the following way: Σ(R) = Σ₁₀₀(100 AU/r)^p. For consistency with previous works (e.g., Williams et al. 2004), we assume a standard millimeter-wavelength mass opacity of 1.8(ν/10¹² Hz)^0.8 cm² g⁻¹ (D'Alessio et al. 2001) evaluated at the 880 μm wavelength of the SMA data.

Because the outer dust belt is not bright enough at short wavelengths to reproduce the 24 and 33 μm Spitzer MIPS fluxes, we add a low-mass 69 K dust belt component to the SED (it would be entirely unresolved and contribute no detectable flux in the SMA data). The temperature is set to match that derived from the IR spectrum in Carpenter et al. (2009), and the mass is then varied to best reproduce the observed SED. We do not expect the properties of this belt to substantially affect the properties of the outer disk that is the focus of this analysis.

We also generate a model image assuming an inclination i = 25° and position angle of the major axis PA = 148° from the scattered light images (Ardila et al. 2004). We compare the spatially resolved model to the SMA observations in the visibility domain, using the MIRIAD task uvmodel to sample the model image at the same spatial frequencies as the SMA data after convolution with a 03 Gaussian seeing kernel.

To find the best-fit model parameters, we perform a χ² comparison between the model and the SED and SMA visibilities. We calculate a χ² value for the SED as a whole as well as for all of the real and imaginary visibilities at each spatial frequency and then simply sum the two values so that χ² = χ²_SED + χ_vis². As discussed in Andrews et al. (2009), due to a fortuitous balance between the quality and quantity of samples in the two data sets, each contributes roughly the same amount to the total χ² value, so that neither the SED nor the visibility data dominates the final fit. We vary the dust grain properties {a,β}, the spatial parameters {r_in,r_out,p}, and the total mass in the disk M_disk and in the belt M_belt. The total disk mass is related to the surface density as $M_\mathrm{disk} = 2\pi \Sigma _{100}(100\,{\rm AU})^p \frac{r_\mathrm{out}^{2-p} - r_\mathrm{in}^{2-p}}{2-p}$.

The six parameters effectively divide into two sets of three, with the spatial parameters {r_in,r_out,p} allowing us to reproduce the morphology of the spatially resolved data and the dust parameters {a,β,M_disk} allowing us to match the SED. The characteristic grain size a affects the dust temperature and shifts the peak wavelength of the blackbody. β has a similar effect, but its influence on Q_λ is far more visible in determining the overall shape of the dust spectrum, particularly the slope with which the tail of the curve falls at long wavelengths. For the optically thin disk models we consider, the total mass M_disk acts as a normalizing factor, shifting the blackbody curve vertically to match the observed integrated millimeter flux. The visibilities are particularly sensitive to the inner radius r_in, which is the dominant factor in determining the location of the first null in the visibility function, i.e., the shortest baseline at which the visibilities pass through zero and become negative. However, the narrower the ring, the greater the influence of the outer radius r_out on the null location, which can also be further moderated by how concentrated the emission is toward the inner radius, described by the surface density power law p. There is a fairly strong degeneracy between r_out and p; however, we include both parameters because the surface density power law is of interest for comparisons with the shorter-wavelength data.

We perform the χ² minimization in multiple steps, beginning with coarsely sampled grids and increasing the grid resolution to refine our estimates and avoid settling on a local rather than a global minimum. The initial coarse grid allows a to vary between 0.01 and 100 μm, β between −1 and 1.5, M_disk between 10⁻¹ and 10⁻⁴ M_⊕, r_in and r_out between 1 and 500 AU, and p from −1 to 4. Eventually the grids are refined so that the smallest (linear) steps are 1 μm for a, 0.1 for β, 10 AU for r_in and r_out, 10⁻³ M_⊕ for M_disk, 0.1 for p, and 2×10⁻⁴ M_⊕ for M_belt. The uncertainties on each parameter are larger than the steps, although they are also frequently correlated with other parameters. Given the low signal-to-noise ratio of the data and the rapidly improving capabilities of millimeter interferometers, our intention is not to provide a detailed characterization of the disk and a robust statistical analysis of the errors, but rather to generate a representative model that can reproduce the gross features of the SED and visibilities (both SMA and CARMA, below) to constrain basic properties like the size, width, mass, and degree of axisymmetry of the dust ring.

The best-fit parameters are listed in Table 2. The debris belt is broad, extending from 50 to 170 AU. The negative surface density power law (p = −0.3) is unexpected, as circumstellar disk density is generally expected to decrease with distance from the star. The power-law index is not tightly constrained, however (estimated uncorrelated uncertainty of ±0.3), and with the relatively coarse spatial resolution, the SMA data could be smoothing over multiple belts of differing masses and radii. The masses of the two dust belts (0.41 M_⊕ for the outer belt and 1.6×10⁻⁴ M_⊕ for the inner belt) are fairly typical for multiple-belt debris systems. Figure 1 shows a comparison of the data and model in both the image and visibility domains, as well as for the SED. The images demonstrate that the model subtracts cleanly, leaving no difference at the >3σ level that would indicate an azimuthal asymmetry.

The best previous spatially resolved long-wavelength imaging of the HD 107146 disk is the 1.3 mm CARMA observation by Corder et al. (2009). As an a posteriori check on the robustness of our model, we predict the spatial distribution of flux at the CARMA frequency from our SED + SMA model, subtract the model in the visibility domain using the MIRIAD task uvmodel, and image the residuals. Figure 3 shows the results of this comparison: the SMA model prediction subtracts cleanly from the CARMA visibilities, leaving no residuals at the >3σ level. This also implies that an axisymmetric model is consistent with the millimeter-wavelength data obtained so far from the HD 107146 system.

Zoom In Zoom Out Reset image size

Because of the very low signal-to-noise ratio in the HR 8799 data, the full multiparameter treatment we apply to HD 107146 is not warranted. Instead, we (1) compare our data with the prediction of the SED-based model of Su et al. (2009), and (2) explore whether the data can discriminate between a narrow or broad dust ring around the star. Throughout the analysis, we assume that the disk is viewed face-on (i=0°). The inclination is likely closer to 20°, based on a range of diagnostics including rotational velocity vsin i, astrometry of the planet orbits, stability analyses, and the marginally resolved Spitzer observations of an extended halo (Lafrenière et al. 2009; Reidemeister et al. 2009; Su et al. 2009; Fabrycky & Murray-Clay 2010; Moro-Martín et al. 2010). However, the data are not sensitive enough to distinguish between the two inclinations, nor are they sensitive enough to favor a particular position angle in any statistically meaningful way. As described below, it is difficult to determine the size of the disk from the visibilities alone to within a factor of two. The difference in radii derived for an inclination of 20° and 0° is only ∼7%, far less than the uncertainty in the data. Therefore, in the absence of a measurable position angle it is simplest to assume a face-on geometry.

Su et al. (2009) model three components to the dust distribution in the HR 8799 system: a warm (T ∼150 K) inner dust belt, a cold (T ∼45 K) outer planetesimal belt, and an extended (R ∼300–1000 AU) halo of small grains. Since only the cold planetesimal belt contributes substantially to the flux at millimeter wavelengths, we compare our data to this component. The model has been revised to take into account the analysis of the system conducted by Moro-Martín et al. (2010). It includes a dust belt with constant surface density (but decreasing surface brightness) as a function of radius, with a sharp inner edge at a distance of 150 AU from the star and an outer edge at a radius of 300 AU. The location of the outer edge is consistent with constraints from single-dish observations: since the 870 μm flux is effectively identical in the 19'' APEX beam (11.2 ± 1.5 mJy; P. Kalas 2011, private communication) and the 14'' JCMT beam (10.3 ± 1.8 mJy Williams & Andrews 2006), the disk radius must be ≲300 AU. The total mass of the dust in the belt is 1.2 × 10⁻¹ M_⊕, with grains ranging in size from 8 μm to 1 mm.

The left panels of Figure 2 show an image-domain comparison between the SMA data and the SED-based model of Su et al. (2009), calculated at the same wavelength as the data. The inset in the central panel shows the predicted spatial distribution of emission for the Su et al. (2009) model observed with the same noise level as the SMA data; the 4σ peak offset by 6'' from the star position and the partial arc of emission match well the observed morphology. The residuals subtract fairly cleanly, with only a slight negative 3σ peak. While such a large residual could potentially point to an azimuthal asymmetry in the dust distribution due to the gravitational influence of the planet HR 8799b, which is located at a similar azimuth to the negative peak, the signal-to-noise ratio is too low to permit a firm conclusion. The separation between the star and emission peak in the central panel inset is comparable to that in the data, certainly within the astrometric uncertainty of the observation. The right panel of Figure 2 makes the same comparison in the visibility domain, where the visibilities have been azimuthally averaged to enhance the signal-to-noise ratio in each bin. Within the uncertainties, the Su et al. (2009) model predicts well the spatially resolved data.

We also perform an exploration of the SMA visibilities and SED using a modified version of the method described in Section 4.1 aimed at determining whether the data can distinguish between a narrow or broad dust ring. To do so, we fix a to 1 μm, β to 0.5, and p to 0, implying a ring of constant surface density. We explore the narrow ring scenario by assuming that the width of the ring is 1 AU and allowing only r_in to vary. The SED is assembled from the literature, including 2MASS, IRAS, and Hipparcos fluxes (Moshir 1989; Moór et al. 2006; Williams & Andrews 2006; Su et al. 2009). To reproduce the SED, we allow M_disk to vary, and add a 150 K belt of variable mass M_belt, as in Su et al. (2009), to reproduce the observed mid-IR flux. We use a Kurucz–Lejune model photosphere (Lejeune et al. 1997) for a star of luminosity 4.92  L_☉, radius 1.4 R_☉, effective temperature 7250 K, and metallicity [Fe/H] = −0.5 (Sadakane 2006; Marois et al. 2008). The χ² minimization results in a best-fit model with a ring radius of 170 AU, with M_disk and M_belt equal to 5.5 × 10⁻⁴ and 4.0 × 10⁻⁶ M_⊕, respectively.

The combination of SED and visibilities can help to distinguish between a narrow ring and a broad belt. The null in the visibility function of the SMA data occurs at roughly 12 kλ; using Equation (A11) from Hughes et al. (2007), if the underlying flux distribution were a thin ring, it would have a radius of about 250 AU. On the other hand, if it were a broad ring like the Su et al. (2009) model predicts, then using Equation A9 of Hughes et al. (2007) the 12 kλ null implies a band with a smaller inner radius by up to a factor of three, depending on how quickly the flux decreases with distance from the star. The 170 AU result for the ring obtained in our SED+visibility analysis is odd in this context, since it is roughly two-third the radius implied by the visibility null and the location of the emission peak. The mismatch is due to the SED fit: the ∼45 K temperature implied by the SED (Su et al. 2009) is warmer than would be predicted for a 250 AU ring given the assumed dust properties. The 170 AU result also depends strongly on the choice of a 1 μm dust grain size. If the grain size is allowed to vary, even smaller grain sizes are preferred, since raising the temperature of the dust grains results in a better match between the SED and visibilities. It is therefore difficult to reconcile the large radius of a thin ring implied by the visibility data with the relatively warm dust implied by the SED for a narrow ring configuration, yet we have already shown that the observed flux distribution is also consistent with that predicted for a broad band of emission with an inner radius of ∼150 AU, as in the Su et al. (2009) model, which contains warmer dust than a narrow ring far from the star. This is consistent with the location of the visibility null for a broad ring configuration. Put another way, the spatially resolved SMA data alone are unable to discriminate between a narrow ring located ∼200 AU from the star and a broad disk with an inner radius of ∼150 AU, but the combination of visibilities and SED favors the latter scenario. Spatially resolved data with a higher signal-to-noise ratio are required to place more stringent constraints on the dust morphology, although the spatially extended Su et al. (2009) model is capable of reproducing both the detailed SED properties and the SMA data.

The location of the inner edge of the dust belt has profound implications for our understanding of the directly imaged planetary system orbiting HR 8799. The visibilities alone constrain the inner radius of the dust belt to between 80 and 250 AU, while the combination of SED and visibilities favors an inner radius of approximately 150 AU. As discussed by Su et al. (2009), the 90 AU inner radius initially estimated from the SED is consistent with the estimated extent of the chaotic zone of HR 8799b. Assuming a semimajor axis equivalent to the projected separation of the planet from the star (68 AU) and the nominal planet mass of 5–11 M_J (Marois et al. 2008), the chaotic zone ranges from approximately 17 to 21 AU from the planet depending on the assumed stellar age. Assuming a nearly face-on circular orbit, this places the chaotic zone edge at a distance of ∼85 AU from the star. It is therefore conceivable that the inner edge of the debris belt is being shaped by the outermost planet if it does indeed fall near a radius of 90 AU. Simulations by Moro-Martín et al. (2010) make it clear that the planetary system parameters, for several of the dynamical configurations described in recent stability analyses (e.g., Fabrycky & Murray-Clay 2010; Goździewski & Migaszewski 2009), are consistent with an inner disk radius up to 150 AU in radius. This brings it into agreement with the larger radius implied by the combination of SMA and SED data. In fact, an exploration of the parameter space of stable planet masses and eccentricities reveals that a semimajor axis of 150 AU effectively serves as an upper limit to the radius at which the disk is plausibly truncated by the planetary chaotic zone for any planetary system configuration with long-term stability. Given the ∼20 M_J upper limit on the mass of the outermost planet (Fabrycky & Murray-Clay 2010), only configurations with the outermost planet in a marginally stable orbit with varying semimajor axis can extend the overlapping resonances of the chaotic zone out to such a large distance from the star. While a nonzero eccentricity of the outermost planet could also extend the radius of the dust belt, the dynamical stability criteria for the system generally allow only a planet at the lower end of the allowable mass range to achieve only a moderate eccentricity, ruling out the case of both a large mass and an eccentric orbit for HR 8799b. More sensitive observations that pinpoint the location of the inner edge of the dust belt more precisely will therefore provide important constraints on the long-term stability of the system. The similarity of the SED shape to that of Fomalhaut suggests that the inner edge of the disk is sharp (Su et al. 2009), which is also consistent with dynamical truncation by interactions with the planetary system. Knowing the geometry of the system would aid in determining whether or not the debris belt is truncated by HR 8799b, since the semimajor axis of both its orbit and the inner edge of the debris belt could be more precisely determined. If the belt is truncated by the outermost planet, this would also provide an independent dynamical estimate of the mass of the planet, which currently depends on the highly uncertain age of the system.

No companion has to date been observed in orbit around HD 107146. It seems clear that the inner edge of the debris ring is not sculpted by a stellar companion, as searches by Metchev & Hillenbrand (2009) and Apai et al. (2008) have ruled out the presence of objects more massive than 10–12 M_J at separations of ∼15–75 AU from the star. However, there are no such limits on the presence of a <12 M_J companion. Corder et al. (2009) discuss the features that a putative planet would possess if it were responsible for sculpting resonances in the debris disk; however, since the SMA observations and modeling presented here do not confirm the clumpy structure suggested by the CARMA data alone, the properties of a planet interacting with the debris would likely be somewhat different. Assuming that the inner edge of the debris belt corresponds with the outer edge of the planet's chaotic zone, a planet with mass between 1 and 10 M_J would most likely be located at a deprojected radial separation of ∼44 AU from the star, given the 50 AU inner radius derived from the SMA data.

Of course, there are physical mechanisms other than dynamical interactions with an embedded giant planet that are capable of generating the observed dust ring morphology. Despite the upper limit on the gas density in the disk, even after most of the original mass has dissipated remnant gas can aid in shaping broad rings and haloes in the dust (e.g., Takeuchi & Artymowicz 2001). The collisional cascade initiated by the formation of smaller planetesimals could also account for the set of broad bands implied by the observations of this system (e.g., Kenyon & Bromley 2004), and may be more consistent with the matching extent of the optical and millimeter emission. However, the timescale necessary to form large (∼1000 km) planetesimals at such a great distance from the star is only marginally consistent with the 160 Myr age of the system, and the width of the emission band is too large to arise from only one such ring of planet formation.