Free Access
Issue
A&A
Volume 537, January 2012
Article Number A110
Number of page(s) 15
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201117731
Published online 17 January 2012

© ESO, 2012

1. Introduction

Over the past few decades, hundreds of debris disks have been discovered by various infrared facilities, e.g., IRAS, ISO, JCMT/SCUBA, Spitzer, AKARI, and now Herschel. These disks were all detected because of their dust components and the excess thermal emission they produce. Their observations have revealed a wide variety of inferred dust masses, disk extents, and even disk morphologies (e.g., Wyatt 2008, and references therein). For example, dust has been detected from well within one AU (Beichman et al. 2005; Absil et al. 2006), out to hundreds of AU (e.g., Su et al. 2005). While the former cases seem to be rare (Wyatt et al. 2007), the latter are more commonly reported, corresponding to significantly higher detection rates at longer wavelengths. In addition, from a purely observational point of view, a huge population of fairly massive, but cold disks could still be hidden below current sensitivity limits. Some of these will be revealed by Herschel (Pilbratt et al. 2010; Eiroa et al. 2011).

Most of the known debris disks have been characterised using spectral and photometric data, with only a few dozen having been spatially resolved. When it comes to modelling, these sparse data sets unavoidably introduce degeneracies between the properties of grains (particularly their typical sizes and distances to the star). A given observed colour temperature could, for example, correspond to larger grains closer to the star or smaller ones further away. The success in finding trends and characterising a typical debris disk depends crucially on the few objects where these degeneracies have been broken. Good knowledge of a set of such archetypes eases the navigation in parameter space and allows for the calibration of simpler models.

In this paper, we discuss the modelling of the HD 207129 system. The host star is nearly sun-like with an estimated age of 1.5–3.2 Gyr (Marshall et al. 2011, and references therein). Excess emission was first detected by IRAS (Walker & Wolstencroft 1988), and later observed by ISO (Jourdain de Muizon et al. 1999), Spitzer (Trilling et al. 2008), and APEX (Nilsson et al. 2010). Extended emission is seen and particularly the observations of scattered light with the Hubble Space Telescope (HST) (Krist et al. 2010) have shown the disk to feature an inclined ring with a reported radius of 163 AU. The presumed underlying belt of planetesimals marks the upper end of the range of radial extents of known debris disks. The best images so far have been obtained with Herschel/PACS (Poglitsch et al. 2010) and SPIRE (Griffin et al. 2010) in the framework of the DUNES (DUst around NEarby Stars) survey (Eiroa et al. 2010). They clearly show a ring-like outer disk, separated from the star by an inner depletion zone. This inner hole and the sheer size of the disk render HD 207129 an interesting object to study.

The next sections address the modelling process along the following steps. Section 2 reviews the observational data and the deconvolution results for the PACS images in order to judge where the emission comes from. In Sect. 3, we employ a power-law fitting to model the distribution of dust that produces this emission. Dust dynamics and the collisional cascade in which the dust is produced from parent bodies are included in Sect. 4. We then discuss these modelling results with respect to the origin and evolutionary history of the system in Sect. 5. Finally, in Sect. 6, a summary is given.

2. Observational input and deconvolution

Our modelling is centred on the observational results presented in Marshall et al. (2011). These data were obtained with Herschel/PACS and SPIRE in scan map and small map modes, respectively. With the photometry from both instruments, the SED is more densely sampled, and PACS images (see Fig. 1) additionally yield spatially resolved brightness information. Photospheric estimates for modelling and star subtraction are based on PHOENIX/NextGen models (Hauschildt et al. 1999), normalised to the short-wavelength, non-excess part of the Spitzer/IRS spectrum.

thumbnail Fig. 1

Herschel/PACS images of HD 207129 at (top to bottom) 70, 100, and 160 μm. Columns a) and b) show the observed and star-subtracted images, respectively. Columns c) to f) show deconvolved images using the following methods: c) van Cittert, d) FFT with a Wiener filter, e) MCS code, f) Richardson-Lucy. Flux steps (in mJy/arcsec2) between solid isolines are indicated in the top-right corner of each panel. The photocentre of the disk, i.e. the assumed stellar position, is indicated by a cross. In every panel, North is up and East is left.

2.1. Spectral energy distribution

A broad set of photometric data is available from different instruments and is shown in Table 1. Shortward of about 25 μm, only the HST detected excess emission (Krist et al. 2010). The Spitzer/IRS spectrum then clearly samples the rise of the SED (see, e.g., Krist et al. 2010). To homogenise the data for our modelling purposes, we extracted a single photometric point at 32 μm.

Table 1

Photometry of HD 207129.

In the far infrared, at wavelengths beyond 100 μm, measurements by Spitzer and ISO are contaminated by a nearby (1′) background galaxy to the north west. At submillimetre wavelengths, JCMT/SCUBA provided upper limits, while APEX/LABOCA detected a faint excess (Nilsson et al. 2010).

2.2. Radial profiles

We used the background-subtracted images at 70, 100, and 160 μm to derive observed radial profiles using the following steps. First, the photocentres were determined for all three wavelengths based on a fitting of an ellipse to the outer isophotes (semi-major axis of 13′′). Second, linearly interpolating pixel values along the major and minor axes of the ellipses, four cuts were made for each ellipse from its centre towards 56°, 146°, 236°, and 326° from north. The stellar contribution was not removed. The error estimates were based on squared sums of the background standard deviations of fluxes within  [20′′,50′′]  from the disk centres and the differences between the two opposing directions per axis (NW versus (vs.) SE and SW vs. NE). In the modelling process, profiles were then truncated at the respective distances where the signal drops below the noise.

Figure 2 shows the extracted profiles and their estimated uncertainties. We note that the profiles at 70 and 100 μm feature a consistent flux discrepancy of about 15% (4σ) between the north-western and the south-eastern ansae. This discrepancy is not directly accounted for in the models because these are axisymmetric. However, we discuss the possible nature of this asymmetry in Sect. 5.

thumbnail Fig. 2

Radial profiles of surface brightness at (top to bottom) 70, 100, and 160 μm: (left) along the long axis and (right) along the short axis. The coloured shaded areas mark the uncertainty, while symbols represent the two respective directions along each axis. In the lower-left of each panel, the stellar contribution is plotted. Azimuthally averaged radial profiles of Wiener deconvolved and deprojected images are given as solid line, and their azimuthal standard deviations are represented as dotted lines, both in arbitrary units.

2.3. Deconvolution

For a first approach to what kind of disk could produce the PSF-blurred images that are observed, we deconvolved them using different methods: (1) a classical reverse Fourier transformation (FT) with Wiener-Kolmogorov filtering (Kolmogorov 1941; Wiener 1949) of the noise, (2) an FT-based van Cittert scheme, (3) the MCS code (Magain et al. 1998), and (4) a Richardson-Lucy (R-L) scheme (Richardson 1972; Lucy 1974). In Cols. (c) to (f) of Fig. 1, the resulting images are shown, and the averaged profiles for the van Cittert scheme are plotted in Fig. 2. All methods predict an outer ring and a gradual depletion towards the disk centre.

The clumpy structure that is visible in the deconvolved images does not necessarily correlate with real structure in the disk. On the one hand, some clumps are apparent at every wavelength and with both deconvolution schemes. The most prominent example is the brightness peak in the western ansa. This feature is clearly visible in the original images and the resulting radial profiles. However, in the MCS-deconvolved image at 100 μm, this feature is arc-like and double-peaked. On the other hand, some clumps appear in only one or two images and are, therefore, probably noise fluctuations.

For comparison, Fig. 3 shows the outcome of an even simpler test: the assumption of constant surface brightness from 70 to 190 AU in all three PACS bands (meaning constant temperature in that region) and zero flux beyond these bounds. Despite the very limited number of parameters, this purely empirical model is able to reproduce the observed radial profiles within the error margins. The surprisingly good agreement of the convolved model profiles with the observations has two immediate consequences. First, from these radial profiles alone, it is impossible to differentiate between a completely empty inner hole, i.e. a sharp inner edge at about 70 AU ( ≈ 4.5′′), and one that contains some material. Second, the disk is not compatible with the standard idea of a surface density that peaks near an inner disk edge and falls off with increasing radial distance. Instead, the cross-section density is highest near the outer edge of the disk. This is examined in greater detail later in this paper.

thumbnail Fig. 3

Comparison of (convolved and unconvolved) rectangular surface brightness profiles with the observations and their deconvolutions. Wavelengths (top to bottom): 70, 100, and 160 μm. The dashed and dotted lines are given in arbitrary units.

3. Parametric modelling of dust properties

To provide a first estimate of the properties of the underlying dust that produces the observed emission, we started with an empirical fitting of the size and radial distributions of dust. For this fitting we used both a scheme based on simulated thermal annealing (Press et al. 1992, SAnD; see S. Ertel et al. 2011) and an extensive Bayesian approach with χ2 mapping of 37 million models across the parameter space (GRaTer; see Augereau et al. (1999) and Augereau et al. (in prep.), on models of the q1 Eridani debris disk).

3.1. Parameters

We adopted a basic model for the surface number density n, i.e. the number of particles per unit size and (normal) disk area, consisting of two separate power laws for size and radial distributions n(s,r)=n0(ss0)κ(rr0)αsκrα,\begin{equation} n(s, r) = n_0 \left(\frac{s}{s_0}\right)^\kappa \left(\frac{r}{r_0}\right)^\alpha \propto s^\kappa r^\alpha,\label{eq:PowerLaw} \end{equation}(1)where the grain sizes s were assumed to be bound to a range  [smin,smax]  and the distances r to  [rmin,rmax] . Beyond these boundaries, we had n = 0. The boundaries and the exponents were free parameters. The constants s0 and r0 were freely chosen such that s0 = smin and r0 = rmax (or, alternatively, r0 = rmin). The normal geometrical optical thickness of such a disk is given by τ(r)=sminsmaxπs2n(s,r)ds=πsmin3n03+κ(rrmax)α[(smaxsmin)3+κ1],\begin{equation} \tau(r) = \int_{s\sbs{min}}^{s\sbs{max}}\limits \pi s^2 n(s,r) \d s = \frac{\pi s\sbs{min}^3 n_0}{3 + \kappa} \left(\frac{r}{r\sbs{max}}\right)^\alpha \left[\left(\frac{s\sbs{max}}{s\sbs{min}}\right)^{3+\kappa} - 1\right], \end{equation}(2)for κ ≠  − 3. The mass for one such dust disk is given by Mdust=8π2ρn0rmax2smin43(2+α)(4+κ)[1(rminrmax)2+α][(smaxsmin)4+κ1]\begin{equation} M\sbs{dust} = \frac{8\pi^2\rho n_0 r\sbs{max}^2 s\sbs{min}^4}{3(2+\alpha)(4+\kappa)} \left[1 - \left(\frac{r\sbs{min}}{r\sbs{max}}\right)^{2+\alpha}\right] \left[\left(\frac{s\sbs{max}}{s\sbs{min}}\right)^{4+\kappa} - 1\right] \end{equation}(3)for κ ≠  − 4 or otherwise Mdust=8π2ρn0rmax2smin43(2+α)[1(rminrmax)2+α]lnsmaxsmin·\begin{equation} M\sbs{dust} = \frac{8\pi^2\rho n_0 r\sbs{max}^2 s\sbs{min}^4}{3(2+\alpha)} \left[1 - \left(\frac{r\sbs{min}}{r\sbs{max}}\right)^{2+\alpha}\right] \ln \frac{s\sbs{max}}{s\sbs{min}}\cdot \end{equation}(4)Note that, for GRaTer, the radial distribution was actually assumed to be a seemless combination of two of the above power laws with a smooth turnover at rbreak. Here, the innermost boundary is set at the size-dependent sublimation radius, rsub, because rmin values turned out not to influence the overall fit to the SED and profiles with GRaTer (rmin values between 10 and 80 AU have been explored).

In addition, for both fitting schemes, the chemical composition was varied in the sense that a mixture of water ice (Li & Greenberg 1998) and astronomical silicate (Draine 2003) was assumed, with the volume fraction of ice, σice = Vice/(Vice + VSi) = 0...1, being another free parameter. Here, Vice and VSi denote the bulk volumes filled by ice and silicate, respectively. The resulting bulk densities of the mixtures are given by ρmix=σiceρice+(1σice)ρSi,\begin{equation} \rho\sbs{mix} = \sigma\sbs{ice} \rho\sbs{ice} + (1-\sigma\sbs{ice}) \rho\sbs{Si}, \end{equation}(5)where ρice = 1.2 g cm-3 and ρSi = 3.5 g cm-3 denote the bulk densities of water ice and silicate, respectively.

The ranges across which we explored the parameter values are given in Table 2.

Table 2

Results of power-law fitting.

3.2. Weighting of the data and uncertainties

The absolute flux calibration enters the total photometric flux as well as each individual image pixel and each data point in the derived radial profiles. However, physically, it is represented by one single number for each wavelength. Therefore, we also modelled the deviations in absolute scaling just once per wavelength, i.e. only for the SED. The radial profiles were then compared only in terms of their relative shapes. The discrepancies in absolute scaling can be expressed in terms of wavelength-dependent correction factors Cλ. We have Fi,λ(rj)=CλFi,λ(rj),\begin{equation} F_{i,\lambda} (r_j) = C_\lambda F'_{i,\lambda} (r_j), \end{equation}(6)where Fi,λ and Fi,λ\hbox{$F'_{i,\lambda}$} are the respective absolute and relative surface brightnesses at wavelengths λ along axis/direction i at distance rj from the centre.

This approach is also motivated by uncertainties in the optical properties. Assume, for example, that our derived ratio between fluxes at 70 and 100 μm is systematically off by 10%. Having such a systematic offset combined with the absolute values of the radial profiles would cause the minimization algorithms to yield fits that preserve the total flux (as desired for the SED) at the cost of distorting the relative shapes of the profiles.

Modelling uncertainties are introduced by the intrinsic limitations of assuming the dust grains to be homogeneous spheres and by the optical constants used (owing to, for instance, simplified composition or mixing). An example of the latter is given in Fig. 4, where the imaginary parts of the refractive indices of astronomical silicate (Draine 2003) are compared with standard references for amorphous and crystalline water ice (Warren 1984; Hudgins et al. 1993; Li & Greenberg 1998). It is obvious that a significant portion of the data in the wavelength region of interest is based on a power-law extrapolation rather than experimental results. However, a discrepancy of 0.5 in terms of the slope of the imaginary part of the refractive index translates into a 20% discrepancy in terms of the ratio between fluxes at 70 and 100 μm (ignoring the change of dust temperature that goes along). At long wavelengths, such discrepancies are easily introduced when the optical data are based on layer thicknesses (Hudgins et al. 1993) or particle sizes (Draine 2003) that are much smaller than the wavelength.

thumbnail Fig. 4

Imaginary part of the refractive indices in the far infrared for (solid line) pure amorphous water ice as compiled by Li & Greenberg (1998), (symbols) water ice at different temperatures from Hudgins et al. (1993, the original data source for Li & Greenberg), (dashed line) crystalline ice as compiled by Warren (1984), and (dotted line) astronomical silicate synthesised by Draine (2003).

For all models, the deviations from the observational data were characterised by means of a reduced χ2. This χred2\hbox{$\chi\sbs{red}^2$} is calculated from the sum of the squares of the deviations at individual data points, namely 17 SED points (λ > 20 μm) listed in Table 1 and 96 points from the six radial profiles shown in Fig. 2. The sum of squares is then divided by the degree of freedom, i.e. the total number of data points reduced by the number of free parameters (including the Cλ factors): 113 − 11 = 102.

We note that the noise contribution from the sky background was assumed to be uncorrelated in this study, whereas the true uncertainty in a given pixel is partly correlated with that of neighbouring pixels. One reason for this correlation is the native instrument pixel size (3.2′′ at 70 and 100 ¯m, 6.4′′ at 160 ¯m) being larger than that of our scan maps (1′′ at 70 and 100 ¯m, 2′′ at 160 ¯m). Therefore, a model with, for example, two adjacent pixels that are each off by  + 0.5σ from the observations would actually fit better than a model with one of the two pixels off by  − 0.5σ and the other by  + 0.5σ. Neglecting the correlation, we tend to slightly overestimate uncertainties and underestimate χred2\hbox{$\chi\sbs{red}^2$}. This is reflected by the observed radial profiles in Fig. 2 being smoother than the corresponding confidence intervals would suggest.

3.3. Results

The outcome of the formal best-fit to the observations is shown in Figs. 5 and 6. The values of the individual parameters are listed in Table 2. The setup includes dust with an ice content σice = 50% and a surface density rising from rmin = 57 AU towards an outer edge at rmax = r0 = 194 AU with a slope α = 2.2.

Such a radial profile that rises outward with a sharp outer edge contrasts with those derived for the archetypal disks Vega (Su et al. 2005; Krivov et al. 2006; Müller et al. 2010; Sibthorpe et al. 2010) and β Pictoris (Golimowski et al. 2006; Thébault & Wu 2008; Krivov et al. 2009; Vandenbussche et al. 2010). In these disks, the derived geometrical optical thickness typically falls off with a moderate r-1.5 relationship beyond the peak value, i.e. they do not feature sharp outer edges. Disks are instead often assumed to have inner edges as sharp as that of Fomalhaut (Kalas et al. 2005, 2008), for example.

In our power-law SAnD model, the best-fit slope of the size distribution is κ ≈  −3.8, which is in good agreement with the value of Krist et al. (2010) and slightly steeper than both the canonical value of  −3.5 (Dohnanyi 1969) and the typical slope of  −3.7 expected for collisional equilibrium in the strength regime (O’Brien & Greenberg 2003). The peak geometrical optical depth is τ = 4 × 10-4 and the dust amounts to a total of 7.2 × 10-3   M in objects with radii smaller than 1 mm. On the basis of visual inspection, the agreement with the observations is good. Quantitatively, it boils down to a reduced χ2 of 0.77. Considering only the PACS and SPIRE wavebands, the agreement is worst at 100 μm, where the SED reveals a total flux deficit of 14% compared to the observations.

thumbnail Fig. 5

Spectral energy distributions of the best-fit power law-based models. The photometry points used for the fitting are overplotted. Note that the horizontal stretch is different for wavelengths smaller than 40 μm. To compare the models to the IRS data from the standard pipeline, we corrected the model SEDs for slit loss and then applied a point-source correction. The uncorrected model SEDs are plotted with light grey lines.

thumbnail Fig. 6

Radial brightness profiles of the best-fit power law model: (top) 70 μm, (middle) 100 μm, (bottom) 160 μm. Cuts along the long and short axes are given in the left and right columns, respectively.

To clarify the influence of the fraction of ice, it is worth comparing this solution to the results of an only slightly worse fit (at a reduced χ2 of 0.78, see Table 2) that involves no ice at all. For this SAnD model, the size distribution is slightly steeper (κ =  −4) and the corresponding total dust mass marginally higher at 8.3 × 10-3   M. The optical depth peaks at 4 × 10-4. The outcome is shown in Figs. 5 and 6 with solid lines. The most notable difference to the overall best fit is the closer agreement with the observed total flux at 100 μm. Without ice, the excess at shorter wavelengths is also stronger. However, considering the uncertainties involved, the differences are marginal.

It should be noted that the resultant slopes of the size distributions partly inherit the uncertainty in the slope of the refractive index used for modelling (see Sect. 3.2). Therefore the differences between κ =  −4.0 for pure silicate and κ =  −3.8 for the mixture with ice are likely due to the absorption efficiencies for silicate (Draine 2003) falling off more steeply with increasing wavelength than those of amorphous water ice (Li & Greenberg 1998).

The GRaTer fit (Table 2, Figs. 5 and 6) shows that a narrower disk and a high content of ice (σice = 90 %) are also compatible with the observations. (Note that SAnD did not explore ice fractions as high as that.) The advantage of this fit is that its location of peak brightness agrees with the ring radius derived from the HST images. However, the lack of emission in the inner regions and at short wavelengths leads to a marginally higher χred2\hbox{$\chi\sbs{red}^2$} of 1.26. Owing to the ring being narrower, the peak geometrical optical thickness is higher at 1.4 × 10-3.

When comparing with the model proposed by Krist et al. (2010) for Spitzer images and spectra, the Herschel images provide stronger contraints than the SED (Figs. 5 and 6). While the SED is not that selective, the radial profiles clearly show that a disk with an inner edge at 148 AU cannot reproduce the emission seen closer to the star.

Summarising our findings, the best-fit model is a disk where most of the material is concentrated near the outer edge, which is located at or beyond the ring observed with the HST by Krist et al. (2010). However, we also probably detect material at least as close as about 60 AU from the star. Models of broader rings fare better than models with very narrow rings. Both GRaTer and SAnD predict minimum grain sizes well above the blowout limit (by factors from four to ten).

4. Modelling from the sources

We now seek to more tightly constrain the properties of the dust disk by modelling the production and loss of material in the collisional cascade that is assumed to act in a debris disk. In this approach, the distributions of sizes and radial distances of the dust are coupled. Parameters such as the radial extent of the dust disk can no longer be directly controlled but are instead replaced as free parameters by the radial extent of the underlying planetesimal disk.

4.1. Basic model

The approach is translated into a kinetic problem, where collision rates and outcomes are described statistically. The numerical solution to this problem is implemented in a parallel C++ code (ACE, Analysis of Collisional Evolution).

The code works on a three-dimensional grid with logarithmic bins for object mass and pericentre and linear bins for the orbital eccentricity. It takes into account material strength, mutual gravity, and different relative orientations of the orbits of the colliding particles. The vertical distribution in the disk is assumed to be uniform within a constant opening angle. The resulting coupled size-radial distribution is then again used as an input to Mie calculations of the thermal emission.

For the collisional outcomes, the size-dependent ratios of stellar radiation pressure to gravitational pull are taken into account: the β ratios. Fragments are launched into orbits that become more eccentric for higher β. When launched from initially circular orbits, grains with β > 0.5 acquire eccentricities e > 1 and become “blowout grains” as they leave the system on unbound orbits (e.g., Burns et al. 1979). For HD 207129, this blowout limit is reached at sizes s ≈ 0.6 μm. We refer to Fig. 7 for a comparison of β for different material compositions and grain sizes. We note that, depending on the orbits that they are launched from, grains can potentially stay bound up to β < 1.0.

thumbnail Fig. 7

The β ratio as a function of grain size for different mixing ratios of ice and silicate, based on Mie theory (see, e.g., Bohren & Huffman 1983; Burns et al. 1979).

A more detailed description of the model can be found in previous papers (Krivov et al. 2005, 2006, 2008).

4.2. Improvements made for this work

The critical specific energy for disruption and dispersal, QD\hbox{$Q\sbs{D}^*$}, is a convenient parameter in the description of collisions. However, it is not a material constant but it is well-known to depend on the object size; this dependence is only mild in the strength regime, where cohesive binding forces need to be overcome, but stronger in the gravity regime that determines the break-up of planetesimals (see, e.g. Benz & Asphaug 1999). Furthermore, from experiments and theoretical considerations, this critical energy is also known to depend on the colliders’ mutual impact velocity. This velocity dependence is not commonly used in models for the collisional disk evolution (Krivov et al. 2006; Thébault & Augereau 2007; Kenyon & Bromley 2008), although relative velocities easily drop by an order of magnitude from 1 to 100 AU. Inspired by the elegant reformulation of QD\hbox{$Q\sbs{D}^*$} of Stewart & Leinhardt (2009)  in terms of  reduced  and total  mass of the

colliders, we adopted their QRD\hbox{$Q\sbs{RD}^*$} for our QD\hbox{$Q\sbs{D}^*$}, including the dependence on the impact velocity vimpQD=[QD,s(s1m)bs+QD,g(s1km)bg](vimpv0)0.5,\begin{equation} Q\sbs{D}^*=\left[Q\sbs{D,s}\left(\frac{s}{1\,\text{m}}\right)^{b\sbs{s}}+Q\sbs{D,g} \left(\frac{s}{1\,\text{km}} \right)^{b\sbs{g}} \right] \left(\frac{v\sbs{imp}}{v_0}\right)^{0.5}, \end{equation}(7)where subscripts “s” and “g” stand for strength and gravity regimes, respectively. Setting v0 = 3 km s-1, we again used base coefficients QD,s, QD,g, bs, and bg as derived by Benz & Asphaug (1999). Since most of the collisions take place within the belt of parent bodies, i.e. with similar impact velocities, the main difference to previous calculations is a slightly reduced average QD\hbox{$Q\sbs{D}^*$}. The effects thereof are discussed in Sect. 4.6.

4.3. Reference run and parameter variations

Following the strategy outlined in Müller et al. (2010), we start from a first-guess setup based on the roughly known belt location (Krist et al. 2010; Marshall et al. 2011). After computation of the resulting SED and the radial profiles, we manually alter the most important (and least constrained) disk parameters and re-run the model until we reach a point from where no further significant improvements can be achieved that way. This set of important parameters is comprised of: (1) the range of parent bodies’ orbital semi-major axes, a, i.e. the belt location and extent; (2) their maximum orbital eccentricities, emax; (3) the dust grains’ critical specific energy for disruption, QD,s\hbox{$Q\sbs{D,s}^*$}, i.e. their collisional strength; and (4) their volume fraction of ice, σ. The orbital inclinations i and the disk semi-opening angles ϵ are set to ϵ = i = emax/2. The parameter space for each run is sampled by 70 grid points in mass (from 2.5 × 10-14 g to 4.2 × 1021 g, corresponding to 0.12 μm to 66 km for pure silicate with ρ = 3.5 g cm-3), 42 grid points in pericentric distance (spanning a range from 13 AU to 610 AU), and 50 to 150 bins in eccentricity (from 0.0 to 5.0, in steps of 0.033 to 0.1).

In agreement with the results of Sect. 3.3, we find a reference setup where most material is concentrated near the outer disk edge. The initial slope of the radial profile of optical depth (or surface mass density) is α = 3.5, i.e. rising steeply outwards. The maximum eccentricity of the larger grains and the parent bodies is emax = 0.05, corresponding to a low dynamical excitation, i.e. low collisional velocities.

The model disk contains a total of 7.6 × 10-3   M ⊕  in dust (with radii below 1 mm). The total rate of mass loss is 3 × 1010 g s-1 with a share of only 0.1% accounted for by drag, i.e. Poynting-Robertson (P-R) drag and stellar wind drag, across the inner edge of the grid at 12 AU. The remaining 99.9% is lost in the form of the very smallest grains: the blowout grains, which are steadily produced in collisions and leave the system on unbound orbits due to the stellar radiation pressure.

Starting from that reference point, we produced alternative solutions by varying the individual parameters. These variations also explore the range to which the parameters can be constrained and the degeneracy among them. Table 3 lists all the runs discussed here. Aside from the more obvious parameters such as the belt location (in terms of semi-major axes), that table also lists whether drag forces were considered.

In general, the effects through which these variations may alter the observed emission can be classified according to three, partly entangled categories: (1) changes in the size distribution; (2) changes in the radial distribution; and (3) direct changes in the optical properties. While the Q and e runs mainly belong to the first category, a and, in particular, d belong to the first and the second category. The c runs can be mainly identified with the third category, although the size distribution is also affected through the change in blow-out radius.

The simulated times of the individual runs do usually not correspond to the physical age of the system (1–3 Gyr) because all simulations were run as long as necessary to ensure that the dust is in collisional equilibrium and to reach a dust mass that is compatible with the observations. For disks with only collisions at work, the change rates of particle numbers are purely quadratic in the particle numbers per bin, n, because collisions are assumed to be a two-body process. In that case, the total initial mass could be scaled after the run to a value where real and simulated times matched (Krivov et al. 2008). However, where drag forces are important, the change rate is no longer quadratic in n, but contains additional linear terms i=j,kGijknjnkjLijninj+jTjnj,\begin{equation} \dot{n}_i = \sum_{j,k}\limits G_{ijk} n_j n_k - \sum_j\limits L_{ij} n_i n_j + \sum_j\limits T_j n_j,\label{eq:GainLoss} \end{equation}(8)where Gijknjnk and Lijninj denote the gain and loss of particles in bin i through collisions among particles in bins j and k as well as i and j, respectively. The transport to and from (neighbouring) bins j is denoted by Tjnj. As a result, the straightforward mass scaling presented in Krivov et al. (2008) does not work, and we would need to iteratively find the proper initial mass to get the desired dust mass after a time that corresponds to the physical age of the system. Given that the initial size distribution of planetesimals is essentially unknown, we refrain from performing this iteration and focus on the equilibrium dust distribution.

Table 3

Description of numerical runs.

We note that the reference run is a good fit within the parameter space and the constraints of this model. However, a more precise knowledge and description of the appropriate collisional physics would likely favour slightly different model parameters. The same is true for a deviation from the crude assumption of homogeneous spheres. Owing to these uncertainties, we refrain from any additional fine tuning of the model parameters since this would yield no additional information about the specific object.

4.4. Reference SEDs and radial brightness profiles

The radial profiles in Fig. 8 and the SED in Fig. 9 illustrate that this approach can produce an overall agreement with the observations that is comparable with the power-law fits. A lack of flux at 100 μm is the only significant discrepancy between the reference run and the observations. In agreement with the results from Sect. 3.3, the flux at 100 μm is underpredicted by about 18%. A correction factor of C100 = 1.23 is needed.

thumbnail Fig. 8

Radial brightness profiles of the reference run.

thumbnail Fig. 9

Same as Fig. 5 but for the spectral energy distribution of the reference run and the contributions from its constituent annuli and the star. For comparison, the dragless case (run d) is dashed.

In terms of flux, the small amount of material that traverses the gap towards the star produces a significant contribution only at shorter wavelengths and close to the star. Figures 9 and 10 distinguish the contributions from within and outside the belt’s inner edge to the modelled SED and the radial profile at 70 μm. In agreement with the deconvolution results, the contribution of the inner regions to the total flux is minor in the PACS 100 and 160 μm and SPIRE wavebands.

thumbnail Fig. 10

Radial brightness profiles of the reference run. Shaded regions mark the contributions from the star and from constituent annuli of the disk.

In addition, when compared to the observations, the model predicts a slightly shallower 70-μm slope around the outer disk edge, i.e. at 10–13 arcsec from the centre along the long axis. At 100 μm, the model disk falls off too steeply at 8–11 arcsec and/or is too narrow overall by about 4%.

At 160 μm, the model features a more pronounced plateau or even a depression close to the disk centre, meaning that there is a marginal lack of long-wavelength emission coming from the drag-dominated regions. With the material dragged into these inner regions being hotter, it is more difficult to fill the gap at longer wavelengths – despite the PSF-blurring being stronger. Additional emission from the neglected inner 20 AU would strengthen this trend, which is also evident from the power-law fits (see Fig. 6).

4.5. Size and radial distribution of the dust

In the left panel of Fig. 11, the profile of the normal optical thickness is shown to deviate from the profile typically expected for collisionally very active debris disks. The models of these disks would predict the optical thickness to be dominated at all radial distances by a halo of small grains that are slightly larger than the blowout grains. The stellar radiation pressure pushes these grains to eccentric and extended orbits, and the resulting outer profile of the optical depth follows a power-law slope of  −1.5 (Krivov et al. 2006; Strubbe & Chiang 2006; Thébault & Wu 2008).

thumbnail Fig. 11

(Left) Radial profiles of optical thickness for two runs: reference and iceless. For comparison, the light shaded region shows the profile of a model of the Vega disk (Müller et al. 2010), radially scaled so that the outer edges match. The dark shaded region marks the extent of the belt of parent bodies. (Right) Size distribution of the reference run at distances of 50, 100, 150, and 200 AU.

thumbnail Fig. 12

Optical thickness as a function of radial distance and grain size for the reference run: (left) 3–23 μm, (right) 0.4–1.8 μm.

Here, we report a slope for this disk that is significantly steeper near the belt’s outer edge. As discussed by Thébault & Wu (2008), this phenomenon is caused by the low dynamical excitation of the system, which results in very different collisional probabilities and impact velocities among the different grain populations. The big grains have their collision rates set by the average “dynamical” eccentricity in the system. In the reference run, we found relatively low eccentricities (emax = 0.05), such that the rate of destruction of big grains is low. The resulting production rate of small grains is also low. Their destruction rate, however, is high because it depends on the dynamical excitation of these small grains, which is always high because it is set by the high eccentricities they get from radiation pressure. As a result of this imbalance between destruction and production rates, their abundance is reduced. In consequence, the belt region is dominated by larger grains, as represented by peaks at about 7 and 4 μm in the size distributions at respective distances of 100 and 150 AU, seen in the right panel of Fig. 11. The halo gets suppressed.

It is only further away from the star and the belt that the model profiles flatten and then approach r-1.5. The size-dependent radial profiles in Fig. 12 show that, in these outer regions, the optical thickness is dominated by an agglomeration of barely bound grains at sub-micron sizes around the maximum in the β function (where still β < 1). The profile overplotted in the left panel of Fig. 11 for the iceless run shows that the exact shape depends on the actual model parameters used.

In contrast to the sharp outer edge, drag forces lead to a rather shallow drop at and within the inner edge of the parent belt. The radial profile of the reference run peaks around the outer edge of the planetesimal disk and slowly slopes towards the inner edge of the belt, where it steepens. Within about 90 AU it then flattens again because of drag forces. However, the optical depth of the reference run does not remain constant in these inner regions, as would be expected for pure transport domination. The profile instead roughly follows a r0.9 profile. This can be explained by a mix of collisions and transport, as illustrated in the size-dependent radial profiles in the right panel of Fig. 12. In the inner regions, grains with radii smaller than 2 μm are the most affected by transport and exhibit a characteristic flat profile, while larger grains suffer from ongoing collisional grinding. As a result of the disk’s low dynamical excitation and the overall depletion of small grains, these larger grains dominate the radial profile even in the inner regions modelled here. However, additional flattening towards a constant profile is expected to occur at distances below 10–20 AU, where only smaller grains ( ≲ 2 μm) remain.

Figure 13 summarises these findings and shows that the relation between size and radial distribution is indeed non-trivial. In the case of the reference run, the peak in the size distribution shifts from 3 μm at 20 AU to about 7 μm (i.e. ten times the blowout radius) at the inner edge of the planetesimal belt. Near the belt’s outer edge, the peak is again at 3–4 μm and it drops to about 1 μm in the outer region, the halo that is filled with barely bound grains.

thumbnail Fig. 13

Colour-coded distribution of geometrical optical depth, i.e. particle cross-section, over grain size and distance to the star given for the reference run. The white solid line gives the peak position in the size distribution as a function of radial distance. The horizontal dotted line indicates the blowout size, where β = 0.5.

4.6. Rate and efficiency of collisions vs. transport

Since the size and radial distributions of the dust are entangled, setups with similar profiles can have different SEDs and vice versa. The reason for this behaviour is that the SEDs are more sensitive to grain temperatures, hence grain sizes than the relative radial profiles are. The influences of these parameters are depicted in Figs. 14 and 15, where we compare the results of the reference run to a set of six runs: two for QD\hbox{$Q\sbs{D}^*$} lowered and raised by a factor of 3, two for emax raised to 0.1 and 0.2, respectively. In addition, the dragless run and the run where QD\hbox{$Q\sbs{D}^*$} is independent of radial distance are shown. The radial profiles show the similarities between the effects of a lower collisional strength (lower QD\hbox{$Q\sbs{D}^*$}) and larger relative velocities (higher emax), and vice versa.

thumbnail Fig. 14

Size distribution in terms of optical depth at 150 AU as a function of grain size. The reference run is compared with runs (left) Q1 and Q2, as well as (right) e1, e2, and d.

thumbnail Fig. 15

Same as Fig. 10 but for a comparison of the (solid line) reference run with runs Q1, Q2, e1, e2, and d. See Table 3 for a description of these runs.

While the combined effects of material strength and dynamical excitation are thought to dominate the rate of destructive collisions (Thébault et al. 2003; Krivov et al. 2006), their individual influences differ. The left panel of Fig. 16 shows that the higher dynamical excitation in run e2 strongly increases the typical temperature by reducing the typical grain size. By contrast, the lower QD\hbox{$Q\sbs{D}^*$} in run Q1 leads to an SED that differs only marginally from the reference, while leading to radial profiles that are almost identical to those in run e2. As discussed in Sect. 4.5, the reason is that the dynamical excitation determines the ratio of the number of large grains, which inherit the excitation of the parent bodies, to the number of small grains that are affected by radiation pressure (Thébault & Wu 2008). When varying the coefficient QD,s\hbox{$Q\sbs{D,s}^*$}, the particle strength influences the collision rate similarly at all sizes. Krivov et al. (2006) and Thébault & Augereau (2007) illustrated the impact of QD,s\hbox{$Q\sbs{D,s}^*$} on the waviness of the size distribution. As a result of the low dynamical excitation and inclusion of cratering collisions, the size distributions in our models for HD 207129 (Fig. 14) do not have strong waviness.

thumbnail Fig. 16

Spectral energy distribution of the reference and runs (top) Q1 and Q2 as well as (bottom) e1 and e2.

By considering both the SED and the radial profiles, we can constrain the dynamical excitation of the disk to values of emax < 0.1, i.e.  ⟨ e ⟩  < 0.05 and  ⟨ i ⟩  < 0.025 = 1.5° (roughly assuming that e and i are distributed uniformly). The excitation is lower than what is typically expected for efficient stirring by a planet (see, e.g. Gomes 2003; Wyatt 2003) and lower than what is inferred from a similar model of the Vega disk (Müller et al. 2010, where emax = 0.2). For the classical, dynamically cold population of the Edgeworth–Kuiper belt, different authors have proposed average eccentricities  ⟨ e ⟩  ≈ 0.08 and inclinations  ⟨ i ⟩  ≈ 2.7° (Brown 2001; Elliot et al. 2005; Vitense et al. 2010), values that are again higher than what we infer for HD 207129. As a result of the low orbital velocities at distances of 160 AU, the differences in terms of random velocities, i.e. typical collisional velocities, are even larger. For our reference run, where  ⟨ e ⟩  ≈ 0.025, we have  ⟨ vcol ⟩  ≈  ⟨ e ⟩ vKepler ≈ 60 m s-1. With the Kuiper belt being located at a distance of only about 40 AU, the resulting relative velocities there are six times as large as for HD 207129.

The relative importances of the collisions and transport in a disk change with time. That collisions are a two-body process makes their typical timescale depend on the amount of material in the system – in contrast to transport, where the timescale n/ is constant (cf. Eq. (8)). Therefore, drag becomes increasingly important the older a system gets and the less dust it contains. In consequence, not only are the absolute scaling of dust mass and flux time-dependent but so are the radial distribution, the size distribution, and the spectral energy distribution. The evolution of the SED and the long-axis profile of brightness are given in the left and right panels of Fig. 17, respectively. While the emission from the belt itself quickly drops with time, the inner regions remain almost unchanged. The flux and optical depth there are limited by the maximum level that is possible before collisions dominate. See also Wyatt (2005).

thumbnail Fig. 17

Evolution in time of run c1: (left) the SED and (right) the radial brightness profile at 70 μm.

4.7. Ice fraction of dust and location of planetesimal belt

In terms of their SEDs, the models most widely differ at shorter wavelengths. At wavelengths beyond 100 μm, grains smaller than 10 μm are almost invisible. In consequence, the contribution of barely bound grains and grains that drift inward to submillimetre fluxes is negligible. Only the larger grains that remain within or near the birth ring matter. Taking into account the similarity in their size distribution, as shown in Fig. 14, the convergence of almost all runs in the submillimetre is a natural outcome.

In common with the power-law approach, models that contain ice fare somewhat better than models without ice. Figures 18 and 19 show that the chemical composition influences the emission at both short and long wavelengths. For example, assuming pure silicate instead of the reference mixture with 70% ice, the SED becomes broader and the results conflict with the IRS data, in particular. Apart from that, considering the uncertainties involved with the optical constants (see Sect. 3.2), no strong constraints can be placed on the exact volume fraction of ice.

thumbnail Fig. 18

Spectral energy distribution of the reference and runs (top) c1 and c2 as well as (bottom) a1 and x.

thumbnail Fig. 19

Same as Fig. 10 but for a comparison of the (solid line) reference run with runs c1, c2, a1, and x. See Table 3 for a description of these runs.

While volatiles contribute significantly to the mass of large transneptunian objects, their role in the chemistry and life cycle of debris dust is less obvious. Erosion by means of photo-sputtering and partly sublimation compete with removal through collisions and drag forces. Following the approach described in Grigorieva et al. (2007), we estimate the rate at which stellar UV photons sputter ice from a grain’s surface and reduce its size. The resulting timescale for grains of pure amorphous water ice (Li & Greenberg 1998) is plotted in Fig. 20, where it is compared to the size-dependent collisional lifetimes for three runs. As a result of the low dynamical excitation and the resulting low collision rates, this rough comparison shows that UV sputtering could be an important removal mechanism for (purely) icy grains smaller than a few tens of micrometres, thereby further increasing the effective grain size. The details of the combined influence of collisions and sputtering on the distributions of grain sizes and radial distances are however beyond the scope of this work.

thumbnail Fig. 20

(Grey solid line) UV-sputtering lifetimes compared with collisional lifetimes of objects in (black solid line) the reference run, (dashed line) run d, and (dotted line) run e1.

Within the set of parameters chosen for this study, the overall location of the planetesimal belt can be assessed more easily because the disk is sufficiently resolved. The radial profiles in the right panel of Fig. 19 show that shifting the whole disk either outward or inward by 5% already results in significant deviations.

However, we note that the disk location as such is actually a poorly defined quality. As in the demonstration in Sect. 2.3, opening up the parameter space further and allowing the initial material distribution to vary more freely would likely bring along a set of equally good solutions. The reference solution presented here is therefore not to be considered as unique. Other setups are conceivable. For example, the flux that comes from regions closer to the star might be emitted from dust that is produced in-situ instead of being dragged there. The underlying distribution of planetesimals might, for example, continue to yet smaller radii than assumed in the reference run. Alternatively, it could be more closely confined to the location of the ring seen at around 163 AU in the HST images (Krist et al. 2010). The latter scenario is simulated in run x and the resulting SED and radial profile are given in the right panels of Figs. 18 and 19.

Both the reference and the alternative run agree with the observations. Therefore, only the location of the outer edge or peak density of the planetesimal belt can really be constrained. However, it is impossible to distinguish between a narrower, ring-like belt and one that extends further in with a steeply declining surface density. Given the low dynamical excitation, transport notably fills the inner gap in both cases.

4.8. Scattered light

The most clearly resolved data available for HD 207129 were obtained with the HST (Krist et al. 2010) in scattered light at a wavelength of 0.6 μm, resulting in a diffraction limit of 0.06′′ (or 1 AU). The images show an inclined ring, 30 AU wide, centred around a radial distance of 163 AU ( ≈ 10′′). No strong constraints can be placed on the brightness and dust content in the inner regions.

For comparison, we calculated the brightness of three of our model disks (reference and runs Q1 and x) in scattered light, again using Mie theory and the algorithm of Bohren & Huffman (1983). The resulting images in Fig. 21 show that all three runs can reproduce the reported location of peak brightness, while only run x has a ring that is as narrow as the one deduced from the HST images.

thumbnail Fig. 21

Synthetic scattered light images of three of the runs: (left) the reference, (middle) Q1, and (right) x. Note that run Q1 is artificially truncated at a radial distance of 400 AU (25′′). All three disks are inclined by 60 degrees with respect to the sky plane, the north-eastern part being closer to the observer.

We note that owing to the strength of transport in the reference run and run x, these models predict the surface brightness in scattered light to rise again in the innermost regions, within about 3′′. If no inner planets or other removal mechanisms are invoked, a zodiacal cloud with an optical thickness of a few times 10-6 is the consequence.

Our models also predict a pronounced asymmetry, caused by grains significantly larger than the wavelength, about 5 μm compared to 0.6 μm. As shown above, these grains provide the dominant contribution to optical thickness. If we were to adopt the simplistic assumptions that they are both spherical and homogeneous, the parts of the disks that are closer to the observer would appear brighter owing to strong forward scattering. However, realistic grains have smaller scale substructure both internally and on their surface. The inclusion of these features would probably produce effective scattering properties that are again closer to those expected for isotropy or even the back-scattering of smaller grains (Kimura et al. 2003; Min et al. 2010), and the asymmetry would be less pronounced.

4.9. Disk mass and long-term evolution

The lifetimes of all the modelled disks are rather long. To illustrate this for the reference run, Fig. 20 shows the object lifetimes as a function of size, based on a weighted average over all distances. The comparison with the lifetimes for the dragless run reveals that transport becomes noticeable for grains smaller than 1 mm. For grains smaller than about 1 μm, transport even becomes the dominant loss channel, acting on timescales much shorter than the purely collisional timescales. For these grains, an increase in the mass throughput of the collisional cascade, as in runs e1 and e2, has no effect. However, for larger, collision-dominated objects, the lifetimes are consistently shorter.

Apart from a lower limit, no strong constraints can be put on the total disk mass. From a comparison of the collisional lifetimes with the estimated system age of 2 Gyr, we have found that the majority of objects larger than about 100 m might still be primordial, i.e. not currently distributed as in collisional equilibrium. Up to s = 100 m, the reference run contains 0.1 Earth masses. Above that size, our approach cannot shed any additional light because collisions barely occur. The mass depends on the assumed initial size distribution. For collisional equilibrium, κ =  −3.0, we would have M(s < 100   km) = 100 Earth masses, a value that is in-between the ones inferred for, e.g., Vega (Müller et al. 2010) and q1 Eridani (Augereau et al., in prep.). For κ =  −4.0, there would be only M(s < 100   km) = 0.7 Earth masses.

Assuming that the random velocities evKepler in the disk are comparable with the larger planetesimals’ escape velocities vesc = s(8πρG/3)1/2, their sizes s would need to be at least on the order of 100 km. These largest objects would then not suffer from destructive collisions at all. On the contrary, they would slowly grow further and their size distribution would depart from that of an equilibrium of destructive collisions.

5. Discussion

Several aspects of the current state and the history of HD 207129 are beyond the scope of this paper. We only discuss some of them briefly in this section.

5.1. Dust production and transport

With the power-law and collisional approach consistently predicting that: (1) grains of about 5 μm in radius dominate the cross-section; (2) there is a rather sharp outer edge; and (3) significant flux is emitted from within a distance of 90 AU from the star, it is very likely that transport by means of Poynting-Robertson and stellar wind drag plays a notable role in the dust dynamics around HD 207129. The low temperature seen in the SED in combination with the resolved location of the dust require the typical grains to be significantly larger than the blowout radius. A dynamically cold disk can provide a very plausible explanation of that situation. In turn, the resulting collisional timescales are sufficiently long to allow enough grains to be dragged inward and produce excess emission in the inner regions.

Nevertheless, the efficiency of transport through the whole disk is probably limited because inward drag would produce too much warm emission in the inner regions. As discussed in Sects. 4.4 and 4.6, this preferentially increases the excess flux at 70 μm relative to that at 100 and 160 μm. The central peak at 70 μm becomes too strong or the depression at 160 μm too deep. In consequence, drag needs to be balanced by collisions and/or another means of grain removal such as gravitational scattering by a planet.

5.2. Disk evolution

In addition to the uncertainties in the derived dust distribution, the precise formation history of the disk remains unclear. One possible scenario is that an extended disk gradually depletes inside-out by means of the formation of terrestrial planets and self-stirring. In their equation (28), Kenyon & Bromley (2008) report timescales t1000 = 475 Myr to grow an object of 1000 km in radius at a distance of 80 AU from a solar-mass star in a disk with surface density Σ = Σ0(r/r0)α, where α =  −3/2 and Σ0 = 30 g cm-2 at a distance r0 = 1 AU. At a distance of 160 AU, the distance that is relevant to HD 207129, this corresponds to t1000 = 3.8 Gyr and Σ = 1.5 × 10-2 g cm-2. For comparison, for 100 Earth masses spread from 130 to 175 AU around HD 207129 with α = 3.5 in run x, we would find that Σ = 7 × 10-2 g cm-2 at 160 AU, which is about five times as high as the value used by Kenyon & Bromley (2008). Consequentially, the timescales would be a factor of five shorter, i.e. tiso ≈ 800 Myr. This disk could have stirred itself sufficiently to the current age of HD 207129.

The depletion of material in the inner region can then be attributed to the strong dependence of collisional timescales on the radial distance. During the inside-out stirring, different annular regions in a disk brighten up and enter their erosive debris disk phase. Their individual fractional luminosities decay more slowly than 1/t (Dominik & Decin 2003; Wyatt et al. 2007; Trilling et al. 2008; Löhne et al. 2008). However, this decay is even more rapid than the outward migration of the stirring front. Since the maximum fractional luminosity, fmax, in Eq. (20) of Wyatt et al. (2007) is proportional to r7/3 for any given system age, the brightest region in a stirred disk will always be located at its outer edge, i.e. close to the stirring front.

This process is potentially also capable of producing the sharp outer edge, which is observed for HD 207129. Since the Hill radius of a Pluto-sized object at 160 AU is about one quarter of an AU, the range of influence of such an object is limited – as is the width of the stirring front.

In an alternative scenario, the disk would be less massive, have failed to grow oligarchs, and be stirred by either an internal or external planetary perturber or a stellar fly-by. This last possibility would require that HD 207129 has had a stellar encounter at some point in its history. The perturber would have needed to approach this star to within three times the disk radius, i.e. about 300–500 AU. In an open cluster with 10 % of the density of Trapezium, these events are likely to occur after about 5 Myr, well within the lifetime of such a cluster (Kobayashi et al. 2005). Indeed, there is a small chance that HD 207129 has recently passed through either the β Pictoris or AB Doradus group or the Tucana-Horologium association (Tetzlaff et al. 2011). It needs to be pointed out, though, that the probability of that scenario is rather low. In addition, this scenario is not exclusive: even if a star can perturb the outer regions of a disk first, a system as dust-rich at 160 AU as HD 207129 is likely to have formed a (potentially second) ring of debris on its own in any case – through self-stirring and/or an inner planet.

Stirring by planets faces the same challenge that it must also produce the narrow belt observed with the HST. While the outer regions are unstirred or void from the start, the inner regions would need to be depleted through mutual collisions or through direct interaction with the planet. The former case could work for secular perturbations from a massive inner planet (Moro-Martín et al. 2007; Mustill & Wyatt 2009) and the latter case would apply to resonant perturbations (Wyatt 2003; Quillen et al. 2007).

When comparing the disk around HD 207129 with the Edgeworth-Kuiper belt (EKB), leaving aside the sheer spatial scale, the dynamical excitation is the potentially biggest difference. Owing to the proposed violent history of the outer regions of our own system (e.g., Gomes et al. 2005), the typical orbital eccentricities in the EKB are rather high, with the bulk of objects having eccentricities around 0.1. This results in typical random velocities of vrel ≈ 500 m s-1. In contrast, for the HD 207129 belt, we found vrel < 120 m s-1 as a result of its being four times as distant with maximum eccentricities well below 0.1.

As far as the collisional lifetime of such a cold disk is concerned, HD 207129 does not necessarily need to harbour a “standard” disk that contains everything from dust to Plutos. Because of the low dynamical excitation, a size range with an upper end at about 100 m would suffice to sustain the dust production over the system’s age. This is consistent with the results of Heng & Tremaine (2010, their Figure 5, Disk F), who give  ≈ 100 m as the minimum size for the largest objects in what they call a “hot” disk at 100 AU that produces debris dust for 3 Gyr.

5.3. Asymmetries in the disk

The brightness asymmetries described in Sect. 2.3 might be based on real structures in the disk. Such azimuthal structures could, for example, be imposed by an inner planet in orbital resonance (e.g., Wyatt 2003). The orchestrated motion of objects trapped in resonance would result in denser regions, observable at far-infrared and submillimetre wavelengths. To produce clumps at about 150 AU, a planet would need to orbit the star, e.g., in a 2:1 resonance at 94 AU or in a 3:2 resonance at 114 AU. However, assuming that the planet trapped the planetesimals during its outward migration, the resulting dynamical excitation of the disk would be significantly higher than the values found here (emax < 0.1). This scenario is therefore inconsistent with our modelling results.

Pericentre glow (e.g., Wyatt et al. 1999) in a slightly off-set disk might be involved to explain the overly bright (by 4σ) north-western ansa. In such a slightly eccentric disk, the material is warmer when close to the pericentre. The according increase in brightness is, however, fully countered by the reduced density there that results from the higher (angular) orbital velocity. In addition, Krist et al. (2010) found the maximum offset that is still consistent with the HST observations is 0.4′′. For a disk radius of 10′′, the corresponding maximum eccentricity of 0.04 does not support this scenario.

6. Conclusions

Our study of HD 207129 has provided another good example of how complex the physics in a debris disk can be. The high angular resolution of the Herschel/PACS images and the well-studied spectral energy distribution have proven to be invaluable input to modelling and, in return, can provide some detailed insights into the chemical composition as well as the transport mechanisms and collisional processes at work.

From deconvolution, power-law fitting, and collisional modelling, we have found that most of the material is concentrated near the outer edge of the underlying distribution of unseen parent bodies. This outer edge is located around or beyond the inclined ring that has been observed at around 163 AU with the HST (Krist et al. 2010). In contrast to the widely spread haloes of small grains often observed around debris rings, the HD 207129 disk features a rather sharp outer edge. Only at distances beyond 400 AU do we expect it to follow the classical slope r-1.5 in terms of optical thickness. Significant amounts of dust are found to be present within the inner edge of the ring, though with a profile that steeply rises outward. If it exists at all, the inner cutoff of the dust distribution is likely located within 60 AU from the central star. This is in contrast to the sharp inner edge that is observed for Fomalhaut and presumably cleared by Fomalhaut b (Kalas et al. 2008).

In combination with the low temperatures and the accordingly large effective grain sizes, the inferred dynamical excitation of the disk is low and transport by P-R drag plays a notable role. At 4–8 μm, the radii of the grains that dominate the cross-section are by about a factor of ten above the blowout limit – and by about a factor of three above what is predicted for collisionally active disks. We therefore deduce that the rate of collisions is lower in accordance with the typical eccentricities emax < 0.1. The resulting depletion of the halo of barely bound grains agrees with the observations. While being dominated by collisions, large grains also get dragged inward and notably fill the inner gap in surface brightness, especially at 70 μm.

The overall effect of the chemical composition on the observables is small. We have tested different mixing ratios of amorphous silicate and water ice to verify that it is possible to produce good fits with both high silicate fractions and high ice fractions. The greatest differences occur at short wavelengths, where the Spitzer/IRS data are most consistent with moderate mixing ratios.

With its low dynamical excitation and its brightness peaking at an unprecedented radial distance of about 160 AU, HD 207129 might represent the limit to the population of coldest and most extended circumstellar debris disks detectable with Herschel. Regarding the potential role of planets in the formation of this disk, no stringent conclusions can be drawn.

Acknowledgments

Parts of our work were funded by the Deutsche Forschungsgemeinschaft (grants Lo 1715/1-1, Kr 2164/9-1, and Wo 857/7-1). C. Eiroa, J. Maldonado, J. P. Marshall, and B. Montesinos are partly supported by Spanish grant AYA 2008-01727. A. Mora acknowledges support from the Portuguese Fundação para a Ciência e a Tecnologia (grant PEst-OE/EEI/UI0066/2011). We thank Harald Mutschke for a discussion about refractive indices.

References

  1. Absil, O., di Folco, E., Mérand, A., et al. 2006, A&A, 452, 237 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Augereau, J.-C., Lagrange, A. M., Mouillet, D., Papaloizou, J. C. B., & Grorod, P. A. 1999, A&A, 348, 557 [NASA ADS] [Google Scholar]
  3. Beichman, C. A., Bryden, G., Rieke, G. H., et al. 2005, ApJ, 622, 1160 [NASA ADS] [CrossRef] [Google Scholar]
  4. Benz, W., & Asphaug, E. 1999, Icarus, 142, 5 [NASA ADS] [CrossRef] [Google Scholar]
  5. Bohren, C. F., & Huffman, D. R. 1983 (New York: Wiley) [Google Scholar]
  6. Brown, M. E. 2001, AJ, 121, 2804 [Google Scholar]
  7. Burns, J. A., Lamy, P. L., & Soter, S. 1979, Icarus, 40, 1 [NASA ADS] [CrossRef] [Google Scholar]
  8. Dohnanyi, J. S. 1969, J. Geophys. Res., 74, 2531 [NASA ADS] [CrossRef] [Google Scholar]
  9. Dominik, C., & Decin, G. 2003, ApJ, 598, 626 [NASA ADS] [CrossRef] [Google Scholar]
  10. Draine, B. T. 2003, ApJ, 598, 1017 [NASA ADS] [CrossRef] [Google Scholar]
  11. Eiroa, C., Fedele, D., Maldonado, J., et al. 2010, A&A, 518, L131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Eiroa, C., Marshall, J. P., Mora, A., et al. 2011, A&A, 536, L4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Elliot, J. L., Kern, S. D., Clancy, K. B., et al. 2005, AJ, 129, 1117 [NASA ADS] [CrossRef] [Google Scholar]
  14. Ertel, S., Wolf, S., Marshall, J. P., et al. 2011, A&A, submitted [Google Scholar]
  15. Golimowski, D. A., Ardila, D. R., Krist, J. E., et al. 2006, AJ, 131, 3109 [NASA ADS] [CrossRef] [Google Scholar]
  16. Gomes, R. S. 2003, Icarus, 161, 404 [Google Scholar]
  17. Gomes, R., Levison, H. F., Tsiganis, K., & Morbidelli, A. 2005, Nature, 435, 466 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  18. Griffin, M. J., Abergel, A., Abreu, A., et al. 2010, A&A, 518, L3 [Google Scholar]
  19. Grigorieva, A., Thébault, P., Artymowicz, P., & Brandeker, A. 2007, A&A, 475, 755 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Hauschildt, P. H., Allard, F., & Baron, E. 1999, ApJ, 512, 377 [NASA ADS] [CrossRef] [Google Scholar]
  21. Heng, K., & Tremaine, S. 2010, MNRAS, 401, 867 [NASA ADS] [CrossRef] [Google Scholar]
  22. Hudgins, D. M., Sandford, S. A., Allamandola, L. J., & Tielens, A. G. G. M. 1993, ApJS, 86, 713 [NASA ADS] [CrossRef] [Google Scholar]
  23. Ishihara, D., Onaka, T., Kataza, H., et al. 2010, A&A, 514, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Jourdain de Muizon, M., Laureijs, R. J., Dominik, C., et al. 1999, A&A, 350, 875 [NASA ADS] [Google Scholar]
  25. Kalas, P., Graham, J. R., Chiang, E., et al. 2008, Science, 322, 1345 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  26. Kalas, P., Graham, J. R., & Clampin, M. 2005, Nature, 435, 1067 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  27. Kimura, H., Mann, I., Jessberger, E. K., & Weber, I. 2003, Meteorit. Planet. Sci. Suppl., 38, 5211 [NASA ADS] [Google Scholar]
  28. Kenyon, S. J., & Bromley, B. C. 2008, ApJS, 179, 451 [NASA ADS] [CrossRef] [Google Scholar]
  29. Kobayashi, H., Ida, S., & Tanaka, H. 2005, Icarus, 177, 246 [NASA ADS] [CrossRef] [Google Scholar]
  30. Kolmogorov, A. N. 1941, Bull. Math. Mosc. Univ., 2, 1 [Google Scholar]
  31. Krist, J. E., Stapelfeldt, K. R., Bryden, G., et al. 2010, AJ, 140, 1051 [NASA ADS] [CrossRef] [Google Scholar]
  32. Krivov, A. V., Sremčević, M., & Spahn, F. 2005, Icarus, 174, 105 [NASA ADS] [CrossRef] [Google Scholar]
  33. Krivov, A. V., Löhne, T., & Sremčević, M. 2006, A&A, 455, 509 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Krivov, A. V., Müller, S., Löhne, T., & Mutschke, H. 2008, ApJ, 687, 608 [NASA ADS] [CrossRef] [Google Scholar]
  35. Krivov, A. V., Herrmann, F., Brandeker, A., & Thébault, P. 2009, A&A, 507, 1503 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Li, A., & Greenberg, J. M. 1998, A&A, 331, 291 [NASA ADS] [Google Scholar]
  37. Löhne, T., Krivov, A. V., & Rodmann, J. 2008, ApJ, 673, 1123 [NASA ADS] [CrossRef] [Google Scholar]
  38. Lucy, L. B. 1974, AJ, 79, 745 [NASA ADS] [CrossRef] [Google Scholar]
  39. Magain, P., Courbin, F., & Sohy, S. 1998, ApJ, 494, 472 [NASA ADS] [CrossRef] [Google Scholar]
  40. Marshall, J. P., Löhne, T., Montesinos, B., et al. 2011, A&A, 529, A117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Min, M., Kama, M., Dominik, C., & Waters, L. B. F. M. 2010, A&A, 509, L6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Moro-Martín, A., Malhotra, R., Carpenter, J. M., et al. 2007, ApJ, 668, 1165 [NASA ADS] [CrossRef] [Google Scholar]
  43. Müller, S., Löhne, T., & Krivov, A. V. 2010, ApJ, 708, 1728 [NASA ADS] [CrossRef] [Google Scholar]
  44. Mustill, A. J., & Wyatt, M. C. 2009, MNRAS, 399, 1403 [NASA ADS] [CrossRef] [Google Scholar]
  45. Nilsson, R., Liseau, R., Brandeker, A., et al. 2010, A&A, 518, A40 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. O’Brien, D. P., & Greenberg, R. 2003, Icarus, 164, 334 [CrossRef] [Google Scholar]
  47. Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1 [CrossRef] [EDP Sciences] [Google Scholar]
  48. Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 1992, Numerical recipes in FORTRAN. The art of scientific computing, 2nd ed. (Cambridge: University Press) [Google Scholar]
  50. Quillen, A. C., Morbidelli, A., & Moore, A. 2007, MNRAS, 380, 1642 [NASA ADS] [CrossRef] [Google Scholar]
  51. Rhee, J. H., Song, I., Zuckerman, B., & McElwain, M. 2007, ApJ, 660, 1556 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  52. Richardson, W. H. 1972, Journal of the Optical Society of America (1917–1983), 62, 55 [Google Scholar]
  53. Schütz, O., Meeus, G., & Sterzik, M. F. 2005, A&A, 431, 175 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Sheret, I., Dent, W. R. F., & Wyatt, M. C. 2004, MNRAS, 348, 1282 [NASA ADS] [CrossRef] [Google Scholar]
  55. Sibthorpe, B., Vandenbussche, B., Greaves, J. S., et al. 2010, A&A, 518, L130 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Stewart, S. T., & Leinhardt, Z. M. 2009, ApJ, 691, L133 [NASA ADS] [CrossRef] [Google Scholar]
  57. Strubbe, L. E., & Chiang, E. I. 2006, ApJ, 648, 652 [NASA ADS] [CrossRef] [Google Scholar]
  58. Su, K. Y. L., Rieke, G. H., Misselt, K. A., et al. 2005, ApJ, 628, 487 [NASA ADS] [CrossRef] [Google Scholar]
  59. Tanner, A., Beichman, C., Bryden, G., Lisse, C., & Lawler, S. 2009, ApJ, 704, 109 [NASA ADS] [CrossRef] [Google Scholar]
  60. Tetzlaff, N., Neuhäuser, R., & Hohle, M. M. 2011, MNRAS, 410, 190 [NASA ADS] [CrossRef] [Google Scholar]
  61. Thébault, P., & Augereau, J.-C. 2007, A&A, 472, 169 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Thébault, P., & Wu, Y. 2008, A&A, 481, 713 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  63. Thébault, P., Augereau, J.-C., & Beust, H. 2003, A&A, 408, 775 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  64. Trilling, D. E., Bryden, G., Beichman, C. A., et al. 2008, ApJ, 674, 1086 [NASA ADS] [CrossRef] [Google Scholar]
  65. Vandenbussche, B., Sibthorpe, B., Acke, B., et al. 2010, A&A, 518, L133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Vitense, C., Krivov, A. V., & Löhne, T. 2010, A&A, 520, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Walker, H. J., & Wolstencroft, R. D. 1988, PASP, 100, 1509 [NASA ADS] [CrossRef] [Google Scholar]
  68. Warren, S. G. 1984, Appl. Opt., 23, 1206 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  69. Wiener, N. 1949, Extrapolation, Interpolation, and Smoothing of Stationary Time Series (New York: Wiley) [Google Scholar]
  70. Wyatt, M. C. 2003, ApJ, 598, 1321 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  71. Wyatt, M. C. 2005, A&A, 433, 1007 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Wyatt, M. C. 2008, ARA&A, 46, 339 [NASA ADS] [CrossRef] [Google Scholar]
  73. Wyatt, M. C., Dermott, S. F., Telesco, C. M., et al. 1999, ApJ, 527, 918 [NASA ADS] [CrossRef] [Google Scholar]
  74. Wyatt, M. C., Smith, R., Greaves, J. S., et al. 2007, ApJ, 658, 569 [NASA ADS] [CrossRef] [Google Scholar]

All Tables

Table 1

Photometry of HD 207129.

Table 2

Results of power-law fitting.

Table 3

Description of numerical runs.

All Figures

thumbnail Fig. 1

Herschel/PACS images of HD 207129 at (top to bottom) 70, 100, and 160 μm. Columns a) and b) show the observed and star-subtracted images, respectively. Columns c) to f) show deconvolved images using the following methods: c) van Cittert, d) FFT with a Wiener filter, e) MCS code, f) Richardson-Lucy. Flux steps (in mJy/arcsec2) between solid isolines are indicated in the top-right corner of each panel. The photocentre of the disk, i.e. the assumed stellar position, is indicated by a cross. In every panel, North is up and East is left.

In the text
thumbnail Fig. 2

Radial profiles of surface brightness at (top to bottom) 70, 100, and 160 μm: (left) along the long axis and (right) along the short axis. The coloured shaded areas mark the uncertainty, while symbols represent the two respective directions along each axis. In the lower-left of each panel, the stellar contribution is plotted. Azimuthally averaged radial profiles of Wiener deconvolved and deprojected images are given as solid line, and their azimuthal standard deviations are represented as dotted lines, both in arbitrary units.

In the text
thumbnail Fig. 3

Comparison of (convolved and unconvolved) rectangular surface brightness profiles with the observations and their deconvolutions. Wavelengths (top to bottom): 70, 100, and 160 μm. The dashed and dotted lines are given in arbitrary units.

In the text
thumbnail Fig. 4

Imaginary part of the refractive indices in the far infrared for (solid line) pure amorphous water ice as compiled by Li & Greenberg (1998), (symbols) water ice at different temperatures from Hudgins et al. (1993, the original data source for Li & Greenberg), (dashed line) crystalline ice as compiled by Warren (1984), and (dotted line) astronomical silicate synthesised by Draine (2003).

In the text
thumbnail Fig. 5

Spectral energy distributions of the best-fit power law-based models. The photometry points used for the fitting are overplotted. Note that the horizontal stretch is different for wavelengths smaller than 40 μm. To compare the models to the IRS data from the standard pipeline, we corrected the model SEDs for slit loss and then applied a point-source correction. The uncorrected model SEDs are plotted with light grey lines.

In the text
thumbnail Fig. 6

Radial brightness profiles of the best-fit power law model: (top) 70 μm, (middle) 100 μm, (bottom) 160 μm. Cuts along the long and short axes are given in the left and right columns, respectively.

In the text
thumbnail Fig. 7

The β ratio as a function of grain size for different mixing ratios of ice and silicate, based on Mie theory (see, e.g., Bohren & Huffman 1983; Burns et al. 1979).

In the text
thumbnail Fig. 8

Radial brightness profiles of the reference run.

In the text
thumbnail Fig. 9

Same as Fig. 5 but for the spectral energy distribution of the reference run and the contributions from its constituent annuli and the star. For comparison, the dragless case (run d) is dashed.

In the text
thumbnail Fig. 10

Radial brightness profiles of the reference run. Shaded regions mark the contributions from the star and from constituent annuli of the disk.

In the text
thumbnail Fig. 11

(Left) Radial profiles of optical thickness for two runs: reference and iceless. For comparison, the light shaded region shows the profile of a model of the Vega disk (Müller et al. 2010), radially scaled so that the outer edges match. The dark shaded region marks the extent of the belt of parent bodies. (Right) Size distribution of the reference run at distances of 50, 100, 150, and 200 AU.

In the text
thumbnail Fig. 12

Optical thickness as a function of radial distance and grain size for the reference run: (left) 3–23 μm, (right) 0.4–1.8 μm.

In the text
thumbnail Fig. 13

Colour-coded distribution of geometrical optical depth, i.e. particle cross-section, over grain size and distance to the star given for the reference run. The white solid line gives the peak position in the size distribution as a function of radial distance. The horizontal dotted line indicates the blowout size, where β = 0.5.

In the text
thumbnail Fig. 14

Size distribution in terms of optical depth at 150 AU as a function of grain size. The reference run is compared with runs (left) Q1 and Q2, as well as (right) e1, e2, and d.

In the text
thumbnail Fig. 15

Same as Fig. 10 but for a comparison of the (solid line) reference run with runs Q1, Q2, e1, e2, and d. See Table 3 for a description of these runs.

In the text
thumbnail Fig. 16

Spectral energy distribution of the reference and runs (top) Q1 and Q2 as well as (bottom) e1 and e2.

In the text
thumbnail Fig. 17

Evolution in time of run c1: (left) the SED and (right) the radial brightness profile at 70 μm.

In the text
thumbnail Fig. 18

Spectral energy distribution of the reference and runs (top) c1 and c2 as well as (bottom) a1 and x.

In the text
thumbnail Fig. 19

Same as Fig. 10 but for a comparison of the (solid line) reference run with runs c1, c2, a1, and x. See Table 3 for a description of these runs.

In the text
thumbnail Fig. 20

(Grey solid line) UV-sputtering lifetimes compared with collisional lifetimes of objects in (black solid line) the reference run, (dashed line) run d, and (dotted line) run e1.

In the text
thumbnail Fig. 21

Synthetic scattered light images of three of the runs: (left) the reference, (middle) Q1, and (right) x. Note that run Q1 is artificially truncated at a radial distance of 400 AU (25′′). All three disks are inclined by 60 degrees with respect to the sky plane, the north-eastern part being closer to the observer.

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.