EDP Sciences
Highlight
Free Access
Issue
A&A
Volume 591, July 2016
Article Number A108
Number of page(s) 22
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201628196
Published online 22 June 2016

© ESO, 2016

1. Introduction

Debris disks are the leftovers of star and planetary formation processes (see Wyatt 2008; Krivov 2010; Matthews et al. 2014, for recent reviews). Departure from photospheric emission at infrared (IR) wavelengths was first discovered around Vega (Aumann et al. 1984) using the Infrared Astronomical Satellite (IRAS). This excess emission was originally thought to be the remnant of the cloud out of which Vega formed. Several decades later, we now know that the dusty “debris” responsible for the IR emission are arranged in a disk comparable to the Edgeworth-Kuiper belt in the solar system. The dust grains, with sizes between a few μm to a few millimeters, located at tens of au from the central star are heated by the stellar radiation and re-emit at mid-IR and mm wavelengths. Since the original discovery, and mostly thanks to space-based missions such as the Infrared Space Observatory, Spitzer, and Herschel, several hundred of main sequence stars are known to harbor debris disks (e.g., Eiroa et al. 2013; Chen et al. 2014).

Table 1

Log for the VLT/SPHERE and ALMA observations.

Recent decades have seen incredible progress in the field of disk observations (e.g., Augereau et al. 1999; Kalas et al. 2005; Buenzli et al. 2010; Lebreton et al. 2012; Millar-Blanchaer et al. 2015). Observations with ever improving spatial resolution have revealed asymmetric disks (e.g., Lagrange et al. 2016; Kalas et al. 2015) as well as complex, moving, small-scale structures (Boccaletti et al. 2015). Nonetheless, spatially resolved observations of debris disks remain relatively rare and even though there are theoretical works focusing on the dynamical evolution of debris disks (e.g., Dominik & Decin 2003; Kenyon & Bromley 2006), they still need to be confronted with the observations. The current paradigm is that the primordial gas-rich proto-planetary disks are thought to have a half-life time of about 2−3 Myr (Hernández et al. 2007). As a disk evolves Pluto-sized planetesimals can form (Johansen et al. 2015) along with giant planets which may accrete their mass from the gas reservoir (core accretion or gravitational instability scenarios). After a few Myr, the gaseous content is quickly removed from the disk by efficient processes such as photo-evaporation (e.g., Alexander et al. 2006; Owen et al. 2011). Only already formed planets, planetesimals and dust grains will thus remain while the disk enters its debris disk phase. After a short phase of runaway growth, terrestrial planets may form in the inner regions of the disk, with Pluto-sized bodies in the outer regions, via chaotic growth of these oligarchs, on a timescale of 10−100 Myr (Kenyon & Bromley 2006, 2008, 2010). The time evolution of the entire system then becomes more regular and less chaotic. The km-sized bodies, arranged in one or more planetesimal belt(s), evolve under their mutual gravitational influence. Through collisions, they continuously release small particles in a collisional cascade. Small dust grains, in turn, are removed from the system either by radiation pressure or Poynting-Robertson drag. Therefore, one can consider that after ~100 Myr, planetary formation has stopped. The system is left with one (or more) planetesimal belt(s) and, quite possibly, with planets of various masses. By observing systems in the range 10−100 Myr, one can therefore study the time evolution of debris disks. Spatially resolved images can provide constraints on the radial and azimuthal distribution of the dust, giving us insight about the dynamics at stake in these systems.

Here, we present Very Large Telescope (VLT) Spectro-Polarimetric High-contrast Exoplanet REsearch (SPHERE) and Atacama Large Millimeter/submillimeter Array (ALMA) observations of the debris disk around the solar type star HD 61005 (G8V), located at a distance of 35.4 ± 1.1 pc (van Leeuwen 2007). The age of the system is believed to be within  Myr old, based on membership of the Argus association (Desidera et al. 2011; De Silva et al. 2013; Elliott et al. 2014). The uncertainties for the age are mostly related to the dispersion in ages reported in the literature for both the Argus association and the IC 321 super cluster. In the last ten years, it has been spatially resolved on multiple occasions with several instruments; Hines et al. (2007, Hubble ; Maness et al. (2009, HST/ACS); Buenzli et al. (2010, VLT/NaCo); Ricarte et al. (2013, Submillimeter Array SMA); and Schneider et al. (2014, HST/STIS). The disk earned its nickname of The Moth because of the swept-back wings first revealed in the HST observations of Hines et al. (2007). The wings may originate from the interaction between the interstellar medium (ISM) and the disk itself; as the star moves through the local ISM, (small) dust grains are set on eccentric orbits, drifting away from the central star.

While Hines et al. (2007) and Maness et al. (2009) mostly studied the wings, the study presented in Buenzli et al. (2010) resolved the debris disk as a ring, thanks to the better angular resolution provided by the NaCo instrument. They found the disk to be almost edge-on (i = 84.3 ± 1°), to have a semi-major axis of 61.25 ± 0.85 au with an eccentricity of e = 0.045 ± 0.015, which translates into an offset of 2.75 ± 0.85 au of the star with respect to the disk. A brightness asymmetry between the two ansae is observed, which cannot be fully explained by the off-centering. No planets were detected by Buenzli et al. (2010), for observations that should, in principle, have detected planets with masses starting from a few Jupiter masses. The integrated luminosity of the disk is large at IR wavelengths (Ldisk/L ~ 3 × 10-3) and is best modeled by two spatially separated dust belts. The main dust belt is the one located at ~60 au presented in Buenzli et al. (2010), but a warm component is usually required to reproduce the spectral energy distribution (SED). Ricarte et al. (2013) argue that this additional belt is mandatory to match the SED at about 20 μm, however its location remains unconstrained. In this paper, we focus on constraining the properties of the debris disk which is resolved at unprecedented angular resolution at both near-IR and mm wavelengths. Studying the properties and origin of the swept-back wings is beyond the scope of this paper, since they are marginally detected with these observations.

2. Observations, data processing, and stellar parameters

Table 1 summarizes the VLT/SPHERE and ALMA observations presented in this study.

2.1. VLT/SPHERE IRDIS observations and data reduction

The star HD 61005 was observed with the VLT/SPHERE (Beuzit et al. 2008), within the guaranteed time consortium. The observations were obtained in different instrumental set-ups in February, March, and May 2015, using the dual-band imager (IRDIS, Dohlen et al. 2008; Vigan et al. 2010), the integral field spectrograph (IFS, Claudi et al. 2008), and the dual-polarization imager (IRDIS DPI, Langlois et al. 2014).

2.1.1. IRDIS dual-band observations

The February observations were conducted in the H2H3 dual band (centered on 1.59 and 1.67μm) for IRDIS and the YJ-band (0.95−1.35μm, at spectral resolution R ~ 54) for IFS. The March observations used the K1K2 filters (centered on 2.11 and 2.25μm) for IRDIS and the YH-band (0.95−1.65μm, R ~ 33) for IFS. All of these observations were performed using an apodized Lyot coronagraph, consisting of a focal mask with a diameter of 185 milli-arcsec (N_ALC_YJH_S) and a corresponding pupil mask. Coronagraphic observations were performed in pupil stabilized mode to use angular differential imaging (ADI) post-processing (Marois et al. 2006) to attenuate residual speckle noise. The observation strategy can be summarized as follows: 1) photometric calibration: imaging of star offset from coronagraph mask to obtain PSF for relative photometric calibration; 2) centering: imaging with star behind mask with four artificially induced satellite spots for centering; 3) science: coronagraphic sequence; 4) centering: same as point two; 5) photometric calibration: same as point one; 6) sky background observation using same DIT as coronagraphic sequence. Finally, true north and plate scale are determined using astrometric calibrators as part of the SPHERE GTO survey for each run (Maire et al. 2016).

Basic reduction of the IRDIS data (background subtraction, flat fielding, centering) was performed using the SPHERE Data Reduction Handling (DRH) pipeline (Pavlov et al. 2008, version 15.0). The output consists of cubes for each filter, re-centered onto a common origin using the satellite spot reference. The cubes are then corrected for the true north position determined from the astrometric calibrations and for distortion. We collapsed the two filters together to increase the signal-to-noise ratio (S/N), and from now on we will refer to the H2H3 and K1K2 datasets as H and K observations, respectively. The datacubes in both bands were then processed using a principal component analysis (PCA, using the implementation of the scikit-learn Python package, Pedregosa et al. 2011) approach from all frames within the data-cube.

2.1.2. IRDIS dual polarization observations

On May 1, 2015, the target was observed using IRDIS in DPI mode in H-band. The same coronagraph was used for these observations. IRDIS DPI splits the light into two perpendicular polarization directions imaged at the same time on the same detector. Full cycles of half-wave plate (HWP) positions (0, 22.5, 45, and 67.5°) were taken to construct the Stokes Q and U vectors. The strategy of the observations was to take as long exposures as possible without saturating the detector just outside the coronagraph to achieve the best possible inner working angle as well as the best S/N for the outer part of the disk. Two integration times (DIT = 64 s with a total integration time of 768 s and DIT = 16 s with a total integration time of 3008 s) were used. IRDIS suffers from a de-polarization effect at certain detector position angles, which depend on the parallactic angle at the time of observation. Because the parallactic angle changed rapidly during the observations, we updated the detector angle at regular intervals.

The DPI data were reduced using a custom pipeline that is different to the DRH pipeline, which closely follows the processes described in Avenhaus et al. (2014), using the double-difference method (see also Canovas et al. 2011) to construct Stokes Q and U vectors from the data. The pipeline has been adapted to suit the IRDIS instrument. The data were centered using the centering frames taken just before and after the science observations. A de-rotation was applied to bring all files to the same orientation (see section above). Furthermore, the frames were corrected to take account of the fact that the IRDIS pixel scale differs slightly (~ 0.6%) in the two principal detector directions. The files were then corrected for true north as determined for IRDIS.

Because scattered light in an optically thin (debris) disk is expected to be polarized perpendicular to the line between the star and the image point in question, we then construct local Stokes Q and U vectors, denoted Qφ and Uφ (see also Benisty et al. 2015). In case of single scattering of the stellar light, Qφ is expected to contain the disk signal, while Uφ is expected to contain no signal (with possible exceptions when the optical depth is large, Canovas et al. 2015) but noise on the same level as the Qφ image and can serve as a noise estimator. Qφ and Uφ can be calculated as: (1)where φ refers to the azimuth in polar coordinates, and Φ is the position angle of the location of interest (x, y), with respect to the stellar location (x0, y0) as: (2)where θ corrects for instrumental effects such as a small misalignment of the half-wave plate.

During the data reduction process, one HWP cycle equivalent to 256 s of data (DIT = 16 s) was taken out because the telescope had lost tracking for a short amount of time, rendering this data unusable. The result are two pairs of Qφ and Uφ images, one for the DIT = 16 s and DIT = 64 s observations each. These were then combined with a weighted average to produce the final Qφ and Uφ images.

thumbnail Fig. 1

Reduced SPHERE observations of HD 61005 used in the analysis. North is up, East is left. From left to right; IRDIS ADI in H and K-bands (PCA with 6 components), and IRDIS DPI Qφ in H-band. Top row shows the data in linear stretch, with a central mask of 0.15′′, and the bottom row shows estimated signal-to-noise maps (see text for detail), with a stretch between [−3σ, 3σ].

Open with DEXTER

thumbnail Fig. 2

From left to right: observations, best-fit model and residuals of the ALMA data. For all panels the color map is a linear stretch between −0.27 and 1.23 mJy/beam. The standard deviation estimated in an empty region of the observations is of 0.09 mJy/beam. For both the observations and the model the contours are set at [3, 5, 7.5, 10]σ, and [σ, σ] for the residuals (no residuals beyond the 2σ level). The beam size is shown in the lower left corner of each panels.

Open with DEXTER

2.1.3. IFS observations

The IFS data proved difficult to be properly reduced at the time of this analysis, mainly because of centering problems. Presenting and analyzing these observations will be postponed for a future study.

2.1.4. Processed images

Figure 1 shows the final reduction for our dataset (with a central mask of radius 0.15′′). The top row displays the IRDIS ADI data in H and K-bands (left and middle panels, respectively), and the IRDIS DPI Qφ image (right), all with a linear stretch. For each image, the bottom row shows a S/N map, with a linear stretch between [−3σ, 3σ]. The noise map is calculated from the reduced images and represents the standard deviation in concentric annuli, centered on the star, with a constant width (2 pixels). We did not mask out the disk when computing the noise maps, hence the uncertainties might be slightly over-estimated. For the DPI observations, the noise map is computed from the Uφ image which does not seem to contain any signal from the disk (Fig. C.1 shows the Uφ image with the same linear stretch as the Qφ image in Fig.1).

We note that the disk is not homogeneously detected at high S/N. The east side is detected at larger S/N in all datasets and it appears brighter than the west side. Such asymmetry, already reported in Buenzli et al. (2010), are further investigated in this paper. For the sake of simplicity, in the rest of the paper, we refer to the ADI and DPI datasets as “scattered” and “polarized” observations, respectively. Strictly speaking, this is not correct as polarized photons must have been scattered by dust grains.

2.2. ALMA observations

HD 61005 was observed with ALMA in Band 6 (PI: David Rodriguez, program 2012.1.00437.S), in the frequency range 211.97−230.99 GHz. The target was observed several times, with precipitable water vapor ranging from 5.19 to 0.79 mm. We only kept the observations performed on the 20 of March 2014 for which the water vapor was minimum. Out of the four spectral windows, three were used to derive the continuum emission (211.97−228.97 GHz, with a 31 250 kHz resolution), while the last one was used to search for CO gas emission (230.05−230.99 GHz, with a 488.28 kHz resolution), which was not detected. Data processing was performed within CASA using the standard scripts provided by the observatory. We only kept the spectral windows used for the continuum observations and averaged the complex visibilities along the 128 different spectral channels (while flagging points with negative weights). Figure C.2 shows the (u, v) plane coverage for the continuum observations, with minimum and maximum baselines of 11.8 and 334.9 m, respectively. The reconstructed image (with so-called briggs weighting and pixel size of 0.13′′) is shown in the left panel of Fig. 2 (sensitivity of 0.09 mJy/beam). The beam size is 1.36′′ × 0.73′′ (48 au × 26 au) with a position angle of −86.5°. We do not attempt to measure the total flux of the disk directly from these observations, but we do estimate it when modeling the complex visibilities (Sect. 3).

2.3. Spectral energy distribution

Table 2

Broadband photometric measurements of HD 61005, and the equivalent widths of the far-IR filters (see text for details).

The star HD 61005 was observed by Herschel (Pilbratt et al. 2010) with the Photodetector Array Camera and Spectrometer instrument (PACS, Poglitsch et al. 2010), within the program OT2_tcurrie_1. The observation numbers (OBSID) are the two pairs 1342270977, 1342270978 and 1342270979, 1342270980 for the 70 μm and 100 μm observations, respectively. The 160 μm map used the four OBSID combined. The data were processed using the HIPE software (build 12.0.2083, Ott 2010), the very same way as described in Olofsson et al. (2013). HD 61005 was also observed with the Spectral and Photometric Imaging Receiver instrument (SPIRE, Griffin et al. 2010) in small scan map mode (OBSID: 1342255147 within the program OT2_kstape01_1). We used the Timeline Fitter task in HIPE to derive SPIRE photometry for our target. Calibration errors (~5.5%, Bendo et al. 2013) are included in the uncertainties. We also gathered photometric observations using VOSA1 (Bayo et al. 2008) and the dataset used to build the SED can be found in Table 2. The meaning of the third column is explained in Sect. 5.

Finally, we downloaded the Spitzer/IRS spectrum from the Cornell Atlas of Spitzer/IRS Sources database2 (Lebouteiller et al. 2011).

2.4. Stellar parameters

The stellar photospheric model is taken from the ATLAS9 Kurucz library (Castelli et al. 1997) with an effective temperature of T = 5500 K (Casagrande et al. 2011). With the dilution factor used to scale the photospheric model to the optical and near-IR photometric measurements, at a distance of 35.4 pc, we find a radius R = 0.84R. We derived a luminosity of L = 0.58L. To derive the stellar mass, which will become important when discussing the dust properties and the effect of radiation pressure on dust grains, we use isochrones from Siess et al. (2000), for an age of 40 Myr and effective temperature of 5500 K. We find that the stellar mass must be of about 1.1M (the corresponding luminosity matching our estimated L). We find a slightly smaller mass (1M) when using the isochrones from Baraffe et al. (2015), but the differences may arise from different model prescriptions (e.g. overshooting). In the following, we adopt a mass of 1.1M. The SED with the broadband photometric measurements, the Spitzer/IRS spectrum as well as the photospheric model are shown in Fig. 7 of Sect. 5.

2.5. Preamble on the modeling strategy

In this study, we aim to model observations from different facilities, at different wavelengths, using different techniques (interferometry and direct imaging). Therefore, before detailing the modeling strategy for each individual dataset, we provide a quick preamble on the methodology.

We first model the ALMA data (Sect. 3) assuming a circular disk, fitting the reference radius (r0, where the dust density peaks), the outer slope for the dust density distribution (αout), the position angle (φ), the inclination (i), and the total flux at 1.3 mm (f1300).

Prior to the modeling of the SPHERE observations, we attempt to constrain some of the dust properties, to limit the number of free parameters. We use the best fit results for the inclination and position angle from the modeling of the ALMA data to derive the polarized intensity as a function of the azimuthal angle from the Qφ image. We constrain the minimum and maximum grain sizes (smin and smax, respectively) as well as the porosity fraction of the dust grains (Sect. 4.1). This enables us to reduce the pool of free parameters when modeling the SPHERE DPI observations. For the SPHERE ADI observations, we use the Henyey-Greenstein approximation for the phase function, which disconnects the modeling process from the aforementioned dust properties. Thanks to the great complementarity between the ADI and DPI observations, we model them simultaneously to best constrain the azimuthal and radial dust density distribution (Sect. 4).

Finally, in Sect. 5, we model the SED of HD 61005, using the results inferred from the modeling of the SPHERE data on the location of the disk to derive stringent constraints on the dust properties (minimum grain size, dust composition, and total dust mass).

thumbnail Fig. 3

One- and two-dimensional (diagonal and lower triangle, respectively) projections of the posterior probability distributions for the results of the modeling of the ALMA observations.

Open with DEXTER

3. The parent planetesimal belt: constraints from the ALMA observations

Roughly speaking, different wavelengths trace different grain sizes. Therefore, we chose to model the ALMA observations independently of the SPHERE ones, the latter probing the small dust grains in the debris disk while the millimeter observations most likely trace a population of larger grains that more closely follow the parent planetesimals’ belt.

3.1. Modeling strategy

The modeling of the ALMA observations is performed in the Fourier space, attempting to reproduce both the real and imaginary parts of the complex visibilities (averaged along the spectral dimension). In Appendix A, we explain how we generate synthetic images and the different notations are summarized in Table C.1. From a synthetic image at the wavelength of 1.3 mm, we first scale the total flux of the image to the free parameter f1300 (in mJy) before computing the Fourier transform of the image. We then interpolate the Fourier transform at the spatial frequencies of the observations. The goodness of fit is the sum of the weights (estimated in CASA3) times the squared difference between the observed and modeled complex visibilities (the weights are proportional to 1 /σ2). We consider the following free parameters: the inclination i, the position angle φ, the total flux of the disk at 1.3 mm f1300, the reference radius4r0, and the outer power-law slope for the dust distribution αout, which is parametrized as (3)where r is the distance from the star and n the number density. The inner power-law slope is set to αin = 5, preliminary tests indicating this parameter is poorly constrained by the observations. To find the most probable solution, we use an affine invariant ensemble sampler Monte-Carlo Markov Chain, implemented in the emcee package, using 200 walkers, a burn-in phase of 500 iterations and a total length of the chains of 2000 iterations after the burn-in phase. At the end of the run, we find that the mean acceptance fraction (the mean fraction of steps accepted for each walker within the chain) is of 0.48 (a good sign of convergence and stability, Gelman & Rubin 1992). The maximum auto-correlation time for all the parameters is of 60 steps, indicating that the chains should have stabilized by the end of the simulations.

3.2. Results

Table 3

Best fit results for the modeling of the ALMA observations.

The projected posterior probability distributions are displayed in Fig. 3, for the different free parameters (using the triangle Python package, Foreman-Mackey et al. 2014). To derive the best-fit values as well as the uncertainties, we smooth the distributions with a kernel density estimator (the width of the Gaussian kernel σkde are reported in Table 3), and the best-fit value is the peak position of the distribution. The confidence intervals (a1, a2) for the parameter a are estimated as follows: (4)where γ = 0.68 and p(a) is the smoothed posterior probability distribution (integral normalized to 1) for parameter a (e.g., Pinte et al. 2008). We note that not all the distributions reach zero on each side of their maximum (especially for αout and i) and, therefore, these uncertainties should be treated carefully. Table 3 summarizes our results for the modeling of the Band 6 observations, and Fig. 2 shows the observations, the best-fit model, and the residuals (from left to right). The synthetic image of the best fit model is processed through CASA (using the ft method with the same antenna configuration as the observations) and the image is reconstructed with the clean algorithm with the same parameters as for the observations. We note that there are some residuals on the east side that may suggest that the disk is brighter on one side, even at mm wavelength. However, these residuals are below 3σ, therefore we cannot conclude they are significant. The apparent brightness asymmetry could be due to the asymmetric (u, v) coverage of the observations.

Overall, we find that most of the parameters are well constrained, except for the outer power-law slope of the dust density distribution, for which we can safely exclude slopes shallower than αout = −4. This is explained by the beam size of the observations which is larger than the debris disk for steep values of αout. Otherwise, we find the reference radius of the disk to be r0 ~ 66 au, the position angle φ ~ 70.7°, the inclination i ~ 84.5°, and the flux at 1.3 mm f1300 ~ 4.6 mJy. These results agree well with the parameters reported in Buenzli et al. (2010, i = 84.3 ± 1°, φ = 70.3 ± 1°, r0 = 61.25 ± 0.85 and Ricarte et al. (2013, r0 = 67 ± 2 au, φ = 71.5 ± 5°. The relatively large beam size of the ALMA observations can explain the slight discrepancy for r0 between the modeling of the ALMA data and the value inferred by Buenzli et al. (2010). This value will be revisited when modeling the SPHERE observations (Sect. 4). Ricarte et al. (2013) obtained a total flux of 7.2 ± 0.3 mJy at 1.3 mm with their SMA observations (Steele et al. 2016 obtained 8.0 ± 0.8 mJy analyzing the same observations), while we find the total flux to be well constrained at 4.6 ± 0.7 mJy (within the SMA and ALMA respective 3σ uncertainties). The shortest baselines being of Bmin ~ 11 and 16 m (for the ALMA and SMA observations, respectively), the largest scales the observations are sensitive to are of the order of 14.2′′ and 10.1′′, respectively (0.6λ/Bmin), much bigger than the disk. It is therefore unlikely that flux from the disk is filtered out by the interferometers. The differences between the SMA and ALMA data may arise from missing frequencies, differences in beam sizes, or calibration uncertainties (the uncertaintites reported by Steele et al. 2016 being more conservative than the ones reported by Ricarte et al. 2013).

Finally, we note that the famous wings responsible for the disk’s nickname are not detected in the ALMA data. Despite a good angular resolution, a disk model (without wings) can successfully reproduce the observations, and we see no trace of the wings in the residuals (with an rms of ~ 0.09 mJy/beam). Our results therefore agree with the ones of Ricarte et al. (2013); only small dust grains are likely present in the wings.

4. Constraining the dust radial distribution from the SPHERE observations

We model the SPHERE images by producing synthetic images at the central wavelength of the H-band observations, λc = 1.63μm. Given the low S/N of the K-band ADI observations, preliminary attempts to model these data showed that the dust distribution cannot be better constrained than with the H-band observations. We therefore focus the modeling effort on the H-band ADI and DPI data.

Figure 1 highlights the complementarity of both the scattered and polarized light images; the ADI and DPI data have very different S/N at the ansae and along the semi-minor axis of the disk. Combining both datasets, therefore, offers the opportunity to study the dust distribution in great detail.

The modeling process has a high dimensionality with many possible free parameters, and regions of intermediate to low S/N. Therefore, to obtain novel yet reliable constraints on the dust distribution we choose to perform a prior analysis on the observations to reduce the number of free parameters. Having a proper description of the polarized phase function prior to the modeling of the DPI observations greatly helps reducing the dimensionality of the modeling (e.g., the minimum and maximum grain sizes as well as the porosity of the dust grains). In this section, we first describe how the polarized phase function is derived and modeled, then we present the modeling strategy and summarize the results we obtain.

4.1. Phase function of the polarized light

To compute the polarized phase function from the DPI dataset, we define an elliptical mask with the following parameters: inner and outer radii (rin and rout), the inclination i, and the position angle φ. For each pixel within the elliptical mask, we compute the scattering angle as the dot product of the line of sight and the location of the pixel with respect to the star. We divide the elliptical mask in two, for the east and west sides and, for each side, we compute the minimum and maximum scattering angles. We then divide each side of the mask into 30 smaller regions corresponding to different bins of the phase function (see Fig. C.3 for an illustration). Each pixel is multiplied by its squared distance to the star, to account for illumination effect. The measured phase function is found by averaging the flux in the observations, in each individual region. Since the observed uncertainties σi are not the same for each pixel within a given intersection, the “average” uncertainty in the intersection is computed as follows: (5)where N is the number of pixels in the considered region. The phase function is then normalized to its maximum value (east and west sides are divided by the same value).

For the elliptical mask, we use the best-fit results of the modeling of the ALMA observations (see Table 3) for i and φ and choose rin and rout small and large enough (50 and 72 au, respectively) so that they encompass the trace of the debris disk. Figure 4 shows the phase function for the DPI H-band observations, the east side being in black, and the west side in red. The uncertainties are shown in a shaded color, and the gray area indicates regions of low S/N, which were estimated from the S/N map derived from the Uφ data. As can be seen in Fig. 1 (bottom right panel), the disk is only detected for a narrow range of azimuthal angles on the west side. This strongly suggests that we underestimated the uncertainties calculated using Eq. (5).

thumbnail Fig. 4

Azimuthal dependency of the polarized intensity in the DPI H-band observations, in black and red for the East and West sides, respectively. The shaded curves are representative of the uncertainties estimated from the noise map and the shaded gray area denotes low S/N regions estimated from the S/N map in Fig. 1. The thin solid lines display examples for different grain sizes for the same dust composition as the best fit.

Open with DEXTER

The difference in brightness between the two sides of the disk is striking in Fig. 4 and the east side appears almost twice as bright as the west side for most of the scattering angles. Even though the S/N for the west side is relatively low it seems that the polarized phase function peaks at an angle compatible with the east side.

To alleviate the number of free parameters, we aim to constrain the grain size distribution for the modeling of the DPI observations directly from the phase function displayed in Fig. 4. We emphasise that we do not have access to the polarization degree () as we do not have an unbiased measure of the total intensity I prior to the modeling. Indeed, the ADI process introduces self-subtraction effects (e.g., Milli et al. 2012), which can eventually be quantified with a model that describes the observations well (which we do not yet have). Since the west side suffers from low S/N, we perform the modeling on the east side’s phase function. To reproduce the phase function, we assume that the signal in the Qφ image is proportional to the size dependent S12(s) element of the Müller matrix. The matrix enables us to compute the Stokes vectors I and Q for the scattered and polarized light, respectively. Assuming single scattering event (reasonable assumption in low density environment such as debris disks), the scattered light will be the product of the first diagonal element of the matrix S11 times the stellar intensity I0. The polarized intensity will be proportional to the second element of the first column of the matrix S12 times I0. We compute S12 for different grain sizes between smin and smax, using the Mie theory, and we average S12 over a grain size distribution with a slope p< 0 (dn(s) ∝ spds) as follows, (6)We then compute the best scaling factor fS12 (to be multiplied to ) that will minimize the difference between the profiles and (with uncertainties σ) (7)We perform a simple grid search over the following parameters: smin, smax, and the optical properties of the dust grains. We fix p = −3.5, as expected for a collisional cascade in a debris disk (Dohnanyi 1969). For the dust composition, we consider a base medium of amorphous silicate with olivine stoichiometry (MgFeSiO4, Dorschner et al. 1995) to which we can add some porosity using the Bruggeman mixing rule. The minimum grain size can vary between 0.01 and 10μm. To constrain the maximum grain size, we vary the quantity Δs (= smaxsmin) between 0.01 and 100μm (in log space), and finally the porosity can change by steps of 10%. We find that the polarized phase function at 1.63μm is best reproduced by small spherical dust grains, with typical sizes in the range 0.3 ≤ s ≤ 0.35μm and a porosity fraction of ~ 80%. The best-fit solution is shown with a thick cyan line in Fig. 4, and reproduces well both the overall shape and the peak position of the observed phase function. Also shown in Fig. 4 are two examples for different grain size distributions around 0.25−0.3μm and 0.35−0.4μm, to illustrate how sensitive the phase function is with respect to the grain sizes. We note that, within such a narrow range of sizes, the value of the slope p of the grain size distribution remains unconstrained. In this section, we assumed that the polarized flux is directly proportional to S12, while it is also related to the dust density in the disk. In the rest of this paper, we try to determine the azimuthal distribution of the dust in the disk. Therefore, this prior analysis must be regarded as a first order approximation. The motivation of this prior analysis was to reduce the dimensionality of the modeling, but it also comes at a slight cost in the interpretation of the grain size distribution. The main conclusion of this section is that it does not seem that we observed large dust grains (which would have a stronger S12 signal at small scattering angles). In Sect. 6.1 we discuss this result further, but it could also be the consequence of low S/N along the semi-minor axis of the disk.

4.2. Modeling strategy

Since both the ADI and DPI datasets were taken at the same wavelength, for a given set of parameters we compute two images: an unpolarized light image using the Henyey-Greenstein (HG) analytical prescription of phase function S11 (Henyey & Greenstein 1941) and a polarized light image using the Mie theory. Using the HG approximation, which is parametrized by the anisotropic scattering factor g (−1 ≤ g ≤ 1), gives us more control when trying to reproduce the observed phase function (hence less free parameters to be considered during the modeling). To reproduce the DPI observations, we use the grain properties derived previously (Sect. 4.1). The absorption and scattering efficiencies, as well as the Müller matrix element S12, are computed with the Mie theory, which is valid for compact spherical grains. This approach may not appear self-consistent, but it is a way to disentangle the modeling of unpolarized and polarized light that may not be well accounted for by spherical grains (e.g., Milli et al. 2015).

The pool of free parameters includes the reference radius r0, the inclination i, the position angle φ, the opening angle ψ, and the outer power-law slope αout for the dust density distribution. Preliminary tests indicated that we can hardly constrain the inner power-law slope, so we fixed αin = 5. This enables us to focus on other parameters that describe the geometry of the debris disk. For instance, Buenzli et al. (2010) conclude that the eccentricity of the disk was not enough to explain the brightness asymmetry between the east and the west sides, and we aim to address this interpretation of previous observations. Therefore, we include the eccentricity e, the rotation angle Φe, the density damping η, and its azimuthal shape (via the width w of the Gaussian profile) and its reference angle Φη (see Appendix A). The azimuthal profile has the shape of a Gaussian profile with a σ = w, a peak of 1 for the azimuthal angle Φη and a minimum value of η ≥ 0 (hence an amplitude of 1−η).

Because we now consider azimuthal variations for the dust density distribution, it is parametrized slightly differently. The semi-major axis of the disk is defined as r0/ (1−e2). Although the formal denomination of r0 is the so-called semi latus rectum, we will simply refer to it as the reference radius in the rest of the analysis. The radius at which the dust density peaks now depends on the azimuthal angle θ. This radius rθ is defined as r0/ (1 + e × cosθ), and the azimuth-dependent dust density distribution n(θ) is parametrized similarly to Eq. (3), replacing r0 with rθ for each azimuthal angle.

The dust mass Mdust is not varied, but the synthetic images are scaled during the fitting process. For the polarized images, the modeled image is the absolute value of the Stokes Qφ parameter (we are assuming that there are no multiple scattering events in the debris disk). We scale the modeled image by a factor fDPI, which is found analytically to minimize the residuals (similarly to Eq. (7)). For the ADI dataset, this approach is not possible because the post-processing of the data cube can introduce self-subtraction effects. We therefore opted for a forward modeling strategy, similar to the one described in Thalmann et al. (2014). For a given set of parameters, we compute one synthetic image at the wavelength 1.63μm. We then produce a cube of 64 images, each one rotated to match the parallactic angles of each frame of the observations (total rotation of 93.3°). Each frame of the cube of synthetic images is multiplied by fADI and is then subtracted from the corresponding observed frame. The PCA process is performed, keeping only the six main components. Another approach would be to perform the PCA on the observations, save the coefficients, and apply them to the modeled cube. But the observations contain signal from the disk that may be accounted for in some of the principal components (even though they should be similarly subtracted in the modeled cube). Therefore, to ensure as little subtraction as possible, we chose to perform the PCA on the model-subtracted cube. The end goal being to minimize the flux in the final image.

The goodness of the fit is the sum of the squared ratio between the final images (residuals for the scattered and polarized data) and the noise map (bottom row of Fig. 1). For the modeling of both the ADI and DPI observations, we therefore have a total of 12 free parameters: r0, i, αout, e, Φe, φ, η, Φη, w, ψ, g, and fADI. Here, we also use the emcee Python package, with 100 so-called walkers in the chain, and first burn in 500 runs for each of these walkers. We then run the chain for 2000 iterations in total (see Millar-Blanchaer et al. 2015 for a similar modeling approach). To speed up the process, we cropped and re-sampled the SPHERE images. The size of one pixel in the new image (or cube) is resampled to be 2.25 times bigger than in the original images, the size of the new image being 200 × 200 pixels, and we use a central mask of radius 0.15′′. We chose not to convolve the synthetic images by a PSF because the cropped and down-sampled images have a pixel size of 0.028′′, while the approximation of the instrumental PSF with a Gaussian profile would have a width of 0.024′′ (approximating the PSF as an Airy disk for an 8 m telescope). At the end of the run, we find that the mean acceptance fraction (the mean fraction of steps accepted for each walker within the chain) is of 0.35, with a maximum auto-correlation time of 79 steps.

4.3. Results

thumbnail Fig. 5

Left to right: observed, residuals, and best-fit models for the ADI and DPI datasets (top and bottom, respectively).

Open with DEXTER

Table 4

Best fit results for the ADI and DPI H-band observations.

The projected posterior probability distributions are displayed in Fig. B.1, for the different free parameters. To derive the best-fit values as well as the uncertainties, we proceed similarly as in Sect. 3.2. The results are presented in Table 4 and Fig. 5 shows the observations, the residuals and the best-fit models (from left to right) for the ADI and the DPI data (top and bottom rows, respectively).

All parameters seem to be well constrained. The most probable solution has a semi-major axis of 60.9 au (r0/ [1−e2]) for an eccentricity of 0.093. With a rotation angle of Φe ~ 127° the pericenter is located slightly toward the observer, on the east side5. An azimuthal density variation seems to be necessary to reproduce the observations with a damping factor η of ~ 0.47 with a reference angle of ~ 138° (hence almost co-located with the pericenter of the eccentric disk). The azimuthal variation of the dust density distribution is a Gaussian profile with a width of ~ 50°. The position angle and inclination are consistent with the results from Buenzli et al. (2010). The aspect ratio of ψ ~ h/r ~ 0.06 agrees well with numerical simulations of the vertical structures of debris disks (e.g., Thebault 2009). Finally, the phase function is anisotropic for low scattering angles with g ~ 0.54. The dust density distribution for the best-fit model, viewed from above the disk, is displayed in Fig. B.2.

Nonetheless, as shown in the residuals of the ADI observations, our best-fit model does not manage to remove all the signal along the semi-minor axis of the disk. It successfully suppresses most of the signal at larger scattering angles, but the best-fit model seems to fail at properly describing the scattering at smaller angles. We tried to implement a phase function with two weighted HG functions, but could not significantly improve the residuals. We included all the “basic” parameters related to the geometry of the disk (i, φ, r0, ψ) and yet failed to perfectly match the observations. Possible explanations can be related to the phase function (the HG phase function remains an approximation), the radial segregation of the grain size distribution, or the azimuthal dust density distribution (this will be further discussed in Sect. 6.1). We crudely assumed a Gaussian profile for the azimuthal distribution, but the actual distribution could be skewed in one or another direction which could explain the brighter region along the semi-minor axis. Another (highly speculative but interesting) explanation could be that we may be seeing an inner disk. Because of the high inclination of the disk, an inner dust belt may appear to the observer as if it was merging with the main belt. The main challenge with this explanation is that we see no indication of this type of a belt in the DPI observations, but the noise is larger in the innermost regions of this dataset. One possible way to address this point would be to detect gas, which could trace velocities compatible with a radius smaller than 61 au. Yet no CO was detected in the ALMA observations (Sect. 6.2).

5. Constraining the dust mineralogy from the SED

The SED of the debris disk around HD 61005 is constructed from the fluxes reported in Table 2 and the Spitzer/IRS spectrum. Modeling an SED from unresolved observations is a degenerate problem. We are basically trying to find the adequate temperature of the grains, which can be changed either by their radial distances, their sizes, or their nature. The modeling of the SPHERE images provide strong constraints on the geometric parameters of the disk, and we can therefore focus on better constraining the dust properties.

5.1. Modeling the thermal emission from the disk

To model the SED of HD 61005, for a given set of parameters, the goodness of fit is computed as (8)where Fobs is the observed flux with its associated uncertainty σi, F the stellar contribution, and Fmodel the modeled thermal emission of the debris disk. The ωi values are weights to the observed data points at wavelength i. This weighting is designed to account for the fact that we simultaneously model broadband photometric observations (e.g., Herschel/PACS) and spectroscopic observations (Spitzer/IRS): two adjacent points in the IRS spectrum do not have the same significance as PACS 70 and 100 μm observations. We therefore opt for a similar strategy to the one described in Ballering et al. (2013). We compute the average spectral resolution of the IRS data and compare it to the equivalent widths of the broadband filters. One IRS point will have a weight ω = 1 and broadband photometric point will have a weight equal to the number of IRS points that would fit in the corresponding equivalent width of the filter. In Table 2, we report the equivalent widths for the various instruments used to model the thermal emission in the far-IR. These values were taken from the Spanish Virtual Observatory filter profile service6. For the ALMA Band 6 observations, we assumed the equivalent width to be equal to the spectral range (105μm).

To account for different mineralogical components, we use amorphous silicate grains of olivine stoichiometry (MgFeSiO4, Dorschner et al. 1995, ρ = 3.5 g cm-3) and amorphous water ices (Li & Greenberg 1998, ρ = 1.2 g cm-3). To mimic porosity, we also add the fraction of vacuum to the pool of free parameters. We mix the optical constants of the different dust components, using the Bruggeman mixing theory, and the Mie theory to compute the absorption coefficients Qabs. Since computing Qabs for large grain sizes can be the bottle-neck when computing thousands of models, we first created a library of opacity files (smin = 0.01 μm and smax = 5 mm, with steps of 10% for the water ices, and porosity fractions). During the modeling process, the appropriate Qabs are drawn and interpolated for the proper grain size distribution and composition.

The free parameters include the grain size distribution (smin and p, we fix smax = 5 mm), the relative fractions for the fractions of amorphous water ices and porosity, and the inner and outer slopes for the dust density (αin and αout, respectively). We fix the reference radius r0 to the best-fit results of the modeling of the SPHERE observations. For each set of parameters the dust mass Mdust is found by scaling the thermal emission of the model to the observed SED (similarly to Eq. (7)). The values for the dust mass are saved for posterior estimation of their corresponding uncertainties. To speed up the modeling of the SED, we neglect the non-azimuthal dust density distribution derived from the analysis of the SPHERE data and assume the disk to be a circular ring. We use the emcee package to search for the most probable model, using 100 walkers, burning the first 500 runs for each walker, and then running the chain for 1000 more iterations. The acceptance fraction at the end of the run is of 0.28 (and a maximum auto-correlation time across all parameters of 72 steps).

5.2. Results

thumbnail Fig. 6

Posterior probability distributions for the modeling of the SED. The solid lines indicate the most probable value, and the vertical dashed lines indicate the derived uncertainties.

Open with DEXTER

thumbnail Fig. 7

Spectral energy distribution of HD 61005. Photometric measurements are shown in red, along with the Spitzer/IRS spectrum (uncertainties, shown in gray, are 3σ, most of them being smaller than the symbol). The dashed gray line is the photospheric model, the black and cyan lines are the best fit models (including the stellar contribution, but with and without the inner component, see text for details).

Open with DEXTER

Table 5

Best-fit results for the modeling of the SED.

Table 5 summarizes the best fit results and the derived uncertainties and projected posterior probability distributions for each parameter are displayed in Fig. 6. The SED of HD 61005, along with the photometric and spectroscopic observations, the photospheric model, and the best-fit model are shown in Fig. 7. The shape of the mid- to far-IR excess is well reproduced by our model, but the turn-off point, where the excess emission starts, is not well matched (near λ ~ 20μm). This hints towards the presence of an additional component, an inner disk that has also been postulated by many authors in the literature (e.g., Morales et al. 2011; Ballering et al. 2013; Ricarte et al. 2013; Chen et al. 2014). Even when using detailed dust properties (instead of Planck functions at single temperatures) to model the excess, the best-fit model still underestimates the flux in the mid-IR beyond 3σ. To constrain the properties of an additional belt, we fit a Planck function to the residuals (the only free parameter being the temperature, the scaling to the residuals being done similarly to Eq. (7)). We find the temperature of the Planck function to best reproduce the residuals to be ~ 220 K. In Fig. 7, the Planck function is shown as a dotted-dashed black line and the final best-fit model (stellar contribution plus the inner and outer belts) as the solid black line. To assess the relevance of adding the inner disk we follow a similar approach as in Moór et al. (2015b, and references therein) and use the Akaike Information Criterion (AIC). We find that the addition of the inner component significantly improves the final model (the AIC including the inner belt is about 80% of the AIC without it), and therefore is deemed necessary to reproduce the whole SED of the disk around HD 61005.

Our modeling results suggest that the outer slope of the dust density distribution is relatively shallow, which would indicate that most of the grains responsible for the emission are located in an extended disk. The minimum grain size is found to be ~ 1.9μm, and the grain size distribution has a steep slope of p ~ −3.84. The composition of the grains would be a mixture of amorphous silicate grains with small fraction of water ices and porosity. The total dust mass (between smin and smax) is of the order of 0.37M.

6. Discussion

In the following, we try to put together the results of the modeling of the ALMA, SPHERE datasets as well as of the SED. We discuss what can be inferred regarding the dust properties, the gas content, the morphology of the disk to better characterize what is happening in this system. Because we modeled several datasets (ALMA, SPHERE, SED), Table 6 summarizes which parameters are constrained with which datasets, to clarify the discussion.

Table 6

Summary of the modeling strategy adopted in this study.

6.1. Dust properties

6.1.1. Self-subtraction corrected phase function

thumbnail Fig. 8

Phase function derived from the ADI H-band observations, corrected for self-subtraction effects (see text for details).

Open with DEXTER

With a satisfying model for the scattered light, we can further investigate the phase function in the SPHERE ADI H-band observations. Indeed, the model enables us to estimate the self-subtraction introduced by the PCA analysis. To do so, we first perform the PCA on the raw data cube and save the resulting coefficients. We then take the best-fit model, produce a data cube at the proper parallactic angles, project it into the new orthogonal basis found for the observations, and multiply this projection with the same coefficients as for the raw observations (e.g., Soummer et al. 2012). This provides us with an estimate of the signal from the disk that is subtracted when doing the PCA. The attenuation is estimated by the ratio rcorr between this estimate and the original model. When deriving the phase function for the ADI dataset (similarly to Sect. 4.1), each pixel is multiplied by the quantity 1/(1−rcorr). Figure 8 shows the phase function derived for the east and west sides (black and red, respectively). The bump at ~ 120° on the west side is an artifact introduced by rcorr (low signal in the disk model leads to an erroneous correction factor). In Fig. 8, we also show the phase function derived during the modeling (for g = 0.55), scaled to the phase function of the east for scattering angles larger than 40°. The inferred HG function appears to be a relatively good representation of the corrected phase function for scattering angles larger than ~ 50°, but severely fails to reproduce it at smaller angles. This would be an explanation for the strong residuals along the semi-minor axis in Fig. 5. Overall, this suggests quite extreme forward scattering in the disk around HD 61005 – a result that remains model-dependent since the PCA attenuation is estimated from our results (see also Milli et al. 2016, for a similar discussion for the disk around HR 4796 A). Nonetheless, the strong peak of scattering for small angles (down to ~ 10°) strongly points towards the presence of large dust grains in the debris disk. However, this seems in contradiction with the prior analysis of the DPI dataset where we had to actually reject the presence of large dust grains (s ~ 0.3μm). For the SED, we need a whole range of sizes (from μm- to mm-sized grains), but the emission is actually dominated by the smaller grains (discussed further in this section). Numerous studies of debris disks around different stars (e.g., Rodigas et al. 2015; Lebreton et al. 2012; Milli et al. 2015, 2016) already reported and discussed similar issues in reconciling the modeling results of different kind of observations. Throughout this paper, for their computational merits, we have either used the Mie theory or the HG prescription while both of them may not be accurate descriptions of the nature of the dust grains.

6.1.2. Reconciling different observations

To try to reconcile the different dust properties inferred from our modeling results, the radial extent of the disk is interesting to discuss. The ALMA observations point towards a narrow belt (αout ≤ −5), while the SED, ADI, and DPI observations suggest that the disk is extended (αout ~ −2 or −3). To assess the typical grain sizes we probe with the SED, we computed the integrated infrared luminosity for each grain size bins between smin and smax, and built the cumulative distribution function. We find that more than 98% of the disk infrared luminosity is accounted for by dust grains smaller than 5μm. The rest of the (very steep) grain size distribution contributes marginally to the total emission. This means that for the DPI observations and the SED, we are observing an extended disk containing small dust grains. We computed the β ratio between the radiation pressure and gravitational forces as in Burns et al. (1979), assuming L = 0.58L and M = 1.1M. For the dust composition found from the SED modeling, we find that grains larger than sblow ~ 0.8μm, for which β ≤ 0.5, should remain on bound orbits (assuming they were released on circular orbits). There is a small discrepancy between sblow and the minimum grain size we found (smin ~ 2μm) but as discussed in Pawellek et al. (2014), this could be caused by the assumption of spherical grains when modeling the SED. Assuming our inferred dust composition is not too far off compared to the real composition of the dust grains, it seems that radiation pressure should be efficient for (sub-) μm-sized grains. Consequently, the parent belt (inferred from the ALMA data) would be fairly narrow, but small grains placed on eccentric orbits (or even unbound) owing to radiation pressure make the disk look more extended at near-IR wavelengths.

The strong forward scattering peak derived from the corrected phase function from the ADI observations suggests the presence of large grains in an extended disk, which seems in contradiction with the aforementioned scenario. Nonetheless, according to Min et al. (2016), large dust aggregates behave like large grains (the size of the whole aggregate) in scattered light and like small grains (the size of the individual monomers) for polarized light. Such aggregates also display mild backward scattering that we do not detect with our observations but may be the consequence of low S/N along the back side of the disk.

Therefore, we hypothesize that the parent planetesimals are arranged in a narrow belt (traced by ALMA) and that the small-size end of the grain size distribution consists of dust aggregates that are pushed away by radiation pressure. This makes the disk appear extended in the SPHERE observations and when modeling the SED. For future studies of this system, it would be interesting to implement a radial segregation of the grain size distribution (varying dn(s) as a function of r, see for instance Stark et al. 2014 for a discussion of the phase function depending on the radial distance).

6.2. Gas mass upper limits

No CO emission was detected in the ALMA data, and we follow a similar procedure to that of Matrà et al. (2015) and Moór et al. (2015a) to derive upper limits for the CO gas mass. We extract the spectrum from the datacube reduced within CASA, centered around 230.5 GHz (12CO 21). The upper limit for the CO mass was determined as (9)where m is the mass of the CO molecule, d the distance of the star, ν2−1 the frequency of the transition, A2−1 the Einstein coefficient7, x2−1 the fractional population of the upper level, and S2−1 the observed integrated line flux. The fractional level was calculated assuming local thermal equilibrium, using the Boltzmann equation, with a temperature of 20 K (as discussed in Moór et al. 2015a, MCO is not strongly dependent on the temperature). S2−1 was computed from the spectrum as , where Δv is the channel velocity width, N the number of channels over a width of 15 km s-1 (centered at the local standard of rest velocity, vLSR), and Srms the standard deviation of the spectrum within that velocity range. From the heliocentric radial velocity of 22.5 km s-1 (Desidera et al. 2011), we obtain a vLSR of 3.68 km s-1. Assuming optically thin gas, this leads to an upper limit of 6.9 × 10-7M. As discussed in Matrà et al. (2015), the LTE hypothesis may not hold in low-density environment and our upper limit could underestimate the CO gas mass. In the ISM the CO/H2 abundance is 10-4 and with our inferred dust mass of ~ 0.37M (for our grain size distribution), this would lead to a gas-to-dust ratio that is much smaller than unity. However, it is unlikely that this ratio will remain the same in a circumstellar disk because of effects such as photo-dissociation. Overall, with the available observations leading to a non-detection, and given the several assumptions made (LTE, CO/H2 ratio) we can hardly conclude on the gas-to-dust ratio on the debris disk. This question should be addressed by future, deeper CO observations (or other species than the CO molecule).

6.3. Stirring of the planetesimals

In the last few years, theoretical works have aimed at characterizing the time evolution of debris disks (e.g., Kenyon & Bromley 2006, 2008). One of the purpose of these studies has been to better understand the growth of planetesimals to sizes of about a thousand km, in the framework of terrestrial planetary formation. It is believed that once these large bodies have been formed, they stir the population of smaller planetesimals. The stirring increases the relative velocities of these planetesimals and increases their chances of collisions, therefore initiating a collisional cascade. According to these models, collisions between planetesimals become destructive only once a Pluto-sized object is formed. Since the timescale for the growth of planetesimals scales with the distance r to the star, self-stirring is thought to be an inside-out process. Here, we follow a similar approach to that described in Moór et al. (2015b), who studied several debris disks spatially resolved with the Herschel observatory. They examined if the disks in their sample are consistent with a self-stirring scenario. To achieve this, they compared the timescale for forming a 1000 km-sized object in a self-stirred debris disk with the age of the stars. To quantify this timescale t1000 (in Myr), we used Eq. (41) from Kenyon & Bromley (2008), which we recall below: (10)where r is the distance of the dust belt, and xm is a scaling factor to parametrize the disk’s initial surface density (the higher, the more massive). Mustill & Wyatt (2009) argue that xm values larger than 10 are highly unlikely since the initial disk would have been gravitationally unstable. Assuming that the timescale for the growth of the planetesimals is the age of the system (i.e., t1000 = 40 Myr), we find that xm should be of the order of about 3.3 (for r = 61 au). This value would suggest that either the primordial disk was ~3 times more massive than the minimum-mass solar nebula, or that other sources of stirring (e.g., induced by a planet) need to be taken into consideration for this system. Debris disks with known planets, such as the ones studied in Moór et al. (2015b, e.g., HR 8799, β, have much larger xm values (≥ 10), while most of the debris disks usually have values for xm smaller than 3. It is therefore reasonable to assume that the disk around HD 61005 is bright because self-stirring by planetesimals has recently started. It does not imply that alternative stirring mechanisms should not be considered, but based on the distance of the dust belt and the age of the system, it is not mandatory to invoke the presence of a massive planet to explain what is currently observed.

6.4. An eccentric and asymmetric disk

thumbnail Fig. 9

From left to right: observations, models including some of the results of the SPHERE modeling, and residuals. Top row is for η = 0.47 and the bottom row for η = 0.6. For all panels the color map is a linear stretch between −0.27 and 1.23 mJy/beam. For both the observations and the model the contours are set at [3, 5, 7.5, 10]σ, and [−2, −1, 1, 2]σ for the residuals. The beam size is shown in the lower left corner of each panels.

Open with DEXTER

With the new SPHERE observations, we confirm the findings of Buenzli et al. (2010) that the disk displays a strong brightness asymmetry. Both the ADI and DPI observations show that the eastern side is brighter than the western side at an unprecedented angular resolution. The asymmetry is detected both in scattered light along the semi-minor axis (ADI dataset) and in scattered and polarized light along the semi-major axis (NACO observations of Buenzli et al. 2010 and DPI observations, respectively). When modeling the observations, we had to include a density damping factor along the azimuthal direction to reproduce the SPHERE observations. From preliminary tests, the DPI dataset (in which the brightness asymmetry is most visible) cannot be properly reproduced without these azimuthal variations. Therefore, we argue it can hardly be a scattering artifact and that it is, indeed, most likely related to a dust density enhancement in the disk.

On the other hand, the ALMA observations do not display such a strong asymmetry at first glance. To further investigate if the best fit model to the SPHERE observations is compatible with the ALMA data, we generated models including eccentricity and azimuthal variations. We took the best fit model to the ALMA data (Table 3), and the following additional parameters from Table 4: e = 0.093, Φe = 127°, w = 52°, Φη = 138°, and two values for η (0.47 and 0.6). Figure 9 shows the comparison between these models and the observations. For η = 0.47 (derived from the SPHERE observations), the brightness of the western side is slightly underestimated but the residuals barely reach the 2σ level. For η = 0.6, the model becomes more similar to the observations and the residuals do not reach the 2σ level. Consequently, we cannot exclude that the disk shows departure from centro-symmetry in the ALMA data, but the sensitivity of the observations does not allow us to put meaningful constraints on the azimuthal variations for the population of large grains. We can exclude damping factors of η ~ 0.5 to the 2σ level, but a value of 0.6 cannot be ruled out based on the available observations.

6.4.1. Possible planet-disk interactions?

Massive planets are often invoked to explain the eccentricity of a debris disk since the planet may shepherd it. The case of Fomalhaut is a good example of these types of studies. The disk shows an eccentricity comparable to the one of HD 61005 (e ~ 0.11 ± 0.1, Kalas et al. 2005) and a candidate companion was detected by Kalas et al. (2008). Studies such as the one of Beust et al. (2014) investigate in detail the possible interactions between the planet and the planetesimal in the disk. In the last two decades, many numerical works have been led to investigate how planets can shape debris disks (e.g., Wyatt et al. 1999; Wyatt 2006; Nesvold et al. 2013; Pearce & Wyatt 2014; Faramaz et al. 2014), or how the planetesimals evolve within their mutual gravitational interactions (e.g., Kral et al. 2015; Jackson et al. 2014).

Via secular interactions, an eccentric planet can shape the debris disk. Because of differential precession, the planet-disk interactions can result in a (short-lived) spiral arm (Mouillet et al. 1997; Wyatt 2005). Then, when taking into account the effect of collisions in the disk, the spiral quickly dissipates (Nesvold et al. 2013), and the disk becomes eccentric (and apse-aligned with the planet’s orbit). Because of the eccentricity, the disk is expected to be brighter at the pericenter (the so-called pericenter glow effect). Nonetheless, this may not be what is happening in the disk around HD 61005. First, we did not detect any point sources in neither the H- nor the K-band ADI datasets, even with their large parallactic rotations. The exact detection limits on the masses of possible planets will be published in the survey paper of the SPHERE consortium, but preliminary tests using the PynPoint package (Amara & Quanz 2012) shows that we should have been able to detect planets with a contrast of about ~10 mag in the near-IR. At an age of about 40 Myr, we can therefore rule out any planets with a mass larger than ~7 Jupiter mass at a distance of 0.2′′. Second, as already discussed in Buenzli et al. (2010, and confirmed in this study), the eccentricity of the disk is not sufficient for the pericenter glow to solely explain the observed brightness asymmetry. Our modeling results suggest that there is a ~50% reduction in density between the pericenter and the apocenter, with a peak density closer to the observer (Φe ~ Φη ~ 130°). Given that both Φe and Φη are free parameters and that their most probable values are so close to each other, it strongly points towards an increase of the dust density distribution at the pericenter of the eccentric disk. Interestingly, the simulations of Pearce & Wyatt (2014) show the opposite effect: the density is lower at the pericenter when including an eccentric planet in the disk (the particles spend more time at the apocenter than at the pericenter).

Wyatt (2006) describes the interactions between a massive planet with the dust belt, in different configurations. Whether the planetesimals are in 2:1 or 3:2 resonances with the planet, we can expect the disk to display asymmetries not only in scattered light but also at millimeter wavelengths. The large grains traced by mm observations are expected to closely follow the planetesimals that are in resonance with the planet. As assessed previously, with the current sensitivity of the mm observations, we cannot confirm nor rule out that large dust grains display departure from centro-symmetry.

Therefore, if a planet were to be detected (at other wavelengths or with other facilities), it would be interesting to confront our observations and modeling results with the orbital parameters of the planet. But in the following, we hypothesize that no massive planet was formed around HD 61005 and we try to explain the inferred morphology of the disk under this assumption.

6.4.2. The impact of gas in a debris disk

Lyra & Kuchner (2013) show that narrow, eccentric debris disks can be the consequence of the presence of gas. While debris disks are generally thought of as depleted of gas, several recent studies have demonstrated it is not necessarily the case, even at ages as old as a few tens of Myr (e.g., Hughes et al. 2008; Dent et al. 2014; Kospal et al. 2013; Moór et al. 2015a). Lyra & Kuchner (2013) suggest that the dust grains can heat up the surrounding gas via photoelectric heating. This will locally increase the gas pressure and, since the grains follow the pressure gradient, the pressure maximum will concentrate them even further. In the end, the disk may appear as narrow and eccentric, even though no planets are shaping it (although the original disk must be narrow). In Sect. 6.2, we concluded that the gas-to-dust ratio remains highly uncertain because of several unconstrained parameters (grain size distribution, CO/H2 abundance ratio, LTE assumption). Overall, this scenario may be a viable solution to explain our findings, but it cannot be verified at the moment (the azimuthal density variation would have to be addressed as well).

6.4.3. The aftermath of a collision between Pluto-sized bodies?

An alternative possibility would be that we are witnessing an impact between two massive bodies in the debris disk. This solution is mostly interesting because we find a peak of density at the pericenter of a slightly eccentric disk. In Sect. 6.3, we concluded that self-stirring by massive planetesimals may be taking place in the 40 Myr old debris disk. At this age, and a radius of 60 au, a few massive oligarchs may have been formed, which are disrupting the orbits of other planetesimals. The collisional cross-sections of these bodies would have increased along with the probabilities of intercepting each others’ orbits. Numerical studies of these types of impact (e.g., Jackson et al. 2014; and Kral et al. 2015) suggest that our results may be compatible with the scenario of an impact between two massive bodies. This type of collision would release significant amount of debris, and all of it would have to pass through the same location where the collision initially occurred. This would enhance the density and therefore the chances of subsequent collisions at the location of this pinch point. In which case, the scattering cross-section will locally increase along the eastern side compared to the western side. Similar scenarios have been postulated for the highly asymmetric disk around HD 15115 (Mazoyer et al. 2014), β Pictoris (Dent et al. 2014), and HD 181327 (Stark et al. 2014). According to the work of Jackson et al. (2014), the debris disk around HD 61005 could maintain this asymmetry for about ~ 0.45 Myr at 60 au (1000 orbital periods). However, in these numerical simulations, the eccentricity and overall complexity of the resulting disk highly depends on the initial conditions of the impact. Furthermore, Leinhardt & Stewart (2012) find that the grain size distribution of the fragments released by the collision of gravity-dominated bodies is best described by a power-law slope in −3.85. Our modeling results of the SED are in excellent agreement with this value, supporting the scenario that the disk mainly consists of fragments from a catastrophic collision.

Nonetheless, a challenging aspect to this postulate is the dust mass inferred from the SED modeling. Simulations of the stochastic phase of planetesimal growth by Stewart & Leinhardt (2012) suggest that collision of massive bodies should release about ~ 15% of the total mass in debris. Assuming “debris” corresponds to grains with sizes between 2μm and 5 mm (used to model the SED), this would correspond to a total mass of about 2.5M for the catastrophic impact, which would make this type of event an outliner in most theoretical works about the formation of planetesimals. For a 40 Myr old solar-mass star, the models of Kenyon & Bromley (2008, 2010) predict a dust mass production of the order of ~ 0.05M (for 1μm s ≤ 1 mm) for a disk that is initially 3 times more massive than the MMSN (to be compared to 0.25M interpolating our results to the same grain size distribution). If the primordial disk around HD 61005 was indeed several times more massive than the MMSN (see Sect. 6.3), it would evolve more rapidly and massive bodies would have had the time to form. Overall, if this interpretation of the multi-wavelength observations of the Moth holds, it would suggest that massive bodies can still be formed at a large distance from the star at an age of ~ 40 Myr. This would provide valuable constraints on the formation and evolution of planetesimals in the first 100 Myr.

The presence of spatially separated belts (inner and outer belts) is often discussed in the context of planet-disk interactions (e.g., Su et al. 2013). Massive planets may scatter away any close-by debris, hence opening a dust-free region in the disk. Planets may also be able to scatter debris from the outer regions and send them inwards, which would result in a bi-modal distribution of the dust (e.g., Bonsor & Wyatt 2012). In the case of HD 61005, we do not have strong constraints on the location of the inner belt. But, given that the evolution of a debris disk is an inside-out process (Kenyon & Bromley 2008), it is unlikely the inner belt is primordial. Given the age of the star, the inner regions should already have been depleted of any detectable amount of small dust grains. Still, assuming that there are no massive planets around HD 61005, is it possible to reconcile the presence of an analog of the asteroid belt with the above interpretation? If it is indeed the case that a few oligarchs have formed in the main belt at 60 au, the disk would be in an active and chaotic state. Some of these massive bodies may have collided in the main belt, but some may have scattered debris inwards. However, simulations, such as the ones presented in Jackson et al. (2014) and Kral et al. (2015), do not display cleared gaps. Consequently, under our assumption that there are no giant planets, it seems that our interpretation of the dynamical evolution of the outer regions is not directly related with what is taking place in the inner regions. It may well be that telluric planets have formed in the inner regions and further investigation remains necessary to better constrain the location of the innermost belt.

7. Conclusion

In this study, we presented new VLT/SPHERE and ALMA observations of the debris disk around the 40 Myr-old solar-type star HD 61005. We modeled both the radial and azimuthal distribution of the dust grains. We find that the disk is slightly eccentric (e ~ 0.1), with a dust density two times smaller at the apocenter compared to the pericenter. The sensitivity of the ALMA observations cannot confirm nor rule out an asymmetric disk for large mm-sized dust grains. We confirm that the eccentricity of the disk is not sufficient to explain the known brightness asymmetry in the near-IR. Thanks to the multi-wavelength analysis, we highlight a possible radial segregation of the grain size distribution since the disk appears more extended at near-IR wavelengths compared to far-IR and mm wavelengths. The parent planetesimals are confined within a narrow belt (with a semi-major axis of about 61 au), while small dust grains (or dust aggregates) would be pushed away by radiation pressure. No CO gas was detected with the ALMA observations and we did not detect any giant planets around the star. For this reason, we attempt to explain the morphology of the disk without invoking planet-disk interactions. We propose that we could be witnessing the aftermath of collision(s) between massive bodies in the main belt. The debris disk may be bright because self-stirring recently started to become significant in the belt, and it triggered a massive collision. The location of the collision became a pinch point where all debris must go through at each orbit, hence locally increasing the dust density and the scattering cross-section, as well as the chances of subsequent collisions. Nonetheless, this scenario would require a massive primordial disk to reconcile numerical simulations of debris disk evolution with the inferred dust mass. Finally, the presence of an inner belt would be consistent with the modeling of the SED, and we conclude its presence and dynamical evolution is unrelated to the evolution of the outer belt.


2

The Cornell Atlas of Spitzer/IRS Sources is a product of the Infrared Science Center at Cornell University, supported by NASA and JPL. http://cassis.sirtf.com/atlas/query.shtml

3

The absolute values of the weights derived within CASA may be inaccurate but their relative values are not.

4

Here we assume the disk is circular.

5

For Φe = 180°, the pericenter would be along the semi-major axis on the east side. Smaller angles would move the pericenter towards the observer, while larger angles would move the pericenter towards the back side of the disk.

7

Taken from the SPLATALOGUE catalog at http://www.cv.nrao.edu/php/splat/

Acknowledgments

We would like to thank the anonymous referee for the valuable feedback we received. The comments improved the discussion of the paper, mostly with respect to the scenario of a massive collision. They also helped with the presentation of the modeling strategy for all the different datasets presented in this paper. This paper makes use of the following ALMA data: ADS/JAO.ALMA 2012.1.00437.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSC and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. J. Olofsson, H. Avenhaus, C. Caceres, H. Canovas, M. R. Schreiber, A. Bayo, and S. Casassus acknowledge support from the Millennium Nucleus RC130007 (Chilean Ministry of Economy). J. Olofsson acknowledges support from ALMA/CONICYT Project 31130027, H. Avenhaus the support from FONDECYT 2015 Postdoctoral Grant 3150643, C. Caceres the support from CONICYT FONDECYT grant 3140592, H. Canovas the support from ALMA/CONICYT (grants 31100025 and 31130027), and A. Bayo the financial support from the Proyecto Fondecyt Iniciación 11140572. This work was supported by the Momentum grant of the MTA CSFK Lendület Disk Research Group. A. Moór acknowledges support from the Bolyai Research Fellowship of the Hungarian Academy of Sciences. J. Olofsson would like to thank Daniel Foreman-Mackey for his help with the parallelization of the code with emcee and Mark Booth for valuable discussions about ALMA data. This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration 2013). This research has made use of the SIMBAD database (operated at CDS, Strasbourg, France). The authors have used the TOPCAT software (Taylor 2005) in this work.

References

Appendix A: Description of the code

Appendix A.1: Computing the spectral energy distribution

thumbnail Fig. A.1

Schematic view of how a disk, seen edge-on, is described in the DDiT code. Different shells (nr = 10) are represented in different fading colors, while dust belts are represented without transparency (see text for details).

Open with DEXTER

The DDiT code is mostly based on the code already described in Olofsson et al. (2012, 2013, which was based on the “Debris disk radiative transfer simulation tool” described in Wolf & Hillenbrand 2005. However, as opposed to the work presented in Olofsson et al. (2012, 2013), only one dust species per model is considered in this study. The dust properties are defined by the optical constants, the density (ρ), the minimum and maximum grain sizes (smin and smax, respectively), and the slope p of the grain size distribution (dn(s) ∝ spds, p< 0). The absorption efficiencies, Qabs(λ,s), are computed using the Mie theory for ng = 100 grain sizes logarithmically spaced between smin and smax.

Concerning the disk model itself, this is defined by the following parameters: a reference radius r0, two power-law indexes αin> 0 and αout< 0, the opening angle of the disk (ψ = arctan(h/r), where h is the height from the midplane), and the total dust mass Mdisk. We first define the inner and outer radii of the disk by finding rin = r0 × C(1 /αin) and rout = r0 × C(1 /αout), where C is a threshold value for the minimum density that we choose to be C = 10-2. The code is 1D and the disk is divided into nr spherical shells centered at rj, logarithmically spaced between rin and rout (see Fig. A.1).

The dust number density distribution for the shell at distance rj and the grain size sk, Ndens(rj,sk), is defined as (A.1)where Ns(sk) is a function depending solely on the grain size (similarly to Eq. (2) of Dullemond & Dominik 2008) (A.2)where Δlog (s)k is the width of each bin in logarithmic space of s. Since the ng grain sizes are logarithmically spaced, each Δlog (s)k is the same, except the first and last ones which are half of that value. The mass density distribution is then given by . Both Ndens(rj,sk) and Mdens(rj,sk) are then normalized so that the integral of Mdens(rj,sk) × Vj (where Vj is the volume of the shell j) over s and r equals the total input dust mass Mdisk.

The radially integrated optical depth at the wavelength of the peak of the stellar emission (τλ peak) is computed at the outer edge rout, assuming a spherical symmetry of the disk. This quantity is then divided by fψ = sin(ψ) to obtain the final optical depth of the model, for a given ψ. This 1 /fψ factor accounts for the fact that the code is 1D, i.e., that the disk is not a sphere but has a maximum vertical height h from the midplane (hence can be seen as dust belts). In Fig. A.1, the shells are shown in fading colors, while the dust belts are represented without transparency. The fψ factor is the ratio of the volume of a sphere of radius R to which we subtract two spherical cones (up and down) of height R × (1−sin(ψ)), divided by the volume of a sphere of the same radius. One spherical cone has a volume of 2π/ 3R3(1−sin(ψ)), which simplifies the volume ratio to sin(ψ). By doing so, we can keep Mdisk constant and leave the output of the code unchanged while simply updating the optical depth. Once the model has been computed, a check is performed to see if τλ peak is smaller than unity for the given dust mass and ψ.

When computing the SED, we consider that the contribution of scattered light is negligible at mid to far-IR wavelengths compared to the thermal emission. Furthermore, since debris disks are optically thin at all wavelengths, geometrical parameters such as the inclination have no impact on the final SED. Finally, as the observed SED is constructed from spatially unresolved observations, disk or aperture sizes do not have to be taken into consideration.

For each shell, the temperature of a dust grain is computed as a function of its size s and distance r to the star. Assuming radiative equilibrium, one can express r as (A.3)where R is the stellar radius, F the stellar surface flux, and Bλ the Planck function at the dust temperature Tdust. By inverting numerically the above equation, we obtain the radial dependency of the temperature. Afterwards, the flux Fth (the thermal radiation of a single dust grain) received by an observer at distance d is computed as (A.4)For a given shell, the total thermal emission is obtained by multiplying Fth(λ,rj,sk) by the normalised Ndens(rj,sk). The final SED is the sum of the thermal emission from all nr shells.

Appendix A.2: Computing synthetic images

thumbnail Fig. A.2

Sketch of a disk (inclined by i) as seen from the side to illustrate the entry and exit points (DEn and DEx) when computing images. The cut shows the y and z in the disk frame axis and is made at x = 0 (the x axis pointing toward the reader).

Open with DEXTER

At optical or near-IR wavelengths, stellar light scattered by dust grains can no longer be considered negligible and has to be accounted for. We follow a parallel approach as for the thermal emission and the stellar light scattered by a single dust grain is computed as follows (A.5)where Qsca is the scattering efficiencies calculated with the Mie theory and S11(θ) is the scattering phase function.

In the model parametrization, we define two quantities for the image itself, the number of pixels np (images are squares) and the size lp of one pixel in units of au. Then, for each of the pixels in the image, we first apply adequate projection and rotation to account for the geometrical parameters i and φ (inclination and position angle, respectively). For φ = 0°, the semi-major axis is along the south-north direction, and the front side of the disk is chosen to be on the east side. For increasing φ the disk rotates counter-clockwise, towards the east.

Since the disk has a non-zero opening angle ψ, it is not infinitely flat in the vertical direction. As in other codes (e.g. GRaTeR, Augereau et al. 1999), we assume a vertical Gaussian profile for the dust distribution, with a σ width equal to ψ. This avoids having blocky images, especially for highly inclined disks. Figure A.2 shows a crude sketch of a cut in the (y,z) reference plane of the disk, for x = 0. The disk has an inclination i and for each pixel, we trace rays along the line of sight and search for intersections with the disk. The disk can be described by four main geometrical shapes: inner and outer shells, and upper and lower planes. Each of them can be either “entry” or “exit” points for the rays (blue and yellow lines in Fig. A.2, respectively). One should note that only slightly different intersection conditions are to be considered for cases where π/ 2−i<ψ (e.g., edge-on disk) and that the code can produce images for such disks. We can therefore obtain the coordinates of both the entry and exit points (DEn and DEx, respectively). Since a pixel has a size lp, the column between DEn and DEx (modulus of | l |) has a volume of . Since this column can probe a significant range of distances r to the star, we divide it into nlos boxes, along the line the sight. For each of these boxes, we evaluate the distance to the star r, the volume of the box Vbox, and the scattering angle (dot product between the line of sight and the position in the disk with respect to the star). We then interpolate Fth, Fscat, and Ndens at the distance r of the center of the box, and interpolate the phase function at the scattering angle. We then sum the product Fth/scat × Ndens × Vbox for each box, to obtain both the scattered light and thermal emission for this given pixel.

Appendix A.2.1: Polarized images

The code can also produce polarized light images on top of the total intensity images. Neglecting multiple scattering events, a reasonable assumption in optically thin disks, we compute Fpol(λ,rj,sk) similarly to Eq. (A.5) simply replacing S11 by S12. The second element of the Müller matrix S12 is calculated using the Mie theory (while S11 can be approximated by the Henyey-Greenstein function). We can therefore produce synthetic images along the Stokes vector Q, which can be used to model polarimetric observations.

thumbnail Fig. A.3

Sketch of an eccentric disk. The color-coding represents the density distribution. The rotation angle is indicated by Φe, the reference angle for the azimuthal variation by Φη, and η = 0.3 in this example.

Open with DEXTER

Appendix A.2.2: Eccentric disks

The code can produce images of eccentric disks. The inner and outer shells as described in Fig. A.2 are modified to inner and outer edges with an eccentricity e parameter and a rotation angle Φe (see Fig. A.3 for a sketch of a face-on disk). For each pixel in the image, at position (x,y) (in the reference frame rotated by Φe), we check if it is in between the inner and outer edges, which are defined as, (A.6)where a = atan2(y,x). For Φe = φ = 0°, the pericenter is located on the south side.

Appendix A.2.3: Azimuthal density variations

It is also possible to introduce azimuthal dust density variations in the synthetic images. The density distribution is parametrized with three parameters: a damping factor η (0 ≤ η ≤ 1), a reference angle Φη, and the azimuthal variation follows a Gaussian profile of σ = w that peaks at 1 for atan2(y,x) = Φη and has an amplitude of 1−η. The attenuation factor is of the form (A.7)For Φη = φ = 0°, the density maximum is located on the south side. For a given azimuthal angle, Ndens is multiplied by Natt before computing the image.

Appendix B: Results

thumbnail Fig. B.1

Projected posterior probability distributions for the combined modeling of the ADI and DPI H-band observations.

Open with DEXTER

thumbnail Fig. B.2

Top view of the dust density distribution (normalized to its maximum) for the best-fit model to the SPHERE observations. The direction to the observer is indicated by the arrow and the spatial scale in the lower right corner.

Open with DEXTER

Appendix C: Miscellaneous

thumbnail Fig. C.1

Uφ component of the H-band DPI observations.

Open with DEXTER

thumbnail Fig. C.2

Coverage of the (u, v) plane for the ALMA observations used in this paper.

Open with DEXTER

thumbnail Fig. C.3

Elliptical mask binned over different scattering angles to depict how the phase function is calculated. The second bin of the phase function would be the mean flux in the region highlighted in white. The color coding shows the scattering angle for a given position in the elliptical mask. The central circle shows the numerical mask to hide the central region.

Open with DEXTER

Table C.1

Notations used in this study.

All Tables

Table 1

Log for the VLT/SPHERE and ALMA observations.

Table 2

Broadband photometric measurements of HD 61005, and the equivalent widths of the far-IR filters (see text for details).

Table 3

Best fit results for the modeling of the ALMA observations.

Table 4

Best fit results for the ADI and DPI H-band observations.

Table 5

Best-fit results for the modeling of the SED.

Table 6

Summary of the modeling strategy adopted in this study.

Table C.1

Notations used in this study.

All Figures

thumbnail Fig. 1

Reduced SPHERE observations of HD 61005 used in the analysis. North is up, East is left. From left to right; IRDIS ADI in H and K-bands (PCA with 6 components), and IRDIS DPI Qφ in H-band. Top row shows the data in linear stretch, with a central mask of 0.15′′, and the bottom row shows estimated signal-to-noise maps (see text for detail), with a stretch between [−3σ, 3σ].

Open with DEXTER
In the text
thumbnail Fig. 2

From left to right: observations, best-fit model and residuals of the ALMA data. For all panels the color map is a linear stretch between −0.27 and 1.23 mJy/beam. The standard deviation estimated in an empty region of the observations is of 0.09 mJy/beam. For both the observations and the model the contours are set at [3, 5, 7.5, 10]σ, and [σ, σ] for the residuals (no residuals beyond the 2σ level). The beam size is shown in the lower left corner of each panels.

Open with DEXTER
In the text
thumbnail Fig. 3

One- and two-dimensional (diagonal and lower triangle, respectively) projections of the posterior probability distributions for the results of the modeling of the ALMA observations.

Open with DEXTER
In the text
thumbnail Fig. 4

Azimuthal dependency of the polarized intensity in the DPI H-band observations, in black and red for the East and West sides, respectively. The shaded curves are representative of the uncertainties estimated from the noise map and the shaded gray area denotes low S/N regions estimated from the S/N map in Fig. 1. The thin solid lines display examples for different grain sizes for the same dust composition as the best fit.

Open with DEXTER
In the text
thumbnail Fig. 5

Left to right: observed, residuals, and best-fit models for the ADI and DPI datasets (top and bottom, respectively).

Open with DEXTER
In the text
thumbnail Fig. 6

Posterior probability distributions for the modeling of the SED. The solid lines indicate the most probable value, and the vertical dashed lines indicate the derived uncertainties.

Open with DEXTER
In the text
thumbnail Fig. 7

Spectral energy distribution of HD 61005. Photometric measurements are shown in red, along with the Spitzer/IRS spectrum (uncertainties, shown in gray, are 3σ, most of them being smaller than the symbol). The dashed gray line is the photospheric model, the black and cyan lines are the best fit models (including the stellar contribution, but with and without the inner component, see text for details).

Open with DEXTER
In the text
thumbnail Fig. 8

Phase function derived from the ADI H-band observations, corrected for self-subtraction effects (see text for details).

Open with DEXTER
In the text
thumbnail Fig. 9

From left to right: observations, models including some of the results of the SPHERE modeling, and residuals. Top row is for η = 0.47 and the bottom row for η = 0.6. For all panels the color map is a linear stretch between −0.27 and 1.23 mJy/beam. For both the observations and the model the contours are set at [3, 5, 7.5, 10]σ, and [−2, −1, 1, 2]σ for the residuals. The beam size is shown in the lower left corner of each panels.

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

Schematic view of how a disk, seen edge-on, is described in the DDiT code. Different shells (nr = 10) are represented in different fading colors, while dust belts are represented without transparency (see text for details).

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

Sketch of a disk (inclined by i) as seen from the side to illustrate the entry and exit points (DEn and DEx) when computing images. The cut shows the y and z in the disk frame axis and is made at x = 0 (the x axis pointing toward the reader).

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

Sketch of an eccentric disk. The color-coding represents the density distribution. The rotation angle is indicated by Φe, the reference angle for the azimuthal variation by Φη, and η = 0.3 in this example.

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

Projected posterior probability distributions for the combined modeling of the ADI and DPI H-band observations.

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

Top view of the dust density distribution (normalized to its maximum) for the best-fit model to the SPHERE observations. The direction to the observer is indicated by the arrow and the spatial scale in the lower right corner.

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

Uφ component of the H-band DPI observations.

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

Coverage of the (u, v) plane for the ALMA observations used in this paper.

Open with DEXTER
In the text
thumbnail Fig. C.3

Elliptical mask binned over different scattering angles to depict how the phase function is calculated. The second bin of the phase function would be the mean flux in the region highlighted in white. The color coding shows the scattering angle for a given position in the elliptical mask. The central circle shows the numerical mask to hide the central region.

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.