Issue 
A&A
Volume 631, November 2019



Article Number  A141  
Number of page(s)  12  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/201935199  
Published online  08 November 2019 
Twisted debris: how differential secular perturbations shape debris disks
Astrophysikalisches Institut, FriedrichSchillerUniversität Jena,
Schillergäßchen 2–3,
07745
Jena,
Germany
email: jan.sende@unijena.de; torsten.loehne@unijena.de
Received:
4
February
2019
Accepted:
13
September
2019
Context. Resolved images suggest that asymmetric structures are a common feature of cold debris disks. While planets close to these disks are rarely detected, their hidden presence and gravitational perturbations provide plausible explanations for some of these features.
Aims. To put constraints on the properties of yet undetected planetary companions, we aim to predict what features such a planet imprints in debris disks undergoing continuous collisional evolution.
Methods. We discuss the basic equations, analytic approximations and timescales governing collisions, radiation pressure and secular perturbations. In addition, we combine our numerical model of the collisional evolution of the size and spatial distributions in debris disks with the gravitational perturbation by a single planet.
Results. We find that the distributions of orbital elements in the disks are strongly dependent on grain sizes. Secular precession is differential with respect to involved semimajor axes and grain sizes. This leads to observable differences between the big grains tracing the parent belt and the small grains in the trailing halo. Observations at different wavelengths can be used to constrain the properties of a possible planet.
Key words: circumstellar matter / planetary systems / planetdisk interactions / celestial mechanics / methods: analytical / methods: numerical
© ESO 2019
1 Introduction
Despite its relatively low total mass, dust in circumstellar debris disks is readily detected because of its high specific cross section. Current farinfrared detection rates are as high as ~20% for nearby stars (Eiroa et al. 2013; Montesinos et al. 2016; Sibthorpe et al. 2018), indicating how ubiquitous these disks are. Akin to the EdgeworthKuiper belt in the Solar system, debris disks are often found at distances of tens or more than one hundred au from their host stars (e.g. Rhee et al. 2007; Pawellek et al. 2014; Matrà et al. 2018). At these distances, they are typically the only detectable components of planetary systems. While thousands of planets have already been detected close to the stars, detections near the cold outer disks are scarce and limited to young, highmass planets and substellar companions (Maldonado et al. 2012). The old Jupiters and Neptunes that dominate the outer Solar System remain hidden around other stars. The most promising way to characterise these is indirect: through their gravitational perturbations of debris disks.
Gravitational perturbation can lead to substantial orbit changes, even far from the perturber. Meanmotion resonances, typically combined with drifting dust or migrating planets, can produce observable clumps (Kuchner & Holman 2003; Wyatt 2003; Ozernoy et al. 2000). Secular perturbation can modify orbital eccentricities and orientations (Murray & Dermott 2000). The latter is a potential explanation for eccentric, offset disks such as those imaged around Fomalhaut (Kalas et al. 2005; Stapelfeldt et al. 2004; Boley et al. 2012; MacGregor et al. 2017), HR 4796A (Schneider et al. 1999; Moerchen et al. 2011; Kennedy et al. 2018), and HD 202628 (Schneider et al. 2016).
The pure gravitational dynamics are well understood. Analytic approximations in the Laplace–Lagrange formalism describe the energyconserving orbital precession caused by mild perturbations (see, e.g. Murray & Dermott 2000). When an orbit is perturbed by a single planet or companion, eccentricity and longitude of pericentre, as well as inclination and longitude of the ascending node oscillate at one frequency. This frequency is determined by the perturber’s distance and mass. The perturber modulates shape and orientation of the orbit, without changing the orbital period.
Observed belt eccentricities and widths can be related to orbit parameters and masses of unseen perturbers (e.g. Rodigas et al. 2014; Thilliez & Maddison 2016), but this approach provides no unique solutions. There is always a degeneracy between perturber orbit orientations, eccentricities, and distances. Due to typical evolution timescales of several Myr it is impossible to observe more than a snapshot of a given system. Further degeneracy comes from the fact that the observed narrowness of some belts can have different reasons. On the one hand, a belt could already be equilibrated, corresponding to orbits that are randomly distributed about the precession cycle. In that case, the precession amplitude, that is, the free elements, must be small. On the other hand, precession could have started recently, being still far from equilibrium. These two cases are discussed in Löhne et al. (2017). The nonequilibrium case is what we focus on throughout the following sections. Figure 1 illustrates the resulting degeneracy for a system with a narrow belt.
What ismore, the observed dust is always affected by radiation pressure. Small dust grains, which dominate the surface brightness at short wavelengths, are placed on wider orbits by stellar radiation pressure and populate a halo that surrounds the belt of bigger parent bodies (Strubbe & Chiang 2006; Krivov et al. 2006; Thébault & Wu 2008). Depending on where in an eccentric belt these small grains are released, they may stay bound or become unbound, stay closer to the parent belt or populate wider orbits. The asymmetry in the parent belt translates to an asymmetry in the halo (Lee & Chiang 2016).
The collisional cascade changes the picture again. In our recent work on these eccentric disks (Löhne et al. 2017), we could show how collisions smoothen or remove observable structure. The parent belt appears broadened at short wavelengths, the halo more symmetric. Kim et al. (2018) show that the amplitude of the azimuthal flux asymmetry in eccentric disks depends on the shape of the dust size distribution, which in turn is affected by details of the underlying collision physics. This dependence on material strength and dynamical excitation adds further degeneracy.
In this paper we show and discuss how differential precession shapes belts and halos. Precession timescales depend on the distance between perturber and perturbee, and even on the particles’ susceptibility to radiation pressure. The timescales can be very different across a disk, resulting in an effective sheer. This sheer can help reduce the degeneracy and constrain the perturber’s orbit. The blue halos in Fig. 1 illustrate this idea.
Section 2 repeats the key points of secular perturbation theory, and amends the formulas to include the dependence on particle size. Furthermore, differential precession is explained, and the relevant timescales are discussed. In Sect. 3 the implementation of perturbations into our existing code ACE (analysis of collisional evolution) is described, and the occurring numerical dispersion is discussed. Moreover, the results of a numerical simulation are shown and explained using differential precession. In Sect. 4 the relevance and implications of our results are discussed. The findings are summarised in Sect. 5.
Fig. 1 Threedifferent interior perturbers in red that can produce the same eccentric belt (thick black) from an initially circular belt. Orbits that deviate from the belt will be sheared by differential precession, as exemplified by orbits of older halo particles (lighter blue) trailing that of freshly produced ones (darker blue). The location of the star is marked by grey cross hairs. The perturber’s periapsis is indicated by red lines. 

Open with DEXTER 
2 Differential secular perturbations
Secular perturbation theory describes the orbital parameter evolution of bodies under mutual gravitational influence, orbiting the same central object. In the following an abridged version of secular perturbation theory is presented. Corrections are made to include the influence of radiation pressure, which introduces a size dependence for smaller particles. A model describing a differentiation in the orbital evolution of particles created by collisions in the debris disk is derived and put into context.
2.1 Secular perturbation theory
While eccentricity e, longitude of pericentre ϖ, inclination i and longitude of the ascending node Ω can change, the semimajor axis a is not modified by secular perturbations. In the following, we focus on the evolution of e and ϖ^{1}, giving a short recapitulation of Murray & Dermott (2000).
In order to simplify the equations, the eccentricities e and the longitudes of pericentres ϖ of all bodies are transformed into the h–k coordinate system: (1)
This transformation simplifies the maths involved and helps interpreting the results. In the h–k space, the complicated changes of e and ϖ become operations on the vectors (h, k).
We are interested in the case of one planet and a lowmass debris disk orbiting around a star. If the planet is much more massive than the debris disk, only the orbital elements of the debris are significantly changed. The orbital elements of the planet are assumed constant because the gravitational influence of the disk on the planet’s orbit is small. However, in a multiplanet system the mutual gravitational influence of the planets would lead to a change in orbital elements over time.
In the case of a single planet, only the differential equations of the debris need to be solved: (2)
The constants A and A_{p} are calculated using the mean motion ^{2}, the mass ratio ɛ = m_{p}∕M_{⋆}, the semimajor axis functions α and , and the Laplace coefficients and : (3)
The semimajor axis functions depend on the position of the debris relative to the planet: if the semimajor axis of the planet a_{p} is smaller than the axis of the debris a we have the case of an inner perturber (a_{p} < a). If the semimajor axis of the planet a_{p} is larger than the axis of the debris a, we have the case of an outer perturber (a_{p} > a). If the semimajor axes of the planet and the debris are similar (a_{p} ≈ a), secular perturbation theory fails to provide a good approximation: (4)
Inserting A and A_{p} from Eq. (3) into Eq. (2), the general solution for the orbit evolution is: (5)
The solution can be interpreted as the orbit of the debris precessing with a constant angular velocity around the forced eccentricity vector (h_{forced}, k_{forced}) = (e_{forced} sinϖ_{p}, e_{forced} cosϖ_{p}) in a constant distance e_{free}.
The forced eccentricity e_{forced} and the precession timescale depend on the planet’s eccentricity e_{p}, and the semimajor axis functions α and : (6)
The free eccentricity e_{free} and the free angle ϕ_{free} are determined from the initial conditions. They can be calculated using the starting position in the h–k space, and inserting t = 0 into Eq. (5), solving for them. If the perturber is far away from the disk (α ≪ 1), the Laplace coefficients can be series expanded and reduced to their leadingorder terms: , and (Murray & Dermott 2000). The forced eccentricity and the precession timescale thus are (7)
In a resolved, extended disk, α can be used to distinguish between possible planets because α increases when the semimajor axis of the debris a is closer to the planets axis a_{p}. Combining Eq. (4) with Eq. (7), one gets the following proportionalities: (8)
For an inner perturber, debris with larger a has a smaller forced eccentricity e_{forced} and a longer precession timescale . The dependency is inverse for an outer perturber – debris with larger a has a larger e_{forced} and a shorter . In both cases, the debris with a semimajor axis closer to the planets axis precesses faster. In the case of an inner perturber, the inner disk edge leads in front of the outer edge, whereas in the case of an outer perturber it drags behind. Furthermore, scales with α^{−3} instead of α^{−2}, leading to aslower precession for an outer planet, compared to an inner planet with the same α. In the following we focus on the case of an inner perturber, but it is possible to derive similar formulas for the outer one.
2.2 Radiation pressure
In the above it is assumed that only gravitational forces act on the debris. In reality however, small dust particles can be heavily influenced by radiation pressure and stellar winds. The above equations are amended, using the Burns et al. (1979) formalism: F = F_{grav} − F_{rad} = (1 − β)F_{grav}, where β is the ratio between radiative push and gravitational pull. (For our analytical model, we ignore stellar winds.) Applying the changes to Eq. (7), it is easy to check, that the forced eccentricity does not change, but the precession timescale now depends on the β value of the dust as well: (9)
Dust with the same semimajor axis (and thus same α), but higher β value precesses faster in the h–k phase space than dust with lower β. Thus, if one knows the position and β value of particles in the h–k space, one can calculate e_{forced} and .
The β value depends on the material and size of the particle, and the luminosity and spectrum of the host star (Burns et al. 1979): (10)
with L_{⋆} and M_{⋆} being the luminosity and the mass of the host star, ρ the dust density,Q_{PR}(s) the radiation pressure efficiency, and s the size of the particle. Usually Q_{PR}(s) is calculated using Mie theory, and depends on the stellar spectral type and the material used.
2.3 Dust release
The detectable dust of a debris disk is created in mutual collisions of bigger objects, and in cometary activity. These source objects are usually not detectable because of their small specific surface area. Yet, the small dust particles they release, have a much larger specific surface area and are therefore detected.
The orbits of fragments released by collisions or cometary activity depend on several parameters, such as the masses of the colliders, the impact velocity and location, etc. However, if a fragment is released from a larger parent body on a roughly circular orbit, the orbit parameters of the fragment only depend on its β value and the semimajor axis of the source a_{src} (Strubbe & Chiang 2006, Eq. (9)): (11)
While these equations are derived for fragment production in a circular parent belt, they still hold for parent belt eccentricities that are small compared to the β values of small halo grains, that is, for e_{src} ≪ β.
As described in Sect. 2.2, the β value depends on the size of the particle, therefore the orbits of collisional fragments depend on the sizes as well. As shown on the left in Fig. 2, larger fragments, which usually have a β close to zero, will stay in nearly circular orbits within the parent belt. Smaller fragments however, due to their larger β, are on elliptic orbits and create the halo of the disk. Particles close to β ≈ 1∕2 are on highly eccentric orbits. They constitute the outermost part of the halo.
When launched from circular orbits, particles with β ≳ 1∕2 have an eccentricity greater than one, and are therefore on unbound orbits, getting ejected from the system. (Due to their anticlockwise motion in Fig. 2, unbound fragments can only be ejected into the red shaded region.) Using Eq. (10) one finds a minimum particle size for each star, below which particles get unbound. This size is called the blowout size s_{blow}. For early type stars it is greater than for late type stars: (12)
where Q_{PR} is the radiation pressure efficiency (Burns et al. 1979). The radiation pressure of low luminosity stars is not strong enough, and thus no small size limit exists (Reidemeister et al. 2011).
2.4 Differential precession
Combining the dependence of Eq. (11) on β with the forced eccentricity and precession timescale of Eq. (9), one gets the full βdependent phase space evolution of fragments produced in a circular belt: (13)
where a_{belt} is the semimajor axis of the initial debris disk belt. According to Eq. (11), fragments of different size will be produced at different locations in the h–k space. These then precess at different speeds around different forced eccentricities, as described in Eq. (13).
Figure 2 illustrates the effect on particles with different β values, produced in the same collision. Small particles, which have higher β values than large particles, begin their dynamic evolution on more eccentric orbits. They start farther from their forced elements in the h–k space and therefore move around wider circles. Small particles have larger semimajor axes as well, leading to lower forced eccentricities and slower precession rates. Although all particles had the same time to precess, large belt particles with β ≈ 0 covered a much longer arc in the h–k space than the small halo particles.
While Fig. 2 shows the case of particles created in a single collision, Fig. 3 shows how the ensemble of collisional fragments created at different longitudes in the parent belt would evolve. For simplicity, the parent belt is assumed to be an azimuthally uniform, narrow, circular ring. The collisional fragments are initially distributed on circles in the h–k space (blackand blue). In realistic belts with finite widths and dynamical excitations, these rings would have a certain width, too.
The parent belt has the smallest semimajor axis, and therefore precesses the fastest in the h–k space. Largeparticles, with low β values, stay close to the parent belt because their eccentricity spread is small. The smaller the grains, the higher are the β values, and the wider is their eccentricity spread. The small particles are on orbits with large semimajor axes and form the disk halo. Their forced eccentricity is low and they precess at a slower rate compared to the parent belt and the large grains.
Figure 3 additionally shows how the phase space evolution is highly dependent on the planets’ semimajor axes. Remembering Eq. (8), we see that for both the planet in the right panel and the planet in the left panel, the forced eccentricity is the same. The precession time of the disks however differs significantly (left: 19.35 Myr; right: 2.15 Myr). After the same evolution time of T = 1.29 Myr both disks look very different. However, with the time T usually unknown, the two cases cannot be distinguished based on their secular precession alone. After T = 11.61 Myr, the distribution on the left would be identical to the one on the right at T = 1.29 Myr. It is only through the combined action with PoyntingRobertson drag and collisional evolution, acting on different timescales, that this degeneracy is broken.
Fig. 2 Left: initial orbits of fragments with β values of 0.0 (orange), 0.1 (black), and 0.3 (blue), emitted in a single collision from a circular orbit at a = 100 au (and ϖ = 30°; particles moving anticlockwise). Particles with β ≥ 0.5 (red) are on unbound orbits. Right: initial positions (circle) of the same fragments in the h–k space and their circles of movement (dashed line) around the forced eccentricity (dot of according colour). The arcs they covered in the h–k space within 5 Myr are shown as solid lines. Given the planet at a_{p} = 20 au with e_{p} = 0.6 (and ϖ_{p} = 0°), they have precession times of (for β = 0.0), (for β = 0.1), and (for β = 0.3). 

Open with DEXTER 
2.5 Comparison of timescales
In addition to the precession timescale, several other timescales need to be considered for a debris disk: the crossing timescale, the timescale of Poynting–Robertson drag, and the average collision time of particles. In order to make comparisons easier, thesetimescales are expressed in terms of the physical properties of the host star, the planet, and the debris disk.
The crossing timescale was derived by Mustill & Wyatt (2009). They showed that particles which are initially separated by different semimajor axes, are forced on crossing orbits by secular perturbations. This leads to an increased destruction of particles in mutual collisions. The time it takes the debris to get on crossing orbits, is the crossing timescale: (14)
Furthermore, we derive a similar form of the Poynting–Robertson timescale for particles created in the initial debris disk. This is the time it takes a dust particle to leave the its initial pericentre position due to Poynting–Robertson drag (15)
In order to allow a better comparison, we rewrite the precession timescale in the same manner: (16)
Lastly, we need the collisional timescales. This is the time it takes a particle on average to be destroyed in another mutual collision. In contrast to the other timescales, the collisional timescale directly depends on the mass of the debris disk. We did not derive an analytic expression, but use the collision timescale calculated in our simulation presented in Sect. 3 instead. (For possible analytical solutions see Wyatt 2005a and Löhne et al. 2008.)
In our later simulation (Sect. 3), we are using an A4Vstar. The relevant timescales are shown in Fig. 4. The dust mass is m_{Dust} = 4.05 × 10^{25}g for particles smaller than1 mm, at a time of 5.75 Myr. The blowout size is s_{blow} = 4 μm. For sizes greaterthan 100 μm, β is nearly zero (β < 0.04). These particles roughly follow the motion of the parent belt^{3}. Because of the small β value, the precession timescale is constant for a wide range of particle sizes. However, near the blowout size, it increases significantly with smaller particle size because these particles have larger semimajor axes. The Poynting–Robertson timescale decreases with smaller particle sizes because the influence of radiation pressure gets larger with smaller particle size. Only for the smallest grains, near the blowout limit, the timescale increases again because of the extremely eccentric, extended orbits the particles areon.
For theparticle sizes shown, the collision timescale is always shorter than the precession timescale. Particles in this size range precess around the forced eccentricity in the h–k space, but will collide with other particles before they can complete a total precession period. In contrast, the larger objects of the debris disk have much longer collisional lifetimes (for particles of 1m already T_{Col} ≥ 10 Myr). They will complete one or more precession periods in the h–k space before they are destroyed in collisions. That is, the larger debris in the parent belt precesses and simultaneously creates a halo of dust, which precesses more slowly, with the precession time depending on the β value, and thus the particle size.
Figure 4 shows the crossing timescale for particles with β ≈ 0 (otherwise it does not apply). For these particles the crossing timescale is shorter than the precession timescale. Thus, large belt particles with β ≈ 0 on noncrossing orbits with marginally different semimajor axis will evolve into a crossing configuration before they can complete a whole precession period, meaning that their collision rate increases. This effect might be significant and even trigger the collisional cascade (Mustill & Wyatt 2009). Yet, in our later numerical example (see Sect. 3.4) we assume a prestirred disk, where it is less important.
The Poynting–Robertson timescale is mostly longer than the precession and crossing timescale, meaning that Poynting–Robertson drag does not play an important role. For particles smaller than 10 μm, the Poynting–Robertson timescale is shorter than the precession timescale, meaning that these particles fall into the star due to Poynting–Robertson drag before they can complete a rotation in the h–k space. Again, in our later example, this effect is not important because the collision timescale is even smaller.
The interplay between these timescales is dependent on the position of the planet a_{p} and the debris disk belt a_{belt}, in particular because the Poynting–Robertson and the collisional timescale only depend on the belt alone. For example, if the belt is placed at 50 au instead of 100 au (other parameters are the same as above), the timescales change as follows: (17)
All these timescales get smaller, however the ratio is much smaller for the crossing and the precession timescales, leading to much more pronounced differential precession.
Fig. 3 h–k space evolution of parent belt particles (orange; a_{belt} = 100 au; β = 0.0; e ≈ 0.0), and particlescreated in collisions of parent belt objects; with fragment β values of β = 0.01;0.03;0.05 (black), and 0.3 (blue). The belt is around a star of M = 1.92 M_{⊙}, and the collisional fragments are perturbed by a m_{p} = 2.5 M_{Jup} mass planet, with a = 20 au and e = 0.6 (left); a = 60 au and e = 0.2 (right). The initial distribution of particles created in the collisions (dashed line), and the position after a time T (solid line) are shown. In between steps are marked with faint colours. 

Open with DEXTER 
Fig. 4 Relevant timescales for fragments created by collisions inside a debris disk (a = 95– 105 au; m_{Dust} = 4.05 × 10^{25}g) around an A4Vstar (M_{⋆} = 1.92 M_{⊙}) with a planet at a_{p} = 20 au with an eccentricity of e_{p} = 0.6 and a mass of m_{p} = 2.5 M_{J}. The β value of the particles (blue, dotted), and the blowout limit (blue zone on the left) are shown. The precession timescale (orange; top), the crossing timescale (black; middle; where applicable), and the Poynting–Robertson timescale (red; crossing the top and the middle plot) are derived analytically. The collision timescale (black, dashed) was taken directly from the numerical results of the simulation. When possible, a range (shaded area) of timescales for a = 95– 105 au is given in addition to the value at a = 100 au (solid line). 

Open with DEXTER 
3 Combining perturbations and collisions
The model presented above provides estimates for individual pieces of the physics. Yet, for an analysis of the interplay between collisions and perturbations a numerical approach is needed. For that purpose we use our existing kinetic code ACE (Krivov et al. 2005, 2006; Reidemeister et al. 2011; Löhne et al. 2017). In ACE the debris disk evolution is described by a fourdimensional phase space distribution. The dimensions are: the masses of the individual particles m, the pericentre distances q, the eccentricities e, and the longitudes of periapsis ϖ. We use the pericentre distance instead of the semimajor axis because the belt and the halo have similar pericentre distances, whereas the semimajor axes differ widely. In order to combine collisions with secular perturbations, we now add the gravitational influence of a single perturbing planet.
3.1 Implementation in ACE
ACE calculates secular perturbations using the transport equation. The phase space density n(x) at the position x = (m, q, e, ϖ) changes depending on the velocity at this point: (18)
For each grid point x in the fourdimensional phase space, the derivatives ė(x) and are calculated using the derivative of Eq. (1): (19)
The derivatives ḣ(x) and are calculated using Eq. (2), thus no additional analytical approximation is necessary. A series expansion of the Laplace coefficients is used.
The derivative of the pericentre distance needed in Eq. (18) is (20)
The contributions to ė, , and ȧ are summed up for secular precession and Poynting–Robertson drag, with ȧ = 0 for the former and for the latter. The transport equation itself is approximated using a variant of the upwind differencing scheme (Press et al. 2007). The parameter velocity v(x) is calculatedfor each phase space bin, and the material shifted into the neighbouring bins depending on the velocity direction.
A sketch of the algorithm is shown in Fig. 5. ACE uses a fourdimensional phase space, however for simplicity we only show a twodimensional sketch of the pericentre q and eccentricity e dimensions. According to our differencing scheme, bin n_{i,j} loses material to its neighbours and gains some from them. The amount lost depends on the parameter velocity v_{i,j}. The material lost in n_{i,j} is moved into its direct neighbours along the q and e axis. In the example, these are the n_{i,j+1} and the n_{i−1,j} bin.
The total phase space change for the n_{i,j} bin is: (21)
where Δe and Δq are the spacing of the phase space grid, and e_{e} and e_{q} are the unit vectors in e and q direction.
One shortcoming of this firstorder method is that no material is moved directly along the direction that the parameter velocity v_{i,j} is actually pointing to. In the example in Fig. 5, no material moves directly into bin n_{i−1,j+1}. However, if we assume that the velocity is similar in the neighbouring bins, the discrepancy is compensated by material moving from bins n_{i,j+1} and n_{i−1,j} into n_{i−1,j+1}. Accuracy is limited by finite bin widths and parameter velocities varying within the grid.
The scheme used in ACE introduces an artificial numerical dispersion into the calculation. Looking at the movement in only the ϖ dimension, and assuming the parameter velocity to be constant, one gets the following differential equation: (22)
Following the notation of Press et al. (2007), and using our differentiation scheme, one gets: (23)
This equation can be rearranged by adding : (24)
Comparing the result with the one derived by Press et al. (2007) for the ForwardTimeCentreSpace scheme, one sees that the discretised equation is the same as the following differential equation: (25)
The equation is the same as Eq. (22), but introduces an additional dispersion term with D = vΔϖ∕2. The dispersion term depends on the grid spacing Δϖ and the parameter velocity v, but is independent of the integration time step. Therefore, the initial phase space distribution always disperses over time.
Fig. 5 Sketch of the upwind scheme, used inside ACE, acting on the pericentre q and eccentricity e dimensions of the fourdimensional phase space. The material is shifted from each bin n_{i,j} into its direct neighbours by the parameter velocity v_{i,j} (calculated with perturbation theory). 

Open with DEXTER 
3.2 Collisionless test case
As an implementation check, we ran a simulation without collisions, but with an inner planet acting on the debris disk. The planet was set on a circular orbit (e_{p} = 0) at a_{p} = 20 au. The debris disk was set between a = 95 au and 110 au, with the eccentricity being distributed in a radius of 0.1 around e = 0.3, ϖ = 0.0 in the h–k space. Inserting the data into Eq. (7), we see that the planet produces a precession of debris in the h–k phase space around e_{forced} = 0. This precession does not modify the eccentricities of the orbits, only the longitudes of pericentre change over time. The precession frequency depends on the slight differences in semimajor axis and on the differences of the β values of the grains. Without collisions, no new particles with different orbital configurations are created and Eq. (13) does not apply. The initial distribution is only moved around. Figure 6 shows a n–ϖ plot of the phase space distribution of objects with s = 52 km and a semimajor axis of a = 95 au and 110 au. We chose the size of s = 52 km to make sure that the differences we see only depend on the semimajor axis and not on the β value. For both distances the distribution moves over time, but because of the numerical dispersion explained in Eq. (25), it gets more spread out along ϖ over time as well.
The simulated distribution was fitted to Eq. (52) of Neuman (1981). The fitted parameter velocities match the analytical expectations. For a semimajor axis of a = 110 au, a parameter velocity of v = 17.94 deg Myr^{−1} was predicted, and a velocity of v = 17.99 deg Myr^{−1} was measured. For a = 95 au a velocity of v = 29.97 deg Myr^{−1} was predicted and v = 31.10 deg Myr^{−1} measured. The same was repeated for the dispersion factor D, using the estimate of D = vΔϖ∕2 from the previous section and the grid spacing of Δϖ = 10 deg. In the case of a = 110 au, a dispersion of D = 89.7 deg^{2} Myr^{−1} was predicted and a D = 95.3 deg^{2} Myr^{−1} was measured. At a = 95 au, D = 149.8 deg^{2} Myr^{−1} was predicted and D = 152.1 deg^{2} Myr^{−1} measured.
These results show that the simulation behaved as expected. The dispersion of the phase space distribution was exaggerated by additional numerical errors. The slight difference in numerical value of v and D were not further investigated, but are probably due to our usage of the upwind scheme, which loses some accuracy in favour of fidelity (Press et al. 2007). The grid spacing of Δϖ = 10 deg was kept forour later example as a compromise between accuracy and simulation runtime. For the eccentricity grid spacing, which was not important in the implementation test, the same rule holds. The phase space distribution disperses with D = vΔe∕2 due to our integration scheme. However, the grid spacing is not linear in order to represent the orbits of particles created in collisions as well as possible. Therefore, the dispersion will depend on the position in the eccentricity grid. The dispersion could be reduced by adding more points to the grid, thereby reducing Δϖ, or by making the now constant grid follow the perturbed motion through e–ϖ space. Both approaches would significantly increase the simulation runtime.
Albeit numerical dispersion limiting the accuracy of the results, the main effects of secular perturbations combined with collisions remained visible in our simulation. The differential secular perturbations caused by differences in semimajor axes and β values equalled or exceeded the numerical dispersion. At the given times, the ϖ distributions at the inner and outer belt edges in Fig. 6 were as far apart as they were wide, which can be considered as a minimum requirement. A comparison between particles with low and high β values will yield larger differences in semimajor axis and therefore larger differences in precession velocity and larger separations in ϖ as well.
Fig. 6 Azimuthal diffusion in the ACE simulation described in Sect. 3.2. The dust distribution of particles at a = 95 au (solid; crosses) and a = 110 au (dashed; circles) along the ϖ dimension is plotted for a simulation time of 0.00 Myr (red), 2.87 Myr (blue), 5.75 Myr (orange), and 8.62 Myr (black). A fit to the analytical solution was made to extract the different parameter velocities and dispersion coefficients. The initial distribution is scaled down by a factor of 5 to fit the plot. 

Open with DEXTER 
3.3 Model parameters
We ran a full collisional simulation including perturbations by a single planet, and Poynting–Robertson drag. The planet orbited a A4V star, with a mass, luminosity, and effective temperature of M = 1.92 M_{⊙}, L = 16.6 L_{⊙}, and T_{eff} = 8600 K. The debris disk extended from a = 95 − 105 au, and had a moderate initial total mass of m_{Disk} = 2 M_{⊕} (cf. Krivov et al. 2018). In contrast to the test case in Sect. 3.2, the debris was initially distributed on nearly circular orbits, spread between eccentricities e = 0 and e = 0.1. We assumed an initial mass distribution index of q = 1.888, a minimum particle radius of s_{min} = 100 μm and a maximum radius of s_{max} = 49 km. The perturbing planet had a mass, eccentricity and semimajor axis of m_{p} = 2.5 M_{Jup}, e_{p} = 0.6, a_{p} = 20 au, causing a forced eccentricity e_{forced} = 0.15 at a = 100 au.
The assumed material was a homogeneous mixture (Bruggeman 1935) of Astrosilicate (Draine 2003) and water ice (Li & Greenberg 1998) in equal volume fractions, with a bulk density of 2.35g cm^{−3}. In order to calculate the β values of the particles, and generate the emission images of the debris disk, the nearest point in the PHOENIX/NextGen grid was used (Hauschildt et al. 1999).
The eccentricity grid of the simulation span [0.01, 1.5] and was logarithmically spaced for e ≤ 0.2, and for 1 < e, with at most a factor of Δe∕e ≈ 0.5 between neighbouring bins. For 0.2 < e < 1, the step size was limited to Δe ≤ 0.1. The transition between the linear and logarithmic regimes are smooth. For a more detailed description of the grid generation, see Eq. (29) in Löhne et al. (2017). The orientation grid had a constant spacing of Δϖ = 10 deg. The pericentre grid had 60 logarithmically spaced grid points, spanning from 40 to 160 au. The mass grid span a range from 7.14 × 10^{−13} to 1.13 × 10^{21}g. Objects unaffected by radiation pressure had a logarithmic spacing of Δm∕m ≈ 12. Small dust grains had a logarithmic spacing of Δm∕m ≈ 3. The transition between the two regimes is smooth, modelled with the number of bins per unit logarithm of mass being (26)
3.4 Numerical results
Figure 7 shows how the initial distribution has evolved after T = 5.75 Myr. The two top rows show the h–k space and the q–e space for grains in bins centred around radii of ≈5, ≈ 12, and ≈ 226 μm. (Single mass bins encompassing a size range are shown.) In the initial distribution grains down to s_{min} = 100 μm were generated. Grains of size ≈5 μm and ≈12 μm are not partof the initial distribution, but were created by collisions within the disk. Particles of size ≈ 226 μm have been initially present, but as shown in Fig. 4, their collisional lifetime is very short. They did not survive until the simulation time of T = 5.75 Myr, but were replenished by collisions of larger bodies.
Grains of size 5, 12, and 226 μm have β values of 0.43, 0.18, and 0.009. According to Eq. (16) the precession timescales are 1990, 42, and 19 Myr (at a distance of a_{belt} = 100 au). Thus, it is clear that at the current simulation time of T = 5.75 Myr no debris in the disk has completely precessed around the forced eccentricity. It is still the initial stage of precession.
The top row of Fig. 7 matches the theory of differential secular perturbations. Grains of all size ranges precess counterclockwise around the forced eccentricity. The ≈ 5 and ≈ 12 μm grains show precession differential to the ≈226 μm grains, which follow the movement of the parent belt. Due to their large β values, small grains are distributed in a larger area of the phase space.
The distribution of ≈5 μm grains seems counter intuitive. According to our theoretical predictions, grains should be distributed in a large circle at e ≲ 1. However, this is only true under the assumption of collisions in a circular belt. Barely bound grains produced by an already elliptic belt are not uniformly distributed in ϖ. Collisions at the pericentre of the disk are more likely to produce unbound grains than collisions at the apocentre, explaining the excess of particles in the bottom right. Furthermore, bins near e = 1 are wider, containing more mass per bin, which results in an apparent density increase at the bottom right. (See Fig. 3 in Löhne et al. 2017 as well.)
In addition, for some of the ≈5 μm grains, the analytic prediction suggests that they should reach eccentricities of e > 1. This is a shortcoming of the linear approximation of the Laplace–Lagrange theory. Secular perturbations cannot make a particle cross the e = 1 boundary because bringing particles on unbound orbits would require energy transfer. (Which is forbidden by ȧ = 0.) In ACE simulations this is handled by setting ė ≤ 0 for the bin closest to e = 1, while still allowing arbitrary . This handling of the boundaries reduces the fidelity of the results for very eccentric grains. However, the phase space asymmetry with respect to the line running from the top left to the bottom right is no pure artefact as these grains have large semimajor axes, long precession timescales and lifetimes, and hence, move slowly through the phase space, allowing pileup near e = 1. A more detailed analysis would require a finer eccentricity grid in this region of the phase space.
The phase space evolution discussed in Sects. 2.4 and 2.5 does not take into account that the initial grains are destroyed over time. Yet, although these particles are destroyed in collisions, new particles are created by other collisions as well. At any moment in time, a new starting distribution of particles is created, which shows differential precession over time. An image of a disk is always a superposition of different, initial distributions. As these initial distributions evolve over time, they occupy different parts of the h–k space, effectively looking like a single, smeared out distribution. At the same time, dust of older initial distributions is lost in collisions, making them fade in intensity. Therefore, observations always show a smeared out snapshot of the most recent initial distributions. In order to extract the correct planet parameters from an image, one needs to simulate the whole phase space evolution to get the smearing right. However, in an approximation, the model described in Sects. 2.4 and 2.5 can still be used.
The second row of Fig. 7 shows the q–e distribution for the same particle sizes. As explained by Eq. (11), particles with a lower β are created on more eccentric orbits, therefore the smaller particles have higher eccentricities. However, the eccentricity is additionally increased by the secular perturbations of the planet, leading to higher eccentricities even for the parent belt: most ≈ 226 μm grains have reached e > 0.2, instead of staying at the initiale = 0 to e = 0.1. This change in eccentricity is reflected in the pericentre distribution as well. With the belts semimajor axis starting between 95 and 105 au, the pericentres have now evolved to approximately 76–84 au.
The same is true for the numerical diffusion. Due to the coupling of pericentre distance and eccentricity, a diffusion in eccentricity leads to a diffusion in pericentre distance as well. This limits the simulation time, as for long enough timescales, the differences in semimajor axis, which cause differential precession, are lost in the diffusion of pericentre distances.
The third row shows the geometrical normal optical thickness (hereafter referred to as optical thickness) of the disk for certain grain size ranges. The optical thickness is the amount of material in the sight line, and while not observable directly, it provides a good illustration of the spatial distribution of material in the disk.
Except for barely bound grains (s = 4–6 μm), the optical thickness maps roughly match the phase space distributions of the top two rows. The orientation of the disk slightlydiffers between different particle sizes, as expected for the differences in orientation shown in the h–k plots. Furthermore, in the maps for both ≈12 and ≈226 μm, we observe the following: the disk is visibly offset from the star; the disk is slightly eccentric; the material density in the disk is the highest south of the star, and not at the closest distance to it. Moreover, the disks get larger for smaller grain sizes.
The offset and the eccentricity of the disk are explained by the planet shaping it. It is more extended for smaller particle sizes because, due to their β value, these particles are on more eccentric orbits with longer semimajor axes. The asymmetry in particle density is a geometric effect. An example of which can be seen in Fig. 1. On the south side, the particle orbits are arranged such that the paths of multiple orbits overlap nearly perfectly, resulting in a higher density. On the other sides, the paths of theorbits diverge significantly, resulting in a lower density. The same effect is discussed in detail by Wyatt (2005b).
The optical thickness of barely bound grains (s = 4–6 μm) shows a different behaviour. It is higher at the apocentre because in collisions of parent bodies at the pericentre, the β values are high enough for these small fragments to become unbound. For collisions at the apocentre, the fragments are still bound. This effect is also visible in Fig. 9 in Löhne et al. (2017) in logarithmic scale. Because this dependence on the collision point restricts the orbits (in addition to the secular perturbations), the axial asymmetry gets diminished.
The fourth row shows the linearly scaled thermal emission surface brightness of the debris disk at λ = 28.3, 70 and 1340 μm. The wavelengths were chosen to closely match astronomical telescopes, and for their associated characteristic grain size to match the particle sizes chosen for the other rows. In the presented case, 70 μm observations could have been made by Herschel/PACS, the 1340 μm observation could be made by ALMA, and the 28.3 μm observations could be made with JWST.
Grains of all sizes contribute to thermal emission and are taken into account in our modelling. Yet, at a given wavelength λ, the emission is dominated by a limited range of grain sizes around the characteristic grain size (Backman & Paresce 1993). In that sense, the wavelengths of the thermal emission images match the grain sizes shown in the phase space plots. We therefore expect the images to roughly reflect the features and differences of the above phase space distributions, albeit with lower significance because of the convolution over a size range. An example is the shortest wavelength (28.3 μm), where the peculiarities of the smallest, barely bound grains (s = 4–6 μm) are not reflected.
All thermal emission images show significant asymmetries. For all wavelengths we observe the following: the disk is visibly offset from the star; the disk is slightly eccentric; the parts of the disk closer to the star are brighter; the parts of increased brightness are not symmetric to the axis between the closest distance of the disk to the star and the star itself.
Again, the observed offset of the disk and its eccentricity is due to the planet shaping the initially circular disk into an ellipse. The disk is brighter at the pericentre^{5}, despite having less material than other regions, because the thermal emission intensity depends on the temperature of the grains. Particles spend more time at the apocentre of their orbit, where they cool down and do not contribute as much to thermal emission, as for the short period they spend close to the star, where they are heated up. The thermal emission therefore depends on the distance to the star. This effect is more pronounced for shorter wavelengths due to the exponential temperature dependence in the Wien regime, leading to a stronger radial dependence (Wyatt et al. 1999).
The axial asymmetry in the thermal emission derives from the optical thickness distribution, which shows more material south of the star. However, due to the pericentre glow the asymmetry is not as strong as for the optical thickness.
Fig. 7 Results of the ACE simulation described in Sect. 3.4. First row: h–k phase space distribution, second row: q–e distribution,and third row: corresponding geometrical optical thickness. Each column shows a different particle size range: ≈ 5 μm (β ≈ 0.43) on the left, ≈12 μm (β ≈ 0.18) in the middle, and ≈226 μm (β ≈ 0.009) on the right. Fourth row: thermal emission over all particle sizes (without applying any instrument function). The characteristic wavelengths shown, match the particle sizes from the plots above. 

Open with DEXTER 
4 Discussion
4.1 Relevance of differential precession
The crucial difference between the simulations shown in Löhne et al. (2017) and the models presented here is how the actual secular perturbations are implemented. Previously, we started with preperturbed belts that had already acquired eccentricity. The disks subsequently evolved only through the collisional cascade and drag forces, without ongoing dynamical perturbation. Now, secular perturbation and the resulting eccentricity evolution are modelled together.
BecauseLöhne et al. (2017) did not account for ongoing, differential precession, their results are applicable only if (a) the perturbations have stopped or (b) collisional replenishing dominates the evolution. Condition (a) could be fulfilled in a multi planet system, where a planet shapes the disk into an ellipse, but is then kicked out of the system in a close encounter with a more massive planet or substellar companion closer to the star. The belt would still be precessing, albeit at a lower rate and with a forced eccentricity closer to zero. If the eccentric belt itself is massive enough (cf. Krivov et al. 2018), it will however induce selfprecession and differential precession of its halo. Condition (b) is met if, for example, the mass of the perturbing planet is low and the resulting precession timescales much longer than collision timescales. Under such circumstances the perturber would still make the belt elliptic, but differential precession could not accumulate because particles would be too shortlived to significantly trail (or lead) freshly produced ones.
Vice versa, differential precession can be significant if precession timescales are sufficiently short because either the perturber or an already eccentric belt are massive enough.
4.2 Multiple perturbers
Our main reason for presenting the case of a single perturbing planet is simplicity. While many multiplanet systems have been observed, the influence of a single, massive planet is much easier to implement and describe because the planet’s orbital parameter do not evolve over time. The basic effects of differential precession can thus be shown more clearly. Our results do thus cover systems where one planet has a dominant influence on the disk, either because of proximity or because of mass.
It is possible to describe secular perturbation for multiple perturbers using the same analytic approximations in the h–k space (Murray & Dermott 2000) and an according implementation in ACE. As a result of the separate and the mutual perturbations of the planets, more complex patterns for precession and differential precession would arise. Given an observed snapshot of an eccentric belt, the already present degeneracy among potential perturber parameters would increase. However, deeper analysis is beyond the scope of this paper.
4.3 Finding hidden planets
If high resolution images at several wavelengths are available, the effects of differential precession on the size distribution can be spatially resolved. Ideally, the wavelengths and corresponding characteristic grain sizes would cover β(s) values from β ≈ 0.5 to β ≈ 0, that is, from barely bound grains in the halo to large grains that stay in the parent belt. Then, for example, the strength of the misalignment between spots of maximum brightness, as discussed in Sect. 3.4, can be used to derive timescales of differential precession relative to the collision timescales for grains in the corresponding size range. If collision timescales are estimatedfrom the relative abundance of halo grains with respect to the belt (Thébault & Wu 2008), the absolute timescales for differential precession can be estimated. Planetary mass and semimajor axis follow from Eq. (25).
In a collisionally active disk such as the one shown in Fig. 7, only the larger grains, and with them the parent belt itself, will have collision timescales long enough for significant differential precession to occur. However, for these large grains the spread in β values is necessarily narrow and size segregation will be unimportant. Precession will then be differential mainly with respect to semimajor axis, with the belt edge that is closer to the perturber precessing at a higher rate.
Given that the dust nearest to the star is brightest and that detection is most significant in bright regions of a disk, the location of the inner disk edge could in principle be measured as a function of observation wavelength. The inner edge is defined by the pericentre distances of the grains of corresponding characteristic size. Assuming that the belt has not evolved far from circularity, we approximate the evolution of the orbital eccentricities and the resulting pericentre distances from Eq. (13) as a function of time t: (27)
where the minus in the exponent is for an inner, the plus for an outer perturber. The precession timescale and the forced eccentricity of the belt are calculated using Eq. (13) with β = 0 and a = a_{belt}.
Equation (27) has effectively two free parameters, the planet’s semimajor axis and eccentricity. Thus, fitting it to observations cannot untangle the planet parameters completely. Competing mechanisms, such as chaotic motion in the zone of overlapping meanmotion resonances near a perturber (Mustill & Wyatt 2012), further complicate the analysis.
In principle, Eq. (27) could provide a simple constraint for possible planets in a debris disk system. However, as has been shown in this work, its practical application is limited by images being convolutions of grain size and orbit distributions. In contrast to thermal emission, scattered light images, which are more dominated by the isotropically scattering small grains, might retain more of the asymmetries. However, in that case, analysis is more sensitive to viewing geometry and grain properties. A more detailed discussion of scattered light is beyond the scope of this work.
4.4 Other approaches to differential precession
The modelling process presented above – dynamics implemented as advection terms in a gridbased collisional code – is computationally expensive. Analytic solutions would be desirable. However, such solutions are usually confined to the purely dynamical case, which is a good approximation only to the big objects with long collision lifetimes. For example, Mustill & Wyatt (2009) discuss the timescale for differential precession by a planet to induce orbit crossing in the belt, thus stirring the belt and triggering the collisional cascade. The actual subsequent collisional evolution cannot be modelled purely analytically.
Another possible modelling approach is to start from an Nbody code that can handle dynamics and add the generation and removal of particles in collisions. SMACK (Nesvold et al. 2013) and LIDTDD (Kral et al. 2013) are two examples that follow this approach. LIPAD (Levison et al. 2012) is a third example, albeit one less focused on dust. In such models, radiation pressure and gravitational interaction with as many perturbers as desired iseasily implemented. Shortterm and resonant perturbation or the aftermaths of single, giant collisions can be followed. While these shortterm phenomena are the domain of the dynamicsbased models, we consider the orbitaveraged statistical approach of ACE to be better suited for smooth longterm evolution and secular perturbations. A compromise between the two approaches could be based on variable instead of the fixed grid points used in ACE. Secular perturbation would then move the grid points along with the material instead of redistributing material in the grid. Compared to our current implementation in ACE this would eliminate the numerical dispersion discussed in Sect. 3.2. As a downside, the moving grid points would require collision properties to be evaluated again for each time step.
5 Summary
The evolution of debris disks under the combined action of a collisional cascade and secular gravitational perturbations were simulated fully consistently for the first time with our code ACE. We came to the following main conclusions:
 1.
The timescales on which secular perturbations create elliptic rings from initially circular distributions are sizedependent through the dust grains’ radiation pressuretogravity ratio β.
 2.
Differential precession, if observed, can help disentangle a perturber’s parameters, such as its mass and semimajor axis.
 3.
When collision lifetimes are long enough, differential precession can notably shear the smallgrain halo with respect to the larger objects in the parent belt. If collision lifetimes are short, shearing due to differential precession can only be observed within the belt itself.
 4.
The detection of planets with this method is pivoting on the availability of images at different wavelengths, covering a wide range of grains sizes, from barely bound halo grains to big debris unaffected by radiation pressure.
Acknowledgements
We thank the anonymous referee for constructive remarks that helped clarify critical sections of the manuscript. We are grateful for support by the Deutsche Forschungsgemeinschaft (project Lo 1715/21, Research Unit FOR 2285 “Debris Disks in Planetary Systems”).
References
 Backman, D. E., & Paresce, F. 1993, MainSequence Stars with Circumstellar Solid Material  The VEGA Phenomenon (Tucson, AZ: University of Arizona Press), 1253 [Google Scholar]
 Boley, A. C., Payne, M. J., Corder, S., et al. 2012, ApJ, 750, L21 [NASA ADS] [CrossRef] [Google Scholar]
 Bruggeman, D. A. G. 1935, Ann. Phys., 416, 636 [NASA ADS] [CrossRef] [Google Scholar]
 Burns, J. A., Lamy, P. L., & Soter, S. 1979, Icarus, 40, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Draine, B. T. 2003, ApJ, 598, 1017 [NASA ADS] [CrossRef] [Google Scholar]
 Eiroa, C., Marshall, J. P., Mora, A., et al. 2013, A&A, 555, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hauschildt, P. H., Allard, F., & Baron, E. 1999, ApJ, 512, 377 [NASA ADS] [CrossRef] [Google Scholar]
 Kalas, P., Graham, J. R., & Clampin, M. 2005, Nature, 435, 1067 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Kennedy, G. M., Marino, S., Matrà, L., et al. 2018, MNRAS, 475, 4924 [NASA ADS] [CrossRef] [Google Scholar]
 Kim, M., Wolf, S., Löhne, T., Kirchschlager, F., & Krivov, A. V. 2018, A&A, 618, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kral, Q., Thébault, P., & Charnoz, S. 2013, A&A, 558, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Krivov, A. V., Sremčević, M., & Spahn, F. 2005, Icarus, 174, 105 [NASA ADS] [CrossRef] [Google Scholar]
 Krivov, A. V., Löhne, T., & Sremčević, M. 2006, A&A, 455, 509 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Krivov, A. V., Ide, A., Löhne, T., Johansen, A., & Blum, J. 2018, MNRAS, 474, 2564 [NASA ADS] [CrossRef] [Google Scholar]
 Kuchner, M. J., & Holman, M. J. 2003, ApJ, 588, 1110 [NASA ADS] [CrossRef] [Google Scholar]
 Lee, E. J., & Chiang, E. 2016, ApJ, 827, 125 [NASA ADS] [CrossRef] [Google Scholar]
 Levison, H. F., Duncan, M. J., & Thommes, E. 2012, AJ, 144, 119 [NASA ADS] [CrossRef] [Google Scholar]
 Löhne, T., Krivov, A. V., & Rodmann, J. 2008, ApJ, 673, 1123 [NASA ADS] [CrossRef] [Google Scholar]
 Löhne, T., Krivov, A. V., Kirchschlager, F., Sende, J. A., & Wolf, S. 2017, A&A, 605, A7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Li, A., & Greenberg, J. M. 1998, A&A, 331, 291 [NASA ADS] [Google Scholar]
 MacGregor, M. A., Matrà, L., Kalas, P., et al. 2017, ApJ, 842, 8 [NASA ADS] [CrossRef] [Google Scholar]
 Maldonado, J., Eiroa, C., Villaver, E., Montesinos, B., & Mora, A. 2012, A&A, 541, A40 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Matrà, L., Marino, S., Kennedy, G. M., et al. 2018, ApJ, 859, 72 [Google Scholar]
 Moerchen, M. M., Churcher, L. J., Telesco, C. M., et al. 2011, A&A, 526, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Montesinos, B., Eiroa, C., Krivov, A. V., et al. 2016, A&A, 593, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Murray, C. D., & Dermott, S. F. 2000, Solar System Dynamics (Cambridge: Cambridge University Press) [CrossRef] [Google Scholar]
 Mustill, A. J., & Wyatt, M. C. 2009, MNRAS, 399, 1403 [NASA ADS] [CrossRef] [Google Scholar]
 Mustill, A. J., & Wyatt, M. C. 2012, MNRAS, 419, 3074 [NASA ADS] [CrossRef] [Google Scholar]
 Nesvold, E. R., Kuchner, M. J., Rein, H., & Pan, M. 2013, ApJ, 777, 144 [NASA ADS] [CrossRef] [Google Scholar]
 Neuman, S. P. 1981, J. Comput. Phys., 41, 270 [NASA ADS] [CrossRef] [Google Scholar]
 Ozernoy, L. M., Gorkavyi, N. N., Mather, J. C., & Taidakova, T. A. 2000, ApJ, 537, L147 [NASA ADS] [CrossRef] [Google Scholar]
 Pan, M., Nesvold, E. R., & Kuchner, M. J. 2016, ApJ, 832, 81 [NASA ADS] [CrossRef] [Google Scholar]
 Pawellek, N., Krivov, A. V., Marshall, J. P., et al. 2014, ApJ, 792, 65 [NASA ADS] [CrossRef] [Google Scholar]
 Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 2007, Numerical Recipes: The Art of Scientific Computing (Cambridge: Cambridge University Press) [Google Scholar]
 Reidemeister, M., Krivov, A. V., Stark, C. C., et al. 2011, A&A, 527, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rhee, J. H., Song, I., Zuckerman, B., & McElwain, M. 2007, ApJ, 660, 1556 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Rodigas, T. J., Malhotra, R., & Hinz, P. M. 2014, ApJ, 780, 65 [NASA ADS] [CrossRef] [Google Scholar]
 Schneider, G., Smith, B. A., Becklin, E. E., et al. 1999, ApJ, 513, L127 [NASA ADS] [CrossRef] [Google Scholar]
 Schneider, G., Grady, C. A., Stark, C. C., et al. 2016, AJ, 152, 64 [NASA ADS] [CrossRef] [Google Scholar]
 Sibthorpe, B., Kennedy, G. M., Wyatt, M. C., et al. 2018, MNRAS, 475, 3046 [NASA ADS] [CrossRef] [Google Scholar]
 Stapelfeldt, K. R., Holmes, E. K., Chen, C., et al. 2004, ApJS, 154, 458 [NASA ADS] [CrossRef] [Google Scholar]
 Strubbe, L. E., & Chiang, E. I. 2006, ApJ, 648, 652 [NASA ADS] [CrossRef] [Google Scholar]
 Thébault, P., & Wu, Y. 2008, A&A, 481, 713 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Thilliez, E., & Maddison, S. T. 2016, MNRAS, 457, 1690 [NASA ADS] [CrossRef] [Google Scholar]
 Wyatt, M. C. 2003, ApJ, 598, 1321 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Wyatt, M. C. 2005a, A&A, 433, 1007 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Wyatt, M. C. 2005b, A&A, 440, 937 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Wyatt, M. C., Dermott, S. F., Telesco, C. M., et al. 1999, ApJ, 527, 918 [NASA ADS] [CrossRef] [Google Scholar]
Contrary to Wyatt (2005b) and Pan et al. (2016) we do not observe apocentre glow because our images are not PSF broadened. See Fig. 15 of Löhne et al. (2017) for the effect of PSF broadening.
All Figures
Fig. 1 Threedifferent interior perturbers in red that can produce the same eccentric belt (thick black) from an initially circular belt. Orbits that deviate from the belt will be sheared by differential precession, as exemplified by orbits of older halo particles (lighter blue) trailing that of freshly produced ones (darker blue). The location of the star is marked by grey cross hairs. The perturber’s periapsis is indicated by red lines. 

Open with DEXTER  
In the text 
Fig. 2 Left: initial orbits of fragments with β values of 0.0 (orange), 0.1 (black), and 0.3 (blue), emitted in a single collision from a circular orbit at a = 100 au (and ϖ = 30°; particles moving anticlockwise). Particles with β ≥ 0.5 (red) are on unbound orbits. Right: initial positions (circle) of the same fragments in the h–k space and their circles of movement (dashed line) around the forced eccentricity (dot of according colour). The arcs they covered in the h–k space within 5 Myr are shown as solid lines. Given the planet at a_{p} = 20 au with e_{p} = 0.6 (and ϖ_{p} = 0°), they have precession times of (for β = 0.0), (for β = 0.1), and (for β = 0.3). 

Open with DEXTER  
In the text 
Fig. 3 h–k space evolution of parent belt particles (orange; a_{belt} = 100 au; β = 0.0; e ≈ 0.0), and particlescreated in collisions of parent belt objects; with fragment β values of β = 0.01;0.03;0.05 (black), and 0.3 (blue). The belt is around a star of M = 1.92 M_{⊙}, and the collisional fragments are perturbed by a m_{p} = 2.5 M_{Jup} mass planet, with a = 20 au and e = 0.6 (left); a = 60 au and e = 0.2 (right). The initial distribution of particles created in the collisions (dashed line), and the position after a time T (solid line) are shown. In between steps are marked with faint colours. 

Open with DEXTER  
In the text 
Fig. 4 Relevant timescales for fragments created by collisions inside a debris disk (a = 95– 105 au; m_{Dust} = 4.05 × 10^{25}g) around an A4Vstar (M_{⋆} = 1.92 M_{⊙}) with a planet at a_{p} = 20 au with an eccentricity of e_{p} = 0.6 and a mass of m_{p} = 2.5 M_{J}. The β value of the particles (blue, dotted), and the blowout limit (blue zone on the left) are shown. The precession timescale (orange; top), the crossing timescale (black; middle; where applicable), and the Poynting–Robertson timescale (red; crossing the top and the middle plot) are derived analytically. The collision timescale (black, dashed) was taken directly from the numerical results of the simulation. When possible, a range (shaded area) of timescales for a = 95– 105 au is given in addition to the value at a = 100 au (solid line). 

Open with DEXTER  
In the text 
Fig. 5 Sketch of the upwind scheme, used inside ACE, acting on the pericentre q and eccentricity e dimensions of the fourdimensional phase space. The material is shifted from each bin n_{i,j} into its direct neighbours by the parameter velocity v_{i,j} (calculated with perturbation theory). 

Open with DEXTER  
In the text 
Fig. 6 Azimuthal diffusion in the ACE simulation described in Sect. 3.2. The dust distribution of particles at a = 95 au (solid; crosses) and a = 110 au (dashed; circles) along the ϖ dimension is plotted for a simulation time of 0.00 Myr (red), 2.87 Myr (blue), 5.75 Myr (orange), and 8.62 Myr (black). A fit to the analytical solution was made to extract the different parameter velocities and dispersion coefficients. The initial distribution is scaled down by a factor of 5 to fit the plot. 

Open with DEXTER  
In the text 
Fig. 7 Results of the ACE simulation described in Sect. 3.4. First row: h–k phase space distribution, second row: q–e distribution,and third row: corresponding geometrical optical thickness. Each column shows a different particle size range: ≈ 5 μm (β ≈ 0.43) on the left, ≈12 μm (β ≈ 0.18) in the middle, and ≈226 μm (β ≈ 0.009) on the right. Fourth row: thermal emission over all particle sizes (without applying any instrument function). The characteristic wavelengths shown, match the particle sizes from the plots above. 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.