# A peculiar class of debris disks from Herschel/DUNES. A steep fall off in the far infrared

S. Ertel, S. Wolf, J. P. Marshall, C. Eiroa, J.-C. Augereau, A. V. Krivov, T. Löhne, O. Absil, D. Ardila, M. Arévalo, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

S. Ertel, S. Wolf, J. P. Marshall, C. Eiroa, J.-C. Augereau, et al. A peculiar class of debris disks from Herschel/DUNES. A steep fall off in the far infrared. Astronomy and Astrophysics - A&A, 2012, 541, 10.1051/0004-6361/201118077. insu-03612435

HAL Id: insu-03612435

https://insu.hal.science/insu-03612435

Submitted on 18 Mar 2022

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A peculiar class of debris disks from Herschel/DUNES*

## A steep fall off in the far infrared

S. Ertel<fnref n="1" />,<fnref n="2" />, S. Wolf<fnref n="2" />, J. P. Marshall<fnref n="3" />, C. Eiroa<fnref n="3" />, J.-C. Augereau<fnref n="1" />, A. V. Krivov<fnref n="4" />, T. Löhne<fnref n="4" />, O. Absil<fnref n="5" />, D. Ardila<fnref n="6" />,
M. Arévalo<fnref n="7" />, A. Bayo<fnref n="8" />, G. Bryden<fnref n="9" />, C. del Burgo<fnref n="10" />, J. Greaves<fnref n="11" />, G. Kennedy<fnref n="12" />, J. Lebreton<fnref n="1" />, R. Liseau<fnref n="13" />,
J. Maldonado<fnref n="3" />, B. Montesinos<fnref n="7" />, A. Mora<fnref n="14" />, G. L. Pilbratt<fnref n="15" />, J. Sanz-Forcada<fnref n="7" />, K. Stapelfeldt<fnref n="16" />, and G. J. White<fnref n="17" />,<fnref n="18" />

Received 12 September 2011 / Accepted 19 March 2012

## ABSTRACT

**Context.** The existence of debris disks around old main sequence stars is usually explained by continuous replenishment of small dust grains through collisions from a reservoir of larger objects.

**Aims.** We present photometric data of debris disks around HIP 103389 (HD 199260), HIP 107350 (HN Peg, HD 206860), and HIP 114948 (HD 219482), obtained in the context of our *Herschel* open time key program DUNES (DUst around NEarby Stars).

**Methods.** We used Herschel/PACS to detect the thermal emission of the three debris disks with a $3\sigma$ sensitivity of a few mJy at 100 µm and 160 µm. In addition, we obtained Herschel/PACS photometric data at 70 µm for HIP 103389. These observations are complemented by a large variety of optical to far-infrared photometric data. Two different approaches are applied to reduce the Herschel data to investigate the impact of data reduction on the photometry. We fit analytical models to the available spectral energy distribution (SED) data using the fitting method of simulated thermal annealing as well as a classical grid search method.

**Results.** The SEDs of the three disks potentially exhibit an unusually steep decrease at wavelengths $\ge 70$ µm. We investigate the significance of the peculiar shape of these SEDs and the impact on models of the disks provided it is real. Using grain compositions that have been applied successfully for modeling of many other debris disks, our modeling reveals that such a steep decrease of the SEDs in the long wavelength regime is inconsistent with a power-law exponent of the grain size distribution –3.5 expected from a standard equilibrium collisional cascade. In contrast, a steep grain size distribution or, alternatively an upper grain size in the range of few tens of micrometers are implied. This suggests that a very distinct range of grain sizes would dominate the thermal emission of such disks. However, we demonstrate that the understanding of the data of faint sources obtained with Herschel is still incomplete and that the significance of our results depends on the version of the data reduction pipeline used.

**Conclusions.** A new mechanism to produce the dust in the presented debris disks, deviations from the conditions required for a standard equilibrium collisional cascade (grain size exponent of –3.5), and/or significantly different dust properties would be necessary to explain the potentially steep SED shape of the three debris disks presented.

**Key words.** circumstellar matter – stars: individual: HIP 103389 – infrared: planetary systems – stars: individual: HIP 107350 – infrared: stars – stars: individual: HIP 114948

## 1. Introduction

Debris disks were first discovered via infrared excess emission associated with main sequence stars such as Vega detected

<footnote n="1">¹ UJF-Grenoble 1 / CNRS-INSU, Institut de Planétologie et d'Astrophysique de Grenoble (IPAG) UMR 5274, 38041 Grenoble, France
e-mail: steve.ertel@obs.ujf-grenoble.fr</footnote>

<footnote n="2">² Institut für Theoretische Physik und Astrophysik, Christian-Albrechts-Universität zu Kiel, Leibnizstraße 15, 24098 Kiel, Germany</footnote>

<footnote n="3">³ Dpt. Física Teórica, Facultad de Ciencias, Universidad Autónoma de Madrid, Cantoblanco, 28049 Madrid, Spain</footnote>

<footnote n="4">⁴ Astrophysikalisches Institut und Universitätssternwarte, Friedrich-Schiller-Universität, Schillergäßchen 2-3, 07745 Jena, Germany</footnote>

<footnote n="5">⁵ Institut d'Astrophysique et de Géophysique, Université de Liège, Allée du Six Août 17, 4000 Sart Tilman, Belgium</footnote>

<footnote n="6">⁶ NASA Herschel Science Center, California Institute of Technology, 1200 E. California Blvd., Pasadena, CA 91125, USA</footnote>

<footnote n="7">⁷ Department of Astrophysics, Centre for Astrobiology (CAB, CSIC-INTA), ESAC Campus, 28691 PO Box 78, Villanueva de la Cañada, Madrid, Spain</footnote>

<footnote n="8">⁸ European Space Observatory, Alonso de Cordova 3107, Vitacura Casilla 19001, Santiago 19, Chile</footnote>

<footnote n="9">⁹ Jet Propulsion Laboratory, California Institute of Technology, Pasadena, CA 91109, USA</footnote>

<footnote n="10">¹⁰ UNINOVA-CA3, Campus da Caparica, Quinta da Torre, Monte de Caparica, 2825-149 Caparica, Portugal</footnote>

<footnote n="11">¹¹ School of Physics and Astronomy, St Andrews University, North Haugh, St Andrews, Fife KY16 9SS, UK</footnote>

<footnote n="12">¹² Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK</footnote>

<footnote n="13">¹³ Onsala Space Observatory, Chalmers University of Technology, 439 92 Onsala, Sweden</footnote>

<footnote n="14">¹⁴ ESA-ESAC Gaia SOC, PO Box 78, 28691 Villanueva de la Cañada, Madrid, Spain</footnote>

<footnote n="15">¹⁵ ESA Astrophysics & Fundamental Physics Missions Division, ESTEC/SRE-SA, Keplerlaan 1, 2201 AZ Noordwijk, The Netherlands</footnote>

<footnote n="16">¹⁶ Code 667, NASA Goddard Space Flight Center, Greenbelt MD 20771, USA</footnote>

<footnote n="17">¹⁷ Department of Physics and Astrophysics, Open University, Walton Hall, Milton Keynes MK7 6AA, UK</footnote>

<footnote n="18">¹⁸ Rutherford Appleton Laboratory, Chilton OX11 0QX, UK</footnote>

<footnote>* Herschel is an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA.</footnote>
by the InfraRed Astronomical Satellite (IRAS; Aumann et al. 1984). The first spatially resolved image of a debris disk was that of β Pictoris (Smith & Terrile 1984) in optical scattered light. In the last few years, the *Spitzer* Space Telescope has revealed that debris disks are common around main sequence stars (e.g., Trilling et al. 2008). Several key programs on the Herschel Space Observatory (Pilbratt et al. 2010) are dedicated to the study of various aspects of the formation and evolution of planetary systems and their attendant circumstellar debris disks (e.g., Augereau et al. 2008). Our Herschel open time key program (OTKP) DUNES (DUst around NEarby Stars; Eiroa et al. 2010) aims to detect debris disks with fractional luminosities similar to the Edgeworth-Kuiper Belt level ($L_d/L_* = 10^{-7}$ to $10^{-6}$; Stern 1996; Vitense et al. 2010) around a volume limited sample ($d \le 20$ pc) of Sun-like stars (F, G, and K spectral type). Some additional sources at $20 \text{ pc} < d \le 25$ pc are included because of their known excesses that have previously been detected with *Spitzer*, or because they are known exoplanet host stars.

As the most readily detectable signposts of other planetary systems, debris disks help us to improve our understanding of the formation and evolution of them as well as of our own solar system (Meyer et al. 2007; Wyatt 2008; Krivov 2010). Studying the spectral energy distribution (SED) of the dust alone usually provides only weak, ambiguous constraints to their properties such as chemical composition, grain size, and spatial distribution (Wolf & Hillenbrand 2003). For example, the location of the inner disk radius of the dust distribution is strongly degenerate with the lower limit of the grain size distribution. However, the grain size distribution can be described by a power-law and the power-law exponent derived from pure SED fitting is usually well constrained and consistent with the analytical value of -3.5 derived by Dohnanyi (1969) under several assumptions: (1) the particles are produced through a collisional cascade from infinitely large to extremely small grains; (2) the grains are not affected by any other influence such as stellar radiation and drag forces; (3) the strength of the particles, i.e., the energy per volume fraction necessary to disrupt them, is independent from the grain size. This is further referred to as a standard equilibrium collisional cascade.

In this paper, the potential Herschel/DUNES discovery of an unusually steep decrease of the SEDs of three spatially unresolved debris disks around Sun-like stars, HIP 103389 (HD 199260), HIP 107350 (HN Peg, HD 206860), and HIP 114948 (HD 219482), is presented<fnref n="1" />. It is demonstrated on the example of the three targets that the analysis of data of faint point sources obtained with Herschel depends very much on the photometric calibration and the exact determination of the uncertainties. In particular, the significance of our results depends on the version of the data reduction pipeline used. Provided it is real, the steep decrease occurs in the range of 70 µm to 160 µm, inconsistent with a dust grain size distribution following a power-law derived from a standard equilibrium collisional cascade. The 100 µm and 160 µm data represent the first photometric measurements of these faint disks at wavelengths >70 µm and, thus, are the first ones to reveal the unusual shape of the SEDs. Results from detailed SED modeling including additional photometric data from the literature are presented as well.

Observations and data reduction are described in Sect. 2 together with the basic observational results. A theoretical

<table><caption>**Table 1.** Summary of Herschel/PACS observations of the three objects.</caption><thead><tr><td>Object</td><td>Observation ID</td><td>OD</td><td>λ [µm]</td><td>T [s]</td></tr></thead><tbody><tr><td>HIP 103389</td><td>1342193157/58</td><td>322</td><td>70/160</td><td>180</td></tr><tr><td>HIP 103389</td><td>1342193159/60</td><td>322</td><td>100/160</td><td>1440</td></tr><tr><td>HIP 107350</td><td>1342195779/80</td><td>355</td><td>100/160</td><td>1440</td></tr><tr><td>HIP 114948</td><td>1342196803/04</td><td>372</td><td>100/160</td><td>1440</td></tr></tbody></table>

**Notes.** OD is the observing date that counts the elapsed time in days since launch on 14 May 2009, λ gives the central wavelength of the two filters used simultaneously, and T is the total on-source integration time.

discussion on how steep the decrease of the SED of a debris disk toward longer wavelengths is expected to be and a characterization of unusually steep SED sources are given in Sect. 3. A detailed description of our modeling of the systems can be found in Sect. 4 and results are discussed in Sect. 5. Conclusions are given in Sect. 6.

## 2. Observations and data reduction

Two Herschel/PACS mini-scan map observations of each target were taken with the 100/160 channel combination at array orientation angles of 70° and 110° providing scan and cross-scan coverage to assist in the removal of noise artifacts from the final composite mosaic. In addition, two scan map observations of HIP 103389 were taken with the 70/160 channel combination with the same array orientations. Each scan map consists of 10 legs of 3' length, with a 4" separation between legs, at the medium slew speed (20" per second). In this way, a region of ≈1 square arc minute around the source position was covered to uniform depth in the resulting mosaic. A summary of the observing log is presented in Table 1.

PACS data reduction was carried out in version 4.2 (Reduction 1) and 7.2 (Reduction 2) of HIPE (Ott et al. 2010) starting from the level 0 products using a modified version of the standard reduction script (provided within HIPE). Two different reductions have been carried out, because changes in the still developing pipeline of HIPE may significantly affect the results. The separate scans at the two position angles of each channel pair were mosaicked to produce a final image at each wavelength. Output scales for the final mosaics are 1" per pixel for the 70 µm and 100 µm images and 2" per pixel for the 160 µm image (smaller than the native detector pixel sizes of 3'2 at 70 µm and 100 µm and 6'4 at 160 µm, respectively).

For the further analysis of the data, different approaches have been used for Reductions 1 and 2. These are described in the following:

**Reduction 1:** A high-pass filter was used to remove large scale background variation from the images, with filter widths of 31 frames at 70 µm and 100 µm and 51 frames at 160 µm (equivalent to ~60" and 100", respectively). A central region of 30" radius in the images was masked from the high pass filter process to prevent the removal of any faint, extended structure near the source position. PACS fluxes and sky noise were measured by aperture photometry using a custom script based on the IDL APER routine, based on the DAOPHOT aperture photometry module and cross-checked with the internal HIPE aperture photometry routines. The aperture radius and sky annulus dimensions were 20" and 30"-40", respectively at all bands. Sky noise for each wavelength was calculated from the rms variance

<footnote n="1">¹ Two of the sources, HIP 103389 and HIP 107350, are shared targets between the DUNES survey in the context of which the analysis is carried out and the DEBRIS survey (Matthews et al. 2010; Phillips et al. 2010).</footnote>
<table><caption>**Table 2.** Observational results obtained in this work<sup>a</sup>.</caption><thead><tr><th rowspan="2">Source</th><th colspan="2">HIP 103389</th><th colspan="2">HIP 107350</th><th colspan="2">HIP 114948</th></tr><tr><th>Reduction 1</th><th>Reduction 2</th><th>Reduction 1</th><th>Reduction 2</th><th>Reduction 1</th><th>Reduction 2</th></tr></thead><tbody><tr><td>F<sub>PACS 70</sub> [mJy]</td><td>47.4 ± 2.7</td><td>44.0 ± 2.3</td><td>...</td><td>...</td><td>...</td><td>...</td></tr><tr><td>F<sub>MIPS 70</sub> [mJy]</td><td>46.6 ± 3.8</td><td></td><td>28.4 ± 2.5</td><td></td><td>68.7 ± 3.0</td><td></td></tr><tr><td>F<sub>PACS 100</sub> [mJy]</td><td>23.7 ± 1.4</td><td>26.3 ± 1.7</td><td>11.0 ± 0.9</td><td>15.1 ± 1.3</td><td>42.5 ± 2.2</td><td>40.8 ± 1.6</td></tr><tr><td>F<sub>PACS 160</sub> [mJy]</td><td>5.0 ± 1.3</td><td>7.7 ± 2.5</td><td>4.4 ± 1.5</td><td>4.4 ± 2.3</td><td>12.7 ± 1.9</td><td>13.3 ± 2.2</td></tr><tr><td>F<sub>*</sub>,70µm<sup>b</sup> [mJy]</td><td></td><td>13.8</td><td></td><td>13.4</td><td></td><td>15.0</td></tr><tr><td>F<sub>*</sub>,100µm<sup>b</sup> [mJy]</td><td></td><td>6.8</td><td></td><td>6.6</td><td></td><td>7.3</td></tr><tr><td>F<sub>*</sub>,160µm<sup>b</sup> [mJy]</td><td></td><td>2.7</td><td></td><td>2.6</td><td></td><td>2.9</td></tr><tr><td>Δ<sub>100,70</sub></td><td>1.94 ± 0.32</td><td>1.44 ± 0.33</td><td>2.66 ± 0.45</td><td>1.77 ± 0.49</td><td>1.35 ± 0.26</td><td>1.46 ± 0.23</td></tr><tr><td>Δ<sub>160,100</sub></td><td>3.31 ± 0.69</td><td>2.61 ± 0.82</td><td>1.95 ± 0.88</td><td>2.62 ± 1.30</td><td>2.57 ± 0.44</td><td>2.38 ± 0.43</td></tr><tr><td>Δ<sub>160,70</sub></td><td>2.72 ± 0.39</td><td>2.11 ± 0.46</td><td>2.25 ± 0.53</td><td>2.26 ± 0.74</td><td>2.04 ± 0.22</td><td>1.99 ± 0.25</td></tr><tr><td>FWHM<sub>PACS 100</sub><sup>c</sup></td><td>6.8 × 7.1</td><td></td><td>7.1 × 6.0</td><td></td><td>7.2 × 7.2</td><td></td></tr><tr><td>RA<sub>PACS 100</sub></td><td>20<sup>h</sup>56<sup>m</sup>47.47</td><td></td><td>21<sup>h</sup>44<sup>m</sup>31.27</td><td></td><td>23<sup>h</sup>16<sup>m</sup>57.46</td><td></td></tr><tr><td>Dec<sub>PACS 100</sub></td><td>-26°17'45.5</td><td></td><td>+14°46'19.4</td><td></td><td>-62°00'05.7</td><td></td></tr><tr><td>RA<sub>optical</sub></td><td>20<sup>h</sup>56<sup>m</sup>47.33</td><td></td><td>21<sup>h</sup>44<sup>m</sup>31.33</td><td></td><td>23<sup>h</sup>16<sup>m</sup>57.69</td><td></td></tr><tr><td>Dec<sub>optical</sub></td><td>-26°17'47.0</td><td></td><td>+14°46'19.0</td><td></td><td>-62°00'04.3</td><td></td></tr></tbody></table>

**Notes.** Uncertainties of the flux measurements are total uncertainties including sky noise and calibration uncertainties. $\Delta v_1, v_2$ is the spectral index measured from $v_1$ and $v_2$ (Sect 3), identified with the corresponding PACS bands. <sup>(a)</sup> The *Spitzer*/MIPS 70 µm photometry is also listed, since these data complement the *Herschel* observations and are important to illustrate the unusual behavior of the SEDs. These data are published by Beichman et al. (2006, HIP 103389 and HIP 114948) and Bryden et al. (2006, HIP 107350) and have been re-reduced in the context of DUNES (Eiroa et al., in prep.). <sup>(b)</sup> Predicted stellar photosphere using a PHOENIX/GAIA synthetic stellar model (Eiroa et al., in prep.). <sup>(c)</sup> FWHM as measured from the PACS images at 100 µm using a 2D Gaussian fit.

of the sky annulus pixel values multiplied by the beam size. Results were scaled by the appropriate aperture correction. Absolute uncertainties are 5% at 70 µm and 100 µm and 160 µm using the calibration files provided along with this HIPE version.

The images of the three sources from Reduction 1 at 100 µm and 160 µm are shown in Fig. 1.

**Reduction 2:** The same filter width as in Reduction 1 has been used for the high-pass filtering, but sources were masked based on a threshold criterion, the value of which was calculated from the standard deviation of all non-zero pixels in the image. PACS fluxes and sky noise were measured by aperture photometry using both the internal HIPE aperture photometry routines and the MIDAS data analysis package to check for consistency. The source flux was measured within a circular aperture of 4", 5", and 8" radius at 70 µm, 100 µm, and 160 µm, respectively, and scaled by the appropriate aperture correction. Sky noise for each wavelength was calculated from the rms variance of ten sky apertures of the same size as the source aperture and randomly distributed across the uniformly covered part of the image. The resulting deviation of these apertures was calibrated by the aperture correction and correlated noise correction factors. Extensive testing has shown that this approach efficiently accounts for the correlated noise in the output images. In fact, this approach provides slightly conservative uncertainties. An extended discussion will be included in the survey overview paper (Eiroa et al., in prep.). Aperture correction factors have also been applied to the rms aperture flux dispersion. Photometric calibration uncertainties in this version of HIPE are 3–5%<fnref n="2" />.

## 2.1. Observational results

The data obtained in this work for the three sources using both reductions as well as the predicted photospheric fluxes are listed in Table 2. Table 3 lists the main stellar parameters

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Images of the three sources at 100 µm (top) and 160 µm (bottom) from Reduction 1 (Sect. 2). The images are displayed in a linear stretch from -0.1 $F_{peak}$ to $F_{peak}$ of the source of interest. The peak flux is 0.38 mJy/pixel, 0.23 mJy/pixel, and 0.54 mJy/pixel at 100 µm and 0.17 mJy/pixel, 0.16 mJy/pixel, and 0.45 mJy/pixel at 160 µm for HIP 103389, HIP 107350, and HIP 114948, respectively. The white circles in the 160 µm images are centered on the brightest pixel of the source, respectively.</figcaption></figure>

of HIP 103389, HIP 107350, and HIP 114948 (Eiroa et al., in prep.). The stellar contribution to the total SED of each star is estimated using a PHOENIX/GAIA synthetic stellar model (Brott & Hauschildt 2005). Specific models for $T_{eff}$, log g*, and [Fe/H] of each star are built by interpolation from the available grid. A normalization is applied to the short wavelength range ($\lambda < 20$ µm) of the long wavelength section of the *Spitzer*/IRS spectrum where available. Otherwise, the normalization is done to the flux measurements at photometric bands from B to K. All three sources are found to have significant excess (≥3 σ) at 70 µm and 100 µm. HIP 114948 also has significant excess at 160 µm. HIP 103389 and HIP 107350 have estimated excesses at 160 µm between 0.8 σ and 2.4 σ, depending on the reduction. Offsets between optical positions and positions in the PACS 70 µm images are below 1 σ of the *Herschel* pointing accuracy.

<footnote n="2">² Release note PICC-ME-TN-037.</footnote>
<table><caption>**Table 3.** Physical properties of the three stars considered in this work (for details and references see Sect. 4.1).</caption><thead><tr><th>Object</th><th>HIP 103389</th><th>HIP 107350</th><th>HIP 114948</th></tr></thead><tbody><tr><td>Name</td><td>HR 8013, HD 199260</td><td>HN Peg, HD 206860</td><td>HR 8843, HD 219482</td></tr><tr><td>Distance (pc)</td><td>21.97</td><td>17.88</td><td>20.54</td></tr><tr><td>Spectral type and luminosity class</td><td>F7 V</td><td>G0 V</td><td>F7 V</td></tr><tr><td>Range of spectral types</td><td>F6V–F8V</td><td>G0IV–V</td><td>F6V–F8V</td></tr><tr><td>V, B – V</td><td>5.70, 0.51</td><td>5.96, 0.59</td><td>5.64, 0.52</td></tr><tr><td>Absolute magnitude M<sub>v</sub>, bolometric correction</td><td>3.99, –0.02</td><td>4.70, –0.05</td><td>4.08, –0.02</td></tr><tr><td>Bolometric luminosity, L<sub>★</sub> [L<sub>☉</sub>]</td><td>2.025</td><td>1.090</td><td>1.867</td></tr><tr><td>Effective temperature (K)</td><td>6257</td><td>5952</td><td>6240</td></tr><tr><td>Surface gravity, log g<sub>★</sub></td><td>4.36</td><td>4.44</td><td>4.31</td></tr><tr><td>Radius, R<sub>★</sub> [R<sub>☉</sub>]</td><td>1.2</td><td>0.99</td><td>1.17</td></tr><tr><td>Metallicity [Fe/H]</td><td>–0.14</td><td>–0.07</td><td>–0.21</td></tr><tr><td>Radial velocity (km s<sup>–1</sup>)</td><td>–16.1</td><td>18.1</td><td>0.47</td></tr><tr><td>Rotational velocity, v sin i (km s<sup>–1</sup>)</td><td>13.7</td><td>12.8</td><td>9.0</td></tr><tr><td>Rotation period, P (days)</td><td>4.4 (P / sin i)</td><td>4.7</td><td>–</td></tr><tr><td>Li I equivalent width (mÅ)</td><td>80.0</td><td>91.9</td><td>73.7</td></tr><tr><td>Space velocities U, V, W (km s<sup>–1</sup>)</td><td>–17.2, –10.2, 1.1</td><td>–13.9, –20.1, –11.1</td><td>14.9, –8.0, –5.7</td></tr><tr><td>Activity, log R'<sub>HK</sub></td><td>–4.402</td><td>–4.48</td><td>–4.434</td></tr><tr><td>X-Ray luminosity, log L<sub>X</sub>/L<sub>★</sub></td><td>–4.71</td><td>–4.37</td><td>–4.43</td></tr><tr><td>Mass, M<sub>★</sub> [M<sub>☉</sub>]</td><td>1.28</td><td>0.98</td><td>1.02</td></tr><tr><td>Age (Gyr)</td><td>0.304, 0.184, 0.748(sin i)</td><td>0.530, 0.162, 0.291</td><td>0.385, 0.115, –</td></tr></tbody></table>

**Notes.** The ages given are derived using the following methods (left to right): log R'_{HK} activity index levels, ROSAT X-Ray luminosities L_X/L_{bol}, rotation period. For more details and references, see Sect. 4.1.

In each case, the observed excess is attributed to the presence of a debris disk associated with the star. The measured full width at half maximum (FWHM) of the sources (derived from a 2D Gaussian fit) are consistent with unresolved objects. One can therefore constrain the extension of the emitting area (adopting the distance of the objects listed in Table 3) to a diameter of less than 7.''1 (156 AU) for HIP 103389, 7.''1 (126 AU) for HIP107350, and 7.''2 (148 AU) for HIP 114948.

The spectral slopes of the three systems form Reduction 1 and Reduction 2 are listed in Table 2. The SEDs of all three sources are found to exhibit an unusually steep decrease at the wavelength range of 70 µm to 160 µm (steeper than a black body radiator in the Rayleigh-Jeans regime) from the results of Reduction 1. The steepness is evaluated based on the spectral index of the SEDs between two different wavelengths (see Sect. 3 for a detailed description). The significance of this steepness is evaluated using error propagation. However, the results from the new (but not necessarily better) version of the data reduction (Reduction 2) give a shallower decrease of the SED toward longer wavelengths. Thus, we emphasize that this is a potential discovery of a steep decrease of the SEDs of the three sources and that a deeper understanding of the data (in particular for very faint sources) than possible until now is necessary to make a final statement. In the subsequent sections, we discuss the consequences of this steepness of the SEDs under the assumption that it is real.

### 3. What is an unusually steep SED?

Infrared excesses observed from debris disk systems are the result of absorption of incident starlight by dust grains in the disk and subsequent re-emission of stellar radiation. Since the dust is colder than the star, its emission peaks at longer wavelengths (we assume that there is only one peak). In this section, a qualitative discussion is given of how steep the SED of a debris disk is expected to be. We first investigate the shape of the excess treating both the star and the disk as single temperature black body radiators. Later, we will discuss the effects of a more realistic

disk model, i.e., allowing for a range of temperatures and more realistic grain properties.

#### 3.1. Treating the dust as a single temperature black body

The spectral index $\Delta$ of an SED is defined as

$$ \Delta = \frac{\partial \log F_v}{\partial \log v} \qquad (1) $$

where $F_v$ is the total flux of the object at the frequency $v$. It is typically measured as the slope of the SED between two data points available, i.e.,

$$ \Delta_{v_1, v_2} = \frac{\log F_{v_2} - \log F_{v_1}}{\log v_2 - \log v_1}. \qquad (2) $$

We further denote the two frequencies $v_1$ and $v_2$ with their corresponding PACS wave bands (e.g., $\Delta_{100,70}$). The spectral index of a black body radiator amounts to $\Delta = 2$ in the Rayleigh-Jeans regime. It is decreasing toward shorter wavelengths where the Rayleigh-Jeans approximation is not valid. A spectral index of $\Delta > 2$ means that the SED falls off toward longer wavelengths steeper than a single temperature black body in its Rayleigh-Jeans regime. In the wavelength range where dust reemission is observed, the emission of the star can be very well approximated by a black body radiator in its Rayleigh-Jeans regime. It will only contribute as an additional component with a spectral index of $\Delta = 2$. Thus, the contribution of the star will – if significant – move the spectral slope of the whole system closer to 2.

#### 3.2. More exact treatment of the disk

Our first step to expand the above discussion is to allow for different dust temperatures, while the dust is still treated as a number of single temperature black bodies (i.e., a radially extended disk). A two component black body has been used successfully to model SEDs of debris disks (e.g., Hillenbrand et al. 2008).
<figure><img src="image_4.png" /><figcaption>Fig. 2. Absorption efficiency $Q_{abs}(\lambda)$ for astronomical silicate grains of different grain radii $a$ (see Sect. 3.2 for details).</figcaption></figure>

This results in a flatter SED ($\Delta_{v_1,v_2} < 2$) than in the above case, since the colder dust component will contribute additional flux in the long wavelength regime of the warmer dust component. Thus, using several components with different temperatures one can only reach $\Delta_{v_1,v_2} = 2$ in the long wavelength regime of the coldest dust component (i.e., at even longer wavelengths compared to the single black body case above). At shorter wavelengths $\Delta_{v_1,v_2}$ will always be <2.

Using Mie theory, one can derive absorption efficiencies $Q_{abs}(\lambda)$ for more realistic dust grains. The emission of a dust grain is then described by the following equation which is a modification of Planck's law:

$$B'_{\lambda}(\lambda, T_{\text{dust}}) = \frac{2hc^2}{\lambda^5} \cdot \frac{1}{\exp[hc/\lambda k T_{\text{dust}}] - 1} \cdot Q_{\text{abs}}(\lambda). \quad (3)$$

The shape of the simple black body SED is modified by the function $Q_{abs}(\lambda)$ which will result in a steeper decrease with wavelength in the Rayleigh-Jeans regime of the black body if $Q_{abs}(\lambda)$ is decreasing with wavelength. Figure 2 shows $Q_{abs}(\lambda)$ for astronomical silicate (Draine 2003) and different grain radii $a$. $Q_{abs}(\lambda)$ is close to 1 at short wavelengths, while it indeed exhibits a break at $\lambda \approx 2\pi a$ and then decreases toward longer wavelengths with $\lambda^{-2}$ (or $\nu^2$). This is analog to the discussion of the opacity index, e.g., by Draine (2006). An SED of a dust disk that is steeper in the long wavelength regime than a black body is obviously possible. However, debris disks are usually expected to be radially extended (e.g., Backman et al. 2009; Ertel et al. 2011; Löhne et al. 2012) and to contain particles of different size. This suggests a broad range of dust temperatures to be present resulting in a flattening of the SED. In addition, the grain size dependence of the break in $Q_{abs}(\lambda)$ means that the grains have to be smaller than $\approx \lambda/2\pi$ to result in an SED that falls off steeper than a black body (Draine 2006). This suggests grains smaller than $\approx 15 \mu m$ for the three disks presented in this paper (Table 2). Thus, one qualitatively expects a narrow dust belt composed of small grains to be present. While a spectral slope steeper than 2 is common for debris disks at (sub-)mm wavelengths, such a behavior is very unusual in the range of $\approx 100 \mu m$ among the few disks with characterized slope in this wavelength range (e.g., Hillenbrand et al. 2008; Roccatagliata et al. 2009). For the debris disks considered there, typical values of $\Delta$ at wavelengths around $100 \mu m$ are in the range of 0 to 2. Note that our Herschel observations are the first to characterize the spectral slope of a large number of debris disks at wavelengths around $100 \mu m$. Earlier

observations in this wavelength range suffer from low sensitivity. Thus they were limited to few bright disks or disks with very flat spectral slope in the far-ir. The completion of our survey and an analysis of the whole data set obtained (Eiroa et al., in prep.) will allow to place the three disks presented in a broader context.

We define an unusually steep SED to be an SED that exhibits a $\Delta_{v_1,v_2} > 2$ in at least one combination of the PACS bands of 70 $\mu m$, 100 $\mu m$, and 160 $\mu m$, and for which a $\Delta_{v_1,v_2} = 2$ would not be possible within the $1\sigma$ uncertainties. These criteria are met by the three debris disks presented in this work (Table 2).

## 4. SED modeling

In this section we perform detailed analytical model fitting to the observed SED data to explore quantitatively the conclusions on the disk properties derived in the above qualitative discussion of the steep slope of the observed SEDs.

### 4.1. Stellar properties

Table 3 lists some of the main stellar parameters and other relevant observational properties of HIP 103389, HIP 107350, and HIP 114948 (Eiroa et al., in prep.). These parameters were obtained using the DUNES Virtual Observatory tool<fnref n="3" />. The spectral types are taken from the HIPPARCOS catalog (Perryman et al. 1997; van Leeuwen 2007), and the range of spectral types from Skiff (2010). Bolometric luminosities and stellar radii have been estimated from the absolute magnitudes and the bolometric corrections using the measurements by Flower (1996); similar values are obtained using the bolometric correction procedure by Masana et al. (2006). Effective temperatures, gravities and metallicities are mean values of spectroscopic and photometric estimates from different works in the literature (HIP 103389 and HIP 114948: Allende Prieto & Lambert 1999; Gray et al. 2006; Holmberg et al. 2009; HIP 107350: Valenti & Fischer 2005; Fuhrmann 2008). Radial velocities are taken from Maldonado et al. (2010) and Kharchenko et al. (2007). Rotational velocities of the stars are taken from Holmberg et al. (2007) and Martínez-Arnáiz et al. (2010). Periods are taken from Reiners (2006) and Messina et al. (2003). The Li I 6707.8 Å line equivalent widths are taken from Maldonado et al. (2010) in the case of HIP 107350 and our own estimates using reduced HARPS<fnref n="4" /> data for the other two stars; for these latest objects the EW(Li I) have been corrected from the FeI 6707.4 Å line contamination as in Maldonado et al. (2010). All three stars show log $R'_{HK}$ activity index levels (Gray et al. 2006; Martínez-Arnáiz et al. 2010) indicating that they are active stars (e.g., Maldonado et al. 2010). The observed ROSAT X-Ray luminosities, $L_X/L_{bol}$, derived in the context of DUNES (Eiroa et al., in prep.) are also given.

The masses of the stars in Table 3 are estimated from their radii and gravities. We estimate the age of the stars (Table 3) from the log $R'_{HK}$ index, and periods using the relationships as given by Mamajek & Hillenbrand (2008) and from the X-Ray luminosity following Garcés et al. (2011). For all three stars these age tracers yield consistent values of less than ~500 Myr. We note that EW(Li I) values are close to those of the Hyades envelope with an age of 600 Myr. Rocha-Pinto et al. (2004) give an age of 370 Myr for HIP 107350 based on the log $R'_{HK}$ activity index, while isochrone-based ages are 2.3 Gyr for HIP 103389,

<footnote n="3">³ http://sdc.cab.inta-csic.es/dunes/searchform.jsp</footnote>

<footnote n="4">⁴ Based on observations made with the European Southern Observatory telescopes obtained from the ESO/ST-ECF Science Archive Facility.</footnote>
with a range of 0.6 Gyr to 3.8 Gyr (Holmberg et al. 2009), 3.1 Gyr for HIP 107350, with a range of 1.1 Gyr to 4.7 Gyr (Valenti & Fischer 2005), and 3.7 Gyr for HIP 114948, with a range of 2.1 Gyr to 5.2 Gyr (Holmberg et al. 2009). Given that the stars are located right on the main-sequence, stellar ages are difficult to estimate on the basis of isochrones, and that they are very sensitive to $T_{eff}$ and metallicity (Holmberg et al. 2009) we would favor in these cases the young values estimated from the activity tracers and rotational periods. This choice is also compatible with the kinematics of the stars and their ascription to the young moving groups, less than 1 Gyr, of HIP 107350 to the local association (~150 Myr) or the Hercules-Lyra (~200 Myr) association (Maldonado et al. 2010; López-Santiago et al. 2006), and HIP 114948 and HIP 103389 to the Castor Moving Group (~200 Myr) and the young disk population (<1 Gyr), respectively.

HIP 107350 hosts a T dwarf companion (HN Peg B, $M = 22 M_{\odot}$) with a projected physical separation of 794 AU (~43′′) from the host star (Luhman et al. 2007). This is sufficiently far away, so that the presence of this companion does not affect our present analysis. A very faint source (~2 mJy) is detected in the PACS 100 µm image close to the predicted position of HN Peg B. However, we do not expect HN Peg B to be detectable by our observations (~0.6 mJy at 8 µm). Thus, we attribute this detection to a coincidental alignment with a background source.

## 4.2. Modeling the excess

In this section, detailed analytical model fitting to the observed SED data is performed to explore quantitatively the conclusions on the disk properties derived in the above qualitative discussion. This is done using the results from Reduction 1, since these results give the steep SEDs. The fitting results are expected to differ less significantly from standard solutions using the photometric data from Reduction 2. For modeling of the three sources, we use SAnD which is part of the DUNES modeling toolbox (Löhne et al. 2012; Augereau et al., in prep.). A description of this tool can be found in the appendix. For the model fitting of the infrared excess SED, the above derived stellar properties are used. It is important to note that uncertainties of the stellar parameters as well as on the distance of the objects are not considered in our modeling of the excess. Thus, the derived uncertainties of the disk parameters are only the formal uncertainties from the fitting procedure. The uncertainties of the stellar parameters and distance would mostly affect the effective temperature and luminosity of the star and, thus, the heating of the dust – slightly increasing or decreasing its temperature. In particular, they will not affect significantly the slope of the derived grain size distribution.

For the comparison between modeled and observed measurements we use a reduced chi-squared

$$ \chi^2_{\text{red}} = \frac{1}{N_{\text{d.o.f.}}} \sum_{i=1}^{N} \left( \frac{F_i - \tilde{F}_i}{\sigma_i} \right)^2; \quad N_{\text{d.o.f.}} = N - N_{\text{free}}, \quad (4) $$

where $N$ is the number of SED data points used, $F_i$, $\tilde{F}_i$, and $\sigma_i$ are the measured flux, the modeled flux, and the uncertainty in a data point, $N_{free}$ is the number of free parameters, and $N_{d.o.f.}$ is the number of degrees of freedom in our fitting. For this reduced $\chi^2$ a value of 1.0 is generally desirable in the fitting if $N_{d.o.f.}$ is large, while in our case $N_{d.o.f.}$ ranges from 1 to 3 where the reduced $\chi^2$ is desirable to be as small as possible.

Only selected data points are used. This results in a total of 8 data points included in the fitting for each source. The SEDs of the sources can be seen in Fig. 3. The references for all photometric data points available for the three sources can be found in Eiroa et al. (in prep.). The selection is done by the following criteria:

- All flux measurements at wavelengths >10 µm are considered, but no upper limits, since they do not give significant additional information in the present case. Measurements at these wavelengths are also included, if they are photospheric, because they exclude significant emission from the disk at these wavelengths.

- Measurements at wavelengths <10 µm are not included, since they are at much shorter wavelengths than the shortest wavelength at which excess is detected (24 µm ... 32 µm).

- To account for the Spitzer/IRS spectrum in a consistent way, an additional photometric point at $\lambda = 32$ µm is extracted from these data following Hillenbrand et al. (2008) and is fitted along with the other photometric values. Therefore, 10 consecutive data points centered at $\lambda = 32$ µm are averaged. The uncertainty is computed by adding in quadrature the standard deviation of these data points and a 5% calibration uncertainty (Spitzer/IRS instrument handbook v4.0; Teplitz et al. 2011). This will prove to be sufficient to get the whole Spitzer/IRS spectrum properly reproduced by all the best-fit models.

The capabilities of SAnD are used to explore a very broad range of parameters. The following model is employed:

- radial surface density distribution $\Sigma(r) \propto r^{-\alpha}$ with inner and outer cut-off radii $r_{in}$ and $r_{out}$;

- differential grain size distribution $dn(a) \propto a^{-\gamma}da$ with lower and upper cut-off size $a_{min}$ and $a_{max}$;

- two possible grain compositions – pure astronomical silicate (Draine 2003) and a 1:1 mixture of astronomical silicate and ice (Löhne et al. 2012; Augereau et al., in prep.) – to explore the possibility that water ice might be a significant constituent of debris disk dust and might be responsible for the peculiar shape of the SEDs.

Since the number of free parameters in the fitting is close to the number of data points used, the fitting is expected to be very degenerate. To explore the parameter space in an efficient way and to find parameters that can be fixed, because they have a unique best-fit result or they have no significant effect on the fitting at all, we use different approaches. These approaches are described below. The parameter space explored for each approach is listed in Table 4. The fitting results are compiled in Table 5. Simulated SEDs from our best-fit models are shown in Fig. 3. The parameter space explored and the parameters fixed in each approach are motivated by the results of the previous approaches as described in the following.

### 4.2.1. Approach 1: Exploring the parameter space

In a first approach, a range of free parameters is considered that includes those used to explain most of the known debris disks. An upper grain size of 1 mm is chosen, large enough to consider any effect of large grains on the fluxes at all wavelengths observed. Seven free parameters are used ($r_{in}$, $r_{out}$, $\alpha$, $a_{min}$, $\gamma$, $M_{dust}$, and dust chemical composition). This is the largest number possible considering a total of 8 SED measurements. Although strong degeneracies are expected, this approach is used
<table><caption>**Table 4.** Explored parameter space for the fitting based on simulated thermal annealing (for details see Sect. 4.2).</caption><thead><tr><th rowspan="2">Parameter</th><th colspan="2">Approach 1</th><th colspan="2">Approach 2</th><th colspan="2">Approach 3</th><th colspan="2">Approach 4</th><th rowspan="2">distribution</th></tr><tr><th>Range</th><th># Values</th><th>Range</th><th># Values</th><th>Range</th><th># Values</th><th>Range</th><th># Values</th></tr></thead><tbody><tr><td>R<sub>in</sub> [AU]</td><td>3–100</td><td>817</td><td>3–100</td><td>817</td><td>3–100</td><td>817</td><td>3–100</td><td>817</td><td>temp</td></tr><tr><td>R<sub>out</sub> [AU]</td><td>5–300</td><td>671</td><td>5–300</td><td>671</td><td>5–300</td><td>671</td><td>5–300</td><td>671</td><td>temp</td></tr><tr><td>α</td><td>-2.0–3.5</td><td>56</td><td>fixed</td><td>1</td><td>fixed</td><td>1</td><td>fixed</td><td>1</td><td>linear</td></tr><tr><td>a<sub>min</sub> [μm]</td><td>0.2–20.0</td><td>448</td><td>0.2–20.0</td><td>448</td><td>0.2–20.0</td><td>448</td><td>0.2–20.0</td><td>448</td><td>logarithmic</td></tr><tr><td>a<sub>max</sub> [μm]</td><td>fixed</td><td>1</td><td>fixed</td><td>1</td><td>fixed</td><td>1</td><td>1.0–1000.0</td><td>674</td><td>logarithmic</td></tr><tr><td>γ</td><td>2.0–5.0</td><td>31</td><td>2.0–10.0</td><td>81</td><td>fixed</td><td>1</td><td>fixed</td><td>1</td><td>linear</td></tr><tr><td>M<sub>dust</sub> [M<sub>⊕</sub>]</td><td>free</td><td>...</td><td>free</td><td>...</td><td>free</td><td>...</td><td>free</td><td>...</td><td>continuous</td></tr><tr><td>Composition</td><td>0%, 50% ice</td><td>2</td><td>0%, 50% ice</td><td>2</td><td>0%, 50% ice</td><td>2</td><td>0%, 50% ice</td><td>2</td><td>...</td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Observed and modeled SEDs of the three sources. The *Spitzer*/IRS spectrum is plotted in red. The data points considered for the fitting are shown in dark blue. The PACS data considered are from Reduction 1, while the results from Reduction 2 are plotted in light blue for comparison. *Top*: final models from the SAnd approaches. The modeled SEDs (solid line) for the three disks are computed from the results of Approach 4 using pure astronomical silicate (Sect. 5). The dashed line in each panel represents the model stellar photosphere. *Bottom*: comparison between the best-fit results from the different approaches. The modeled SEDs are computed from the best-fit result of each approach using pure astronomical silicate. The results from Approach 1 are not included, since the best-fit results are found not to be in the explored range of parameters in this approach.</figcaption></figure>

to explore the parameter space without strong initial constraints that might bias the results. Parameters that have a unique solution or that have no significant effect on the fit can be fixed in subsequent approaches.

A very narrow ring structure is found to be the best-fit for all three disks, although the parameters describing the spatial dust distribution are not constrained very well. Furthermore, a very steep grain size distribution with a value for the exponent $\gamma$ very close to the edge of the explored parameter space (5.0) is found as best-fit for all three SEDs (Table 5). In the following, it is referred to this as an underabundance of large grains. Since this is an atypical, but not unexpected result (Sect. 3), more effort is put on the evaluation of the significance of this result in the subsequent approaches. As the radial extent of the disks has been

found to be narrow, the exponent of the radial surface density distribution is fixed in the further fits to $\alpha = 0.0$ (constant surface density), decreasing the explored parameter space by one dimension. In a narrow ring, this will not have any effect on the SED.

### 4.2.2. Approach 2: A larger range of possible values for $\gamma$

Now, the aim is to obtain a best-fit value for $\gamma$ that is included in the considered parameter space. A range of possible values for $\gamma$ of 2.0 to 10.0 is explored. All other parameters (beside $\alpha$ which is now fixed as described above) have the same ranges as in Approach 1. This results in a total of 6 free parameters.
<table><caption>Table 5. Results from the fitting based on simulated thermal annealing (for details see Sect. 4.2).</caption><thead><tr><th rowspan="3">Parameter</th><th colspan="10">Best-fit value [3σ confidence levels]</th></tr><tr><th colspan="10">HIP 103389</th></tr><tr><th colspan="2">Approach 1</th><th colspan="2">Approach 2</th><th colspan="2">Approach 3</th><th colspan="2">Approach 4</th></tr><tr><th></th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th></tr></thead><tbody><tr><td>R<sub>in</sub> [AU]</td><td>20.7 [8.5–32.2]</td><td>20.8 [9.1–35.9]</td><td>18.2 [7.9–24.5]</td><td>20.9 [8.0–26.6]</td><td>11.9 [4.2–16.0]</td><td>13.6 [4.4–19.1]</td><td>42.3 [12.9–63.3]</td><td>22.5 [8.6–43.0]</td></tr><tr><td>R<sub>out</sub> [AU]</td><td>20.7 [15.3–53.1]</td><td>20.9 [15.4–123.3]</td><td>20.0 [16.3–65.1]</td><td>20.9 [17.3–77.3]</td><td>12.0 [8.8–37.4]</td><td>13.6 [10.2–38.5]</td><td>46.0 [21.3–138.5]</td><td>22.5 [17.9–110.9]</td></tr><tr><td>α</td><td>0.1 [–2.0–3.5]</td><td>2.6 [–2.0–3.5]</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td></tr><tr><td>a<sub>min</sub> [μm]</td><td>6.6 [4.7–8.6]</td><td>9.2 [4.4–11.4]</td><td>9.5 [7.8–10.4]</td><td>12.8 [10.7–13.6]</td><td>6.1 [2.8–9.8]</td><td>6.6 [2.6–12.0]</td><td>4.2 [3.1–7.8]</td><td>9.3 [3.7–14.9]</td></tr><tr><td>a<sub>max</sub> [μm]</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>14.3 [12.7–18.1]</td><td>22.7 [14.3–28.5]</td></tr><tr><td>γ</td><td>5.0 [4.6–5.0]</td><td>5.0 [4.5–5.0]</td><td>7.4 [6.3–10.0]</td><td>9.0 [6.5–10.0]</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td></tr><tr><td>M<sub>dust</sub> [M<sub>⊕</sub>]</td><td>1.98e-5</td><td>1.67e-5</td><td>1.30e-5</td><td>1.29e-5</td><td>4.69e-5</td><td>4.03e-5</td><td>4.86e-5</td><td>1.45e-5</td></tr><tr><td>reduced χ<sup>2</sup></td><td>4.112</td><td>3.664</td><td>0.816</td><td>0.776</td><td>10.277</td><td>8.488</td><td>0.628</td><td>0.760</td></tr></tbody></table><table><thead><tr><th rowspan="3"></th><th colspan="10">HIP 107350</th></tr><tr><th colspan="2">Approach 1</th><th colspan="2">Approach 2</th><th colspan="2">Approach 3</th><th colspan="2">Approach 4</th></tr><tr><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th></tr></thead><tbody><tr><td>R<sub>in</sub> [AU]</td><td>15.3 [3.8–28.5]</td><td>21.5 [7.4–93.0]</td><td>29.1 [7.9–47.2]</td><td>30.6 [5.6–44.1]</td><td>9.6 [3.0–15.6]</td><td>10.9 [3.0–16.5]</td><td>37.1 [4.3–54.4]</td><td>35.2 [7.2–54.2]</td></tr><tr><td>R<sub>out</sub> [AU]</td><td>15.9 [9.3–227.3]</td><td>21.8 [10.3–276.6]</td><td>31.3 [13.9–113.9]</td><td>32.3 [16.0–138.5]</td><td>9.6 [5.8–33.0]</td><td>11.0 [9.1–35.7]</td><td>37.4 [19.2–187.5]</td><td>35.2 [17.4–145.3]</td></tr><tr><td>α</td><td>–0.9 [–2.0–3.5]</td><td>1.4 [–2.0–3.5]</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td></tr><tr><td>a<sub>min</sub> [μm]</td><td>5.9 [4.7–9.0]</td><td>5.9 [4.3–9.1]</td><td>6.9 [2.7–10.9]</td><td>8.2 [3.6–10.7]</td><td>5.7 [4.8–11.0]</td><td>5.8 [4.9–11.7]</td><td>7.8 [1.6–10.4]</td><td>9.6 [2.2–13.5]</td></tr><tr><td>a<sub>max</sub> [μm]</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>7.8 [6.3–13.3]</td><td>9.6 [6.4–17.4]</td></tr><tr><td>γ</td><td>5.0 [4.0–5.0]</td><td>5.0 [4.0–5.0]</td><td>10.0 [6.0–10.0]</td><td>10.0 [5.9–10.0]</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td></tr><tr><td>M<sub>dust</sub> [M<sub>⊕</sub>]</td><td>4.59e-6</td><td>6.03e-6</td><td>1.11e-5</td><td>9.29e-6</td><td>1.05e-5</td><td>9.09e-6</td><td>1.54e-5</td><td>1.08e-5</td></tr><tr><td>reduced χ<sup>2</sup></td><td>7.048</td><td>5.672</td><td>1.652</td><td>1.568</td><td>4.229</td><td>3.869</td><td>1.528</td><td>1.488</td></tr></tbody></table><table><thead><tr><th rowspan="3"></th><th colspan="10">HIP 114948</th></tr><tr><th colspan="2">Approach 1</th><th colspan="2">Approach 2</th><th colspan="2">Approach 3</th><th colspan="2">Approach 4</th></tr><tr><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th><th>Silicate</th><th>Mixture</th></tr></thead><tbody><tr><td>R<sub>in</sub> [AU]</td><td>12.8 [6.8–13.9]</td><td>13.6 [6.8–14.5]</td><td>12.8 [7.1–13.7]</td><td>13.5 [7.8–14.9]</td><td>12.7 [5.5–16.6]</td><td>14.0 [7.9–19.1]</td><td>32.5 [9.0–40.1]</td><td>13.3 [8.4–14.5]</td></tr><tr><td>R<sub>out</sub> [AU]</td><td>12.8 [12.1–39.8]</td><td>13.8 [12.9–34.7]</td><td>12.8 [12.1–23.2]</td><td>13.8 [12.9–24.9]</td><td>12.7 [10.5–25.6]</td><td>14.1 [11.0–30.8]</td><td>34.5 [26.9–81.4]</td><td>13.8 [12.9–22.9]</td></tr><tr><td>α</td><td>2.9 [–2.0–3.5]</td><td>0.0 [–2.0–3.5]</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td><td>0.0 (fixed)</td></tr><tr><td>a<sub>min</sub> [μm]</td><td>9.6 [8.9–10.2]</td><td>12.9 [11.9–13.8]</td><td>9.6 [8.9–10.2]</td><td>12.9 [11.9–13.8]</td><td>6.4 [3.6–9.9]</td><td>7.4 [3.0–11.7]</td><td>3.2 [2.6–4.6]</td><td>11.8 [10.7–13.1]</td></tr><tr><td>a<sub>max</sub> [μm]</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>1000.0 (fixed)</td><td>24.2 [10.0–27.3]</td><td>43.8 [10.0–53.2]</td></tr><tr><td>γ</td><td>4.7 [4.5–5.0]</td><td>4.7 [4.4–5.0]</td><td>4.7 [4.4–5.2]</td><td>4.7 [4.4–5.1]</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td><td>3.5 (fixed)</td></tr><tr><td>M<sub>dust</sub> [M<sub>⊕</sub>]</td><td>2.19e-5</td><td>2.16e-5</td><td>2.19e-5</td><td>2.16e-5</td><td>8.82e-5</td><td>7.39e-5</td><td>5.93e-5</td><td>1.61e-5</td></tr><tr><td>reduced χ<sup>2</sup></td><td>0.568</td><td>0.464</td><td>0.284</td><td>0.232</td><td>3.291</td><td>2.333</td><td>0.132</td><td>0.232</td></tr></tbody></table>
For HIP 107350, the range of explored values of $\gamma$ seems to be still too small. This is ignored, since such a very large value is not significant (confidence levels: 5.9 ... 10.0, Table 5). For the other two objects, the value of $\gamma$ is well within the explored range. The confidence levels of the inner and outer radius suggest that the ring-like shape of the disk is not very significant, even for a constant surface density.

### 4.2.3. Approach 3: Fixing $\gamma$ to 3.5

Fixing the value of $\gamma$ to 3.5 and having a fixed upper grain size of $a_{max} = 1.0$ mm, one can force a fit where the parameters of the grain size distribution are consistent with an equilibrium collisional cascade (Dohnanyi 1969). The explored ranges of all other parameters are unchanged compared to Approach 2. This results in a total of 5 free parameters.

The resulting $\chi^2$ is much worse (by a factor of 2.5 to 12.6, Table 5) compared to the results from Approach 2. The changes of the radial position of the dust ring are in line with the expectations from the much lower abundance of small (warm) particles.

### 4.2.4. Approach 4: Fixed value of $\gamma = 3.5$, but free upper grain size

Another possibility to produce an underabundance of large grains in the model is to let the upper grain size be a free parameter. Since a very steep size distribution and a small upper grain size have comparable results (removing large grains from the model), one can fix the value of $\gamma$ to 3.5 to be consistent with an equilibrium collisional cascade. The resulting upper grain size is expected to be sufficiently small. This results in a total of 6 free parameters.

This approach gives in general the best $\chi^2$ (Table 5). Upper grain sizes of few tens of micron are found. The low radial extend of the debris ring is again not very significant, but still the best-fit.

### 4.2.5. Grid search and inclusion of spatial information for HIP 114948

From our discussion in Sect. 3 we are not able to constrain our parameter space by any reasonable physical assumptions. Since SAnD is designed to explore a large, high-dimensional parameter space using a statistical approach, it was particularly useful to find the global best-fit to the data in the above approaches. The over-all result of our fitting (a lack of large grains) is clearly visible for HIP 103389 and HIP 107350. For HIP 114948 the result is not that stringent. The value of $\gamma = 4.7^{+0.5}_{-0.3}$ from Approach 2 of the SAnD fitting is significantly smaller than the values derived from the other two disks (although still significantly larger than the value of 3.5 expected from a standard equilibrium collisional cascade, Dohnanyi 1969). Thus, we perform a grid search of the parameter space using GRaTer (Augereau et al. 1999; Lebreton et al. 2012) to confirm our results from the statistical approach of SAnD. Furthermore, the constraint that the source is spatially unresolved is included in this fitting. This is not expected to significantly change the best-fit results (since these results are consistent with an unresolved disk) but to further constrain the range of values possible within the uncertainties of the derived values (since very extended configurations that would clearly be resolved are also included in the confidence levels from pure SED fitting). The results from fitting with and without the spatial constraints are compared.

Since the general result from the SAnD fitting (steep grain size distribution, large lower grain size) is independent from the approach used, we decide to use the standard approach of GRaTer. This is a power-law approach ($dn(a) \propto a^{-\gamma}da$) with upper and lower limit for the differential grain size distribution and a two-power-law approach for the radial surface density distribution:

$$ \Sigma(R) \propto \sqrt{2 \left( (R/R_0)^{-2\alpha_{in}} + (R/R_0)^{2\alpha} \right)^{-1}}, \quad (5) $$

where $\alpha$ and $\alpha_{in}$ are slopes of the radial surface density distribution outside and inside a peak position. We fix $\alpha_{in} = 10$ to get a sharp inner edge close to $R_0$ comparable to the single power-law approach used for the SAnD fitting.

The spatial constraints are included in our fitting using a radial brightness profile derived by radial averaging of the sources image in our 100 µm data. After a Gaussian is fitted to the image to derive the source center, the average flux and distance from the center is derived for all pixels in one pixel wide radial bins. The uncertainty on the flux in each bin is derived using the standard deviation of the pixel values in each radial bin. Images are simulated from the models, convolved with the PACS PSF measured on the bright standard star $\alpha$ Boo and profiles are extracted. These profiles are compared to the one derived from the data. The same number of SED data points as before is considered. In addition, a total of 10 profile points are used. An over-all reduced $\chi^2$ including SED and profile data is computed following Eq. (A.1) (using only one profile). The explored parameter space and best-fit values for the fit with and without the profile included are given in Table 6. These results should be compared to the results from Approach 2 of the SAnD fitting, since these two approaches are the most similar. They are found to be consistent. Maps of the distribution of the reduced $\chi^2$ in the parameter space for the fits including and not including the profile have been produced as cuts through the parameter space in three planes. The remaining parameters are fixed to their best-fit values, respectively. These maps are displayed in Fig. 4. We find that the best-fit values do not change significantly when including the spatial information. This confirms that the results from the pure SED fitting are fully consistent with a spatially unresolved disk. Figure 4 illustrates a number of results from the mapping of the parameter space:

- there is a correlation between the lower dust grain size $a_{min}$ and the exponent $\gamma$ of the size distribution. However, from our SED fitting one can rule out that any combination of significantly smaller grains and smaller exponent of the size distribution results in reasonable fits to the data (which would be closer to the expectations from a collisional cascade and grains close to the blow-out size of the system of $a < 1.0$ µm);

- there is nearly no correlation between the radial position of the dust belt and the exponent of the grain size distribution in both the case without and with the spatial information included;

- the degeneracy between the lower dust grain size and the radial position of the dust belt known from SED fitting is clearly visible. It is also obvious that the inclusion of the spatial information in part breaks this degeneracy since it excludes a large (≥70 AU) radial distance of the belt from the central star;

- in particular, the results of $a_{min}$ and $\gamma$ seem to be less constrained when including the spatial information. This is in part because the spatial information gives no constraints on
<table><caption>**Table 6.** Fitting result for HIP 114948 using the grid search method with GRaTer.</caption><thead><tr><th>Parameter</th><th>Range</th><th>Spacing</th><th># Values</th><th>Fit w/o spatial constraints</th><th>Fit w/ spatial constraints</th></tr></thead><tbody><tr><td>R<sub>0</sub> [AU]</td><td>1.0–250.0</td><td>logarithmic</td><td>60</td><td>13.7</td><td>13.7</td></tr><tr><td>α</td><td>0.0–10.0</td><td>linear</td><td>11</td><td>8.0</td><td>9.0</td></tr><tr><td>a<sub>min</sub> [μm]</td><td>0.05–43.5</td><td>logarithmic</td><td>55</td><td>14.1</td><td>14.1</td></tr><tr><td>γ</td><td>2.7–5.9</td><td>linear</td><td>17</td><td>4.7</td><td>4.7</td></tr><tr><td>V<sub>ice</sub>/(V<sub>si</sub> + V<sub>ice</sub>)</td><td>0.0–0.9</td><td>linear</td><td>10</td><td>0.6</td><td>0.6</td></tr><tr><td>M<sub>dust</sub> [M<sub>⊕</sub>]</td><td>1.0 × 10<sup>-6</sup>–4.0 × 10<sup>-3</sup></td><td>logarithmic</td><td>100</td><td>2.36 × 10<sup>-5</sup></td><td>2.20 × 10<sup>-5</sup></td></tr><tr><td>reduced χ<sup>2</sup></td><td>—</td><td>—</td><td>—</td><td>0.19</td><td>0.12</td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>**Fig. 4.** χ² maps for HIP 114948. The maps are created as two dimensional cuts through the χ² distribution in the parameter space searched with GRaTer without (top) and with (bottom) the spatial information included. The parameters not plotted are fixed to their best-fit value, respectively (Table 6). The correlations between the three most relevant parameters (γ, a<sub>min</sub>, and R<sub>0</sub>) are shown. The white cross denotes the position of the best-fit (lowest χ²). The dashed vertical line denotes the value of γ = 3.5 expected from a standard equilibrium collisional cascade.</figcaption></figure>

the parameters of the grain size distribution (beside the partly broken degeneracy of the lower dust grain size and the radial position of the dust which results in an indirect constraint on a<sub>min</sub>). Thus, including the radial profile we include a number of measurements that are not sensitive to changes in γ, which results in a less significant increase of the reduced χ² with changes of γ in this case. In addition, the data points in the profile can not be considered to be independent from each other. Although this is mitigated in part by the weighting in Eqs. (A.1) to (A.3), the combination of χ² derived from SED and profile is not fully consistent, in particular in the case of unresolved sources. To measure the quantity and uncertainty of a model parameter through model fitting using reduced χ², one should in general only include measurements in the fitting process that are sensitive to the parameter (as far as one can tell in advance) and independent from each other. Thus, the uncertainties on γ derived from pure SED fitting should be taken more reliable here.

## 5. Results

The general results and conclusions that can be drawn from the fitting approaches are summarized in the following.

- Values of γ = 4.7 to 10.0 or a small upper grain size are found rather than the expected value of γ ≈ 3.5 and a large a<sub>max</sub>. Models with γ = 3.5 do not reproduce the observed SEDs in a reasonable way. This result is not strongly affected by any modeling degeneracies found (that are included in the estimate of the uncertainties).

- Evidence is found for a large lower grain size of the dust compared to the expected blow-out size of the systems (0.5 µm to 0.9 µm). However, the lower grain size depends very much on the approach used. Earlier studies of debris disk SEDs had similar results (e.g., HD 107146; Roccatagliata et al. 2009) that have been at least mitigated and partially attributed to modeling degeneracies in pure SED fitting after including resolved data in the fitting process (Ertel et al. 2011). Thus, this result has to be treated with caution.

- The dust in all three disks must be located at a distance from the star of no more than a few tens of AU as can be concluded from the modeling results and the fact that the disks are spatially unresolved in the PACS images. Stronger constraints on this value are not possible due to degeneracies in the modeling.

- The disks appear to be narrow rings. However, the uncertainties – in particular of the outer disk radius – are very large and also very broad disks can result in fits on the SED that fall within the derived 3σ confidence levels. This can only partly be excluded by the additional constraint that the disks are spatially unresolved.

It is important to note that these results depend in part on the model used and on the parameter space explored. For example, the composition and shape of the dust grains have not been explored in detail in our modeling, but might allow one to reproduce the SEDs with models that are more in line with a standard equilibrium collisional cascade.

Furthermore, the results from Reduction 2 have not been considered in the fitting. Considering these photometric measurements and uncertainties instead of those from Reduction 1 would lower the significance of the underabundance of large grains. As a result, values are expected to be found for the model parameters that are more in line with those expected from an equilibrium collisional cascade. Without a deeper understanding of the data obtained with Herschel, it is not possible to make a final statement about the significance of our modeling results.

Each approach presented results in a different set of best-fit parameters and a different χ², while the results can be interpreted
<table><caption>**Table 7.** Simulated observational properties of the disks derived from our final models.</caption><thead><tr><th rowspan="2">Source</th><th rowspan="2">L<sub>dust</sub>/L<sub>*</sub></th><th colspan="2">(F<sub>dust</sub>/F<sub>*</sub>)<sub>0.6µm</sub></th><th colspan="2">(F<sub>dust</sub>/F<sub>*</sub>)<sub>0.8µm</sub></th><th colspan="2">(F<sub>dust</sub>/F<sub>*</sub>)<sub>1.1µm</sub></th><th colspan="2">(F<sub>dust</sub>/F<sub>*</sub>)<sub>2.2µm</sub></th></tr><tr><th>face-on</th><th>edge-on</th><th>face-on</th><th>edge-on</th><th>face-on</th><th>edge-on</th><th>face-on</th><th>edge-on</th></tr></thead><tbody><tr><td>HIP 103389</td><td>1.5 × 10<sup>-5</sup></td><td>1.1 × 10<sup>-6</sup></td><td>1.9 × 10<sup>-4</sup></td><td>1.1 × 10<sup>-6</sup></td><td>1.7 × 10<sup>-4</sup></td><td>1.3 × 10<sup>-6</sup></td><td>1.8 × 10<sup>-4</sup></td><td>1.8 × 10<sup>-6</sup></td><td>1.5 × 10<sup>-4</sup></td></tr><tr><td>HIP 107350</td><td>0.6 × 10<sup>-5</sup></td><td>4.9 × 10<sup>-7</sup></td><td>7.4 × 10<sup>-6</sup></td><td>4.6 × 10<sup>-7</sup></td><td>8.5 × 10<sup>-6</sup></td><td>5.5 × 10<sup>-7</sup></td><td>1.3 × 10<sup>-5</sup></td><td>6.9 × 10<sup>-7</sup></td><td>3.9 × 10<sup>-5</sup></td></tr><tr><td>HIP 114948</td><td>2.5 × 10<sup>-5</sup></td><td>2.1 × 10<sup>-6</sup></td><td>3.4 × 10<sup>-4</sup></td><td>2.0 × 10<sup>-6</sup></td><td>3.0 × 10<sup>-4</sup></td><td>2.6 × 10<sup>-6</sup></td><td>3.3 × 10<sup>-4</sup></td><td>3.7 × 10<sup>-6</sup></td><td>2.6 × 10<sup>-4</sup></td></tr></tbody></table>

in a consistent way. To find from the fitting results a final model of each disk, the following selection is applied:

- The results from Approach 1 have to be ruled out, because the best-fit parameters from this approach are close to the edge of the explored parameter space.

- In general, the SEDs cannot be fitted by a grain size distribution that includes a significant amount of large grains (Approach 3).

- The results from Approach 2 (free $\gamma$, $a_{max} = 1.0$ mm) and Approach 4 (free $a_{max}$, $\gamma = 3.5$) are consistent, while the results from Approach 4 and pure silicate represent in general the best fit to the data. Thus, these results are considered as the best-fit models from the fitting.

From these final models, one finds dominating dust temperatures (the temperature of the smallest grains that dominate the emission) of 63 K, 50 K, and 73 K for HIP 103389, HIP 107350, and HIP 114948.

## 5.1. The origin of the dust

The existence of micron-sized dust grains around main sequence stars is usually explained by collisions of larger bodies such as planetesimals (e.g., Krivov 2010). In such a scenario one would expect a significant amount of larger grains (several tens of micron to millimeter sized particles) to be present that are produced by the collisional cascade also producing the smallest grains (Dohnanyi 1969). These grains significantly contribute to the SED in the observed wavelength regime. Assuming that the steep slopes of the SEDs are real, we find from our modeling strong evidence for an underabundance of these larger grains. To our knowledge, these are the first debris disks discovered, that potentially exhibit such a peculiar shape of the SED in the wavelength regime of 70 µm to 160 µm. All this would suggest that these objects are exceptions from the common understanding of dust creation in debris disks. On the other hand, six more disks with similar shape of the SED have recently been identified in our DUNES survey, which might suggest that this phenomenon is common among low-mass (i.e., low luminosity, Tables 5 and 7) debris disks<fnref n="5" />. The fact that they are very faint (Table 7) might imply that we are faced with a new class of debris disks that were not revealed earlier due to limited sensitivity in the relevant wavelength range. In the following, we discuss scenarios that might be capable to explain the presence of a disk as it has been modeled for the three systems.

### Scenario 1: Significant deviation from the conditions required for a standard equilibrium collisional cascade

For the standard equilibrium collisional cascade, a number of assumptions are made that are not necessarily valid in debris disks due to the effects of radiation pressure and Poynting-Robertson (PR) drag (Sect. 1). For the massive debris disks known so far,

PR drag is expected to be negligible (Wyatt 2005). Following Backman & Paresce (1993), one can estimate the ratio between collisional time scale and PR time scale for the three disks to $10^{-2}...10^{-3}$. However, the disk models used (best-fit from Approach 4, pure silicate) represent very narrow rings, which is not a significant fitting result. Assuming a ring with a width of 10 AU starting at the inner radius found from the fitting, one finds $t_{coll}/t_{PR} = 3 \times 10^{-2}...5 \times 10^{-2}$. It is not possible to put strong constraints on the dust dynamics due to the crude estimate of the time scales. Thus, transport mechanisms may play a role in the dust dynamics of the small grains in these disks. The effects of these processes have been modeled, e.g., by Krivov et al. (2006); Thébault & Augereau (2007); Wyatt et al. (2011). In particular, radiation pressure leads to a wavy structure of the grain size distribution increasing the slope of the grain size distribution in certain ranges. PR drag results in a depletion of small grains, increasing the dominating grain size, but mitigating the wavy structure of the size distribution. Pan & Schlichting (2012) found that accounting self-consistently for the velocity evolution in a collisional cascade can increase the derived value of $\gamma$ from the standard value of 3.5 to up to 4.0.

### Scenario 2: Different chemical composition or physical shape of the dust grains than the assumed one

One might imagine a dust composition emitting significantly more efficiently than astronomical silicate or ice in the wavelength range of 70 µm to 100 µm compared to longer wavelengths. Therefore, one would need grains with a break in $Q_{abs}(\lambda)$ (Sect. 3) at shorter wavelengths (in particular for large grains). Also, if the large grains were significantly colder than expected, their emission would be reduced. This would both allow a significantly larger amount of large grains to be present in the systems than modeled. Furthermore, Mie calculations might result in an incorrect approximation of $Q_{abs}$ particularly for large grains (tens of microns or larger), e.g., due to incorrect assumptions on the physical shape of the grains.

Voshchinnikov et al. (2006) simulated the effect of porosity on the absorption efficiency of dust grains. They found that the temperature of dust grains at a given distance from a star decreases significantly with increasing porosity. If one now has small, compact grains and larger, porous grains, the difference in temperatures of the smaller and larger grains can be increased significantly. This would be the case, if the large grains producing the smaller ones through collisions were porous, composed of smaller, compact units in the order of the lower grain size derived from the modeling of the three systems (i.e., 5 µm to 10 µm). Such a scenario would result in an overabundance of grains of this size compared to the larger, weaker grains. Such an overabundance would shift the lower cut-off of our adopted power-law distribution toward this grain size, although smaller grains are present in the system.

### Scenario 3: A shepherding planet

Planets can trap dust particles into mean-motion resonances (e.g., Wyatt 2006). This results in a barrier against particles moving inward due to PR drag. This barrier is less efficient for very

<footnote n="5">⁵ For an overview of the DUNES results we refer to Eiroa et al. (in prep.).</footnote>
small grains for which PR drag is very strong (Reidemeister et al. 2011). On the other hand, very large grains and planetesimals are not significantly affected by PR drag.

The dust seen in these disks might be produced in a faint, transport dominated (Krivov 2010) debris disk further away from the star, too cold and too faint to be detected in the available data. The small particles (few tens of micron and smaller) would then be dragged inward by PR drag. A possible planet present further inward in the system should trap the particles into resonance. This would result in an accumulation of particles with a very distinct range of sizes, which would explain both the lack of large grains (at least in an abundance and at a position where they would be detectable) and the lack of small grains. Such a “dust trail” has been observed, e.g., to be associated with our Earth and has been modeled to predominantly consist of grains of ≈12 µm and larger (Dermott et al. 1994). It will have to be evaluated, whether an outer planetesimal belt capable to produce enough dust through collisions to replenish the dust in the disk can be faint enough not to be detected by the Herschel observations, and under which conditions the planet can trap enough dust into resonance to produce the disk observed. This has to be done through detailed dynamical and collisional modeling (e.g., Reidemeister et al. 2011), which is not within the scope of the present work.

## 5.2. Observational perspectives

Further observations with present and near future instruments can help to reasonably increase our understanding of these debris disks. In the following, we discuss a number of observations that might result in very valuable data.

- Since the most extended models of the three disks predict radial extents of only ≈2′′2, these objects are at the edge of what is resolvable with coronagraphy. Furthermore, the disks are very faint in scattered light (Table 7). It is important to note that in the case of edge-on orientation most of the flux comes from forward scattering and will then be concentrated close to the star (projected physical separation) with no contribution to the signal in coronagraphic observations. A contrast ratio of <10⁻⁶ is not accessible to present instruments.

- Successful optical/near-infrared imaging of planetary companions would give strong evidence that Scenario 3 (shepherding planet) is responsible for the peculiarities of the disk. Determining the position of this planet would also help to further constrain the position of the dust in this scenario. The youth of the stars means good chances to directly image giant planets at separations of ≥1" from the star with present methods (e.g., Marois et al. 2010).

- The predicted extents of the disks are in a range easily resolvable with ALMA. However, the disks are already very faint at PACS wavelengths and the expected surface brightness is decreasing particularly steep toward wavelengths accessible with ALMA. It is not clear, whether ALMA observations are sensitive enough to detect these disks. On the other hand, at (sub-)mm wavelengths even upper limits in the sensitivity range reachable are expected to provide useful further constraints on the shape of the SEDs. From the results of our fitting Approach 4, predicted fluxes of star and disk at 350 µm (the shortest wave band accessible with ALMA) are 0.55 mJy and 0.11 mJy for HIP 103389, 0.53 mJy and 0.05 mJy for HIP 107350, and 0.60 mJy and 0.20 mJy for HIP 114948, respectively.

- Further photometry and spectroscopy in thermal reemission would help to constrain the radial distribution of the dust and the size distribution of the small grains (wavelengths between 40 µm and 70 µm) and the steep slope and the shape of the SEDs at longer wavelengths (λ > 70 µm). The Stratospheric Observatory for Infrared Astronomy (SOFIA) is most promising to provide observational capabilities in the relevant wavelength regime in the near future.

## 6. Conclusions

The first data at wavelengths >70 µm for the three debris disks discussed in this paper have been presented and modeled. All three sources potentially exhibit an unusually steep decrease of the SED in the wavelength range between 70 µm and 160 µm. In a general discussion, it has been shown that this peculiar shape of the SED is an indicator for a deviation from the case of a standard equilibrium collisional cascade assuming standard grain composition (astronomical silicate or ice) and shape (compact, spherical grains). However, it has also been shown that the understanding of faint source photometry obtained with Herschel is still incomplete and that the results presented here depend very much on the actual version of the data reduction pipeline. Provided that the steep decrease of the SEDs is real, modeling implies that the thermal emission from these disks is dominated by a very distinct grain size regime of several micron to few tens of micron. The disks have been modeled as narrow rings with a significant underabundance of large grains. A number of possible explanations for such an unusually steep shape of the SEDs have been discussed. Six more candidates for this new class of debris disks have been identified so far from the ongoing DUNES survey. This is the first published discovery of debris disks that potentially exhibit such a peculiar shape of the SED in this wavelength regime.

*Acknowledgements.* We thank Kate Su for the re-reduction of the Spitzer/MIPS photometry used in this work. Furthermore, we thank the whole DUNES team for valuable discussion. S. Ertel thanks for financial support from DFG under contract WO 857/7-1 and for general support from K. Ertel. C. Eiroa, J. Maldonado, J. P. Marshall, and B. Montesinos are partially supported by Spanish grant AYA 2008/01727. J.-C. Augereau and J. Lebreton thank financial support through PNP-CNES. A. V. Krivov ans T. Löhne thank for financial support from DFG under contracts KR 2164/9-1 and LO 1715/1-1. O. Absil is supported by an F.R.S.-FNRS Postdoctoral Fellowship. S. Ertel, J.-C. Augereau and J. Lebreton thank the French National Research Agency (ANR) for financial support through contract ANR-2010 BLAN-0505-01 (EXOZODI). This work was partly funded by the Fundação para a Ciência e a Tecnologia (FCT) through the project PEst-OE/EEI/UI0066/2011.

## Appendix A: Description of SAnD

When performing analytical model fitting of debris disks, one is faced with a number of challenges. On the one hand, recent results (e.g., Backman et al. 2009; Ertel et al. 2011) reveal more and more the complexity of the known debris disks and the inability to reproduce the available data with simple models. Furthermore, standard assumptions (e.g., lower grain size consistent with the blow-out size of the system, grain size distribution exponent of -3.5) have been proven to be inadequate (e.g., Ertel et al. 2011; this work), and unexpected results like an outward increasing surface density distribution (Löhne et al. 2012) have been found. The increasing number and quality of the data available (in particular, spatially resolved images from Herschel and ALMA) allow one to break modeling degeneracies through simultaneous multi-wavelength modeling and to increase the complexity of the models. On the other hand, finding a
best-fit model in a complex parameter space results in significant challenges on the fitting method and the simulation of the observed data in order to reduce the computational effort necessary.

To rise to these challenges we, developed the tool SAnd. It is able to fit SED data and radial profiles in thermal reemission simultaneously using a simulated annealing approach (Press et al. 1992) on a grid of possible values in the parameter space. It uses an analytical approach for the radial density distribution of a rotationally symmetrical disk as well as for the grain size distribution. The SED is computed in an analogous way to Wolf & Hillenbrand (2005). Resolved images are computed with sufficient resolution and convolved with the telescope PSF. Radial profiles are then extracted along the major and/or minor axis. Each image (i.e., each pair of major and minor axis profiles or single major axis profile, if no minor axis profile is extracted, e.g., in azimuthally symmetrical disk images) is scaled by a factor $x_i$ to minimize the $\chi^2_{i,prof}$ (Eq. (A.3)). This way, only the shape of a radial profile is fitted. This is done because flux calibration uncertainties have to be considered only once per wavelength, but would be included for each data point in the radial profiles if absolute profiles were considered. Furthermore, deficits of our models to reproduce the absolute flux at one wavelength at which also radial profiles are fitted (e.g., due to uncertainties in the optical properties of the dust) would be considered in each profile data point and the SED data point which would result in an overweighting of the data at this wavelength. Simulated and observed multi-wavelength data are then compared using the following reduced $\chi^2$:

$$ \chi_r^2 = \frac{N_{\text{tot}} \left(\sum_{i=0}^{M} w_i\right)^{-1}}{N_{\text{tot}} - N_{\text{free}}} \left( w_0 \frac{\chi_{\text{SED}}^2}{N_0} + \sum_{i=1}^{M} w_i \frac{\chi_{i,\text{prof}}^2}{N_{i,1} + N_{i,2}} \right), \quad (A.1) $$

$$ \text{with} \quad \chi_{\text{SED}}^2 = \sum_{k=1}^{N_0} \left( \frac{F_k - \tilde{F}_k}{\sigma_k} \right)^2, \quad (A.2) $$

$$ \text{and} \quad \chi_{i,\text{prof}}^2 = \sum_{j=1}^{2} \sum_{k=1}^{N_{i,j}} \left( \frac{x_i S_{i,j,k} - \tilde{S}_{i,j,k}}{\sigma_{i,j,k}} \right)^2, \quad (A.3) $$

where $i=0$ is for the SED and $i \ge 1$ for the profiles, $j=1$ for major axis and $j=2$ for minor axis. In Eqs. (A.1) to (A.3), the notations used are:

- $M$ - number of images;

- $N_0$ - number of SED data points;

- $N_{i,1}$ - number of profile data points along major axis;

- $N_{i,2}$ - number of profile data points along minor axis;

- $N_{\text{tot}}$ - total number of data points, $N_0 + \sum_{i=1}^{M} \sum_{j=1}^{2} N_{i,j}$;

- $N_{\text{free}}$ - number of free parameters;

- $F_k, \tilde{F}_k$ - modeled and observed SED, respectively;

- $S_{i,j,k}, \tilde{S}_{i,j,k}$ - modeled and observed surface brightness profiles;

- $\sigma_k, \sigma_{i,j,k}$ - observed uncertainties;

- $x_i$ - surface brightness profile scaling factor;

- $w_i$ - weight.

For the present work, only one profile is used and $w_0$ and $w_1$ are set to 1.0, which means that the profile is given the same weight as the SED.

In order to allow one to find the best-fit model within a given range of parameters of a high dimensional parameter space in reasonable time, one has to compute the value of $\chi^2$ for a given set of parameters very fast. Creating the SED and three profiles

at different wavelengths for one particular model and comparing them with the observations usually takes ~0.1 s on an Intel Xeon E5410 CPU (2.3 GHz). The exact time depends mostly on the particular set of parameters and the number of wavelengths, at which radial profiles are provided. It hardly depends on the number of data points in the SED and the profiles.

The simulated annealing approach enables us to handle a very large and high dimensional parameter space without sampling the whole range of parameters. Therefore, a random walk on the grid of possible parameters is started at an arbitrary position. The probability $p$ to go a certain step depends on $\chi_1^2$ of the actual and $\chi_0^2$ of the previous step. It is 1 for $\chi_0^2 \ge \chi_1^2$ (fit becomes better) and follows a Boltzmann distribution for $\chi_0^2 < \chi_1^2$ (fit becomes worse):

$$ p = \exp \left[ - \frac{\chi_1^2 - \chi_0^2}{T} \right]. \quad (A.4) $$

In this equation, the quantity $T$ is a parameter comparable to the Temperature of an “annealing” physical system and controls the ability to reach areas in the parameter space that give worse fits. It is chosen at the beginning of each run to be a large value based on experiences from previous runs and is lowered with each successful step of the random walk. The width of each step is chosen from a probability distribution that prefers short steps, while the maximum step width is 1/2 of the whole range of values for each parameter. The probability to go long steps is lowered with each successful step. With this approach the code is able to reach each position in the parameter space, but prefers regions with low $\chi^2$. The run stops, when a position in the parameter space with the following properties is reached:

- with the actual value of $T$ the code was not able to leave the position after computing 500 further models;

- the actual $\chi^2$ is lower than or equal to the best reached $\chi^2$ during the whole run;

- the actual $\chi^2$ is lower than or equal to a maximum $\chi^2$ based on experiences from previous runs (input parameter).

A set of model parameters that satisfies all of these conditions is considered as a likely global fit. To increase confidence in this result, one can start several runs at different regions of the parameter space (4 in the recent case for the final fitting). Only if all runs reach the same region in the parameter space, we consider the result to be really a global best-fit. If only the first and not more than one other of these stopping criteria is satisfied, the value of $T$ is increased and the code goes on trying to leave the actual position in the parameter space.

Error estimates are done in the context of the simulated annealing approach by starting a new random walk at the best-fit position with a fixed value of $T$ (usually 10% of the best-fit $\chi^2$ in our case). The code counts how often each value of a parameter is reached. This way, the code samples the projection of the probability distribution of the random walk (Eq. (A.4)) on the axis of each parameter. With the knowledge of the value of $T$ one can then compute the probability distribution

$$ p' = \exp \left[ - (\chi^2 - \chi_{\text{best}}^2) \right] \quad (A.5) $$

and corresponding levels of confidence.

As a feature of the simulated annealing approach the time it takes to find a best-fit (and uncertainties) does not vary very much with the number of free parameters or the total size of the parameter space to be searched, but depends on the shape of the $\chi^2$ distribution in the parameter space. Weak variations
of the $\chi^2$ over a broad range of values of the parameters as well as very degenerate fitting problems require more time, i.e. more models to be computed. For the recent work (SED fitting only) we find the best-fit in a parameter space composed of ~$10^{10}$ to $10^{11}$ knots typically in less than 24 h including the estimates of the uncertainties. Therefore the code computes about 100.000 to 1.000.000 models and accepts (goes) about 10.000 and 100.000 steps for the fitting and computes ~1.000.000 models and accepts ~300.000 steps for the error estimates.

## References

Allende Prieto, C., & Lambert, D. L. 1999, A&A, 352, 555

Augereau, J. C., Lagrange, A. M., Mouillet, D., Papaloizou, J. C. B., & Grorod, P. A. 1999, A&A, 348, 557

Augereau, J.-C., Absil, O., Bouvier, J., et al. 2008, in SF2A-2008, ed. C. Charbonnel, F. Combes, & R. Samadi, 443

Aumann, H. H., Beichman, C. A., Gillett, F. C., et al. 1984, ApJ, 278, L23

Backman, D. E., & Paresce, F. 1993, in Protostars and Planets III, ed. E. H. Levy, & J. I. Lunine, 1253

Backman, D., Marengo, M., Stapelfeldt, K., et al. 2009, ApJ, 690, 1522

Beichman, C. A., Bryden, G., Stapelfeldt, K. R., et al. 2006, ApJ, 652, 1674

Brott, I., & Hauschildt, P. H. 2005, in The Three-Dimensional Universe with Gaia, ed. C. Turon, K. S. O'Flaherty, & M. A. C. Perryman, ESA Spec. Publ., 576, 565

Bryden, G., Beichman, C. A., Trilling, D. E., et al. 2006, ApJ, 636, 1098

Dermott, S. F., Jayaraman, S., Xu, Y. L., Gustafson, B. Å. S., & Liou, J. C. 1994, Nature, 369, 719

Dohnanyi, J. S. 1969, J. Geophys. Res., 74, 2531

Draine, B. T. 2003, ApJ, 598, 1017

Draine, B. T. 2006, ApJ, 636, 1114

Eiroa, C., Fedele, D., Maldonado, J., et al. 2010, A&A, 518, L131

Ertel, S., Wolf, S., Metchev, S., et al. 2011, A&A, 533, A132

Flower, P. J. 1996, ApJ, 469, 355

Fuhrmann, K. 2008, MNRAS, 384, 173

Garcés, A., Catalán, S., & Ribas, I. 2011, A&A, 531, A7

Gray, R. O., Corbally, C. J., Garrison, R. F., et al. 2006, AJ, 132, 161

Hillenbrand, L. A., Carpenter, J. M., Kim, J. S., et al. 2008, ApJ, 677, 630

Holmberg, J., Nordström, B., & Andersen, J. 2007, A&A, 475, 519

Holmberg, J., Nordström, B., & Andersen, J. 2009, A&A, 501, 941

Kharchenko, N. V., Scholz, R., Piskunov, A. E., Roeser, S., & Schilbach, E. 2007, VizieR Online Data Catalog, 3254

Krivov, A. V. 2010, Res. Astron. Astrophys., 10, 383

Krivov, A. V., Löhne, T., & Sremčević, M. 2006, A&A, 455, 509

Lebreton, J., Augereau, J.-C., Thi, W.-F., et al. 2012, A&A, 539, A17

Löhne, T., Augereau, J., Ertel, S., et al. 2012, A&A, 537, A110

López-Santiago, J., Montes, D., Crespo-Chacón, I., & Fernández-Figueroa, M. J. 2006, ApJ, 643, 1160

Luhman, K. L., Patten, B. M., Marengo, M., et al. 2007, ApJ, 654, 570

Maldonado, J., Martínez-Arnáiz, R. M., Eiroa, C., Montes, D., & Montesinos, B. 2010, A&A, 521, A12

Mamajek, E. E., & Hillenbrand, L. A. 2008, ApJ, 687, 1264

Marois, C., Zuckerman, B., Konopacky, Q. M., Macintosh, B., & Barman, T. 2010, Nature, 468, 1080

Martínez-Arnáiz, R., Maldonado, J., Montes, D., Eiroa, C., & Montesinos, B. 2010, A&A, 520, A79

Masana, E., Jordi, C., & Ribas, I. 2006, A&A, 450, 735

Matthews, B. C., Sibthorpe, B., Kennedy, G., et al. 2010, A&A, 518, L135

Messina, S., Pizzolato, N., Guinan, E. F., & Rodonò, M. 2003, A&A, 410, 671

Meyer, M. R., Backman, D. E., Weinberger, A. J., & Wyatt, M. C. 2007, Protostars and Planets V, 573

Ott, Stephan, Herschel Science Centre, European Space Agency, E. 2010, ASP Conf. Ser., 434, 139

Pan, M., & Schlichting, H. E. 2012, ApJ, 747, 113

Perryman, M. A. C., Lindegren, L., Kovalevsky, J., et al. 1997, A&A, 323, L49

Phillips, N. M., Greaves, J. S., Dent, W. R. F., et al. 2010, MNRAS, 403, 1089

Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1

Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 1992, Numerical recipes in FORTRAN, The art of scientific computing (Cambridge, UK: Cambridge University Press)

Reidemeister, M., Krivov, A. V., Stark, C. C., et al. 2011, A&A, 527, A57

Reiners, A. 2006, A&A, 446, 267

Roccatagliata, V., Henning, T., Wolf, S., et al. 2009, A&A, 497, 409

Rocha-Pinto, H. J., Flynn, C., Scalo, J., et al. 2004, A&A, 423, 517

Skiff, B. A. 2010, VizieR Online Data Catalog, 1, 2023

Smith, B. A., & Terrile, R. J. 1984, Science, 226, 1421

Stern, S. A. 1996, A&A, 310, 999

Teplitz, H. I., Chary, R., Elbaz, D., et al. 2011, AJ, 141, 1

Thébault, P., & Augereau, J. 2007, A&A, 472, 169

Trilling, D. E., Bryden, G., Beichman, C. A., et al. 2008, ApJ, 674, 1086

Valenti, J. A., & Fischer, D. A. 2005, VizieR Online Data Catalog, 215, 90141

van Leeuwen, F. 2007, A&A, 474, 653

Vitense, C., Krivov, A. V., & Löhne, T. 2010, A&A, 520, A32

Voshchinnikov, N. V., Il'in, V. B., Henning, T., & Dubkova, D. N. 2006, A&A, 445, 167

Wolf, S., & Hillenbrand, L. A. 2003, ApJ, 596, 603

Wolf, S., & Hillenbrand, L. A. 2005, Comp. Phys. Comm., 171, 208

Wyatt, M. C. 2005, A&A, 433, 1007

Wyatt, M. C. 2006, ApJ, 639, 1153

Wyatt, M. C. 2008, ARA&A, 46, 339

Wyatt, M. C., Clarke, C. J., & Booth, M. 2011, Cel. Mech. Dyn. Astron., 39