Capture and evolution of dust in planetary mean-motion resonances: a fast, semi-analytic method for generating resonantly trapped disk images

8.4 MB·English

by  Andrew Shannon

Mon.Not.R.Astron.Soc.000,000–000(0000) Printed9January2015 (MNLATEXstylefilev2.2)
Capture and evolution of dust in planetary mean-motion
resonances: a fast, semi-analytic method for generating
resonantly trapped disk images
5 Andrew Shannon1, Alexander J Mustill2,3, & Mark Wyatt1
1
1Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge, UK, CB3 0HA
0
2Departamento de F´ısica Te´orica, Universidad Aut´onoma de Madrid, Cantoblanco, 28049 Madrid, Spain
2
3Lund Observatory, Department of Astronomy and Theoretical Physics, Lund University, Box 43, SE-221 00 Lund, Sweden

n
a
J
9January2015

7

]
ABSTRACT
P
DustgrainsmigratingunderPoynting-Robertsondragmaybetrappedinmean-motion
E
resonanceswithplanets.Suchresonantlytrappedgrainsareobservedinthesolarsys-
.
h tem. In extrasolar systems, the exozodiacal light produced by dust grains is expected
p tobeamajorobstacletofuturemissionsattemptingtodirectlyimageterrestrialplan-
- ets. The patterns made by resonantly trapped dust, however, can be used to infer the
o
presenceofplanets,andthepropertiesofthoseplanets,ifthecaptureandevolutionof
r
t thegrainscanbemodelled.ThishasbeendonewithN-bodymethods,butsuchmeth-
s
a ods are computationally expensive, limiting their usefulness when considering large,
[ slowly evolving grains, and for extrasolar systems with unknown planets and parent

bodies, where the possible parameter space for investigation is large. In this work,
1
we present a semi-analytic method for calculating the capture and evolution of dust
v
grains in resonance, which can be orders of magnitude faster than N-body methods.
1
3 We calibrate the model against N-body simulations, finding excellent agreement for
6 Earth to Neptune mass planets, for a variety of grain sizes, initial eccentricities, and
1 initial semimajor axes. We then apply the model to observations of dust resonantly
0 trapped by the Earth. We find that resonantly trapped, asteroidally produced grains
. naturally produce the ‘trailing blob’ structure in the zodiacal cloud, while to match
1
0 the intensity of the blob, most of the cloud must be composed of cometary grains,
5 which owing to their high eccentricity are not captured, but produce a smooth disk.
1
1 INTRODUCTION with Neptune (Vitense, Krivov & Lo¨hne 2014). As similar
:
v processesmightbeexpectedtooccurinextrasolarplanetary
Xi SmalldustgrainsarereleasedintotheSolarsystemincolli- systems,thisraisestheinterestingpossibilitythatstructures
sions between asteroids and by the outgassing of comets. If in extrasolar debris disks are created by this phenomenon,
r
a they are not so small as to be unbound, they will spiral to- and can be used to infer the presence of planets and their
wardstheSunduetoPoynting-Robertson(PR)drag.These properties.
dust grains reflect, absorb, and re-radiate the solar light,
In the outer part (>> 10 au) of a stellar system, de-
producingthezodiacallightandzodiacalthermalemission.
bris disks can be spatially resolved by facilities and instru-
Aroundotherstars,similarprocessesmayproduceexozodi-
mentssuchasHST(Kalas,Graham&Clampin2005),Sub-
acal light. Such exozodiacal light is expected to be a signif-
aru (Thalmann et al. 2011), SMA (Hughes et al. 2011),
icant barrier to future missions attempting to directly im-
VLT (Buenzli et al. 2010), Herschel (Duchˆene et al. 2014;
age terrestrial planets if the level of dust is greater than
Matthews et al. 2014), Gemini (Wahhaj et al. 2014), and
∼ 10 times that of the Solar system (Backman et al. 1998;
ALMA(Boleyetal.2012).Thiscouldallowtheinferenceof
Beichmanetal.2006),alevelthatmaybepresentinmany,
resonantlytrappeddustfromthespatialdistribution(Kuch-
perhaps even most, systems (Kennedy & Wyatt 2013).
ner & Holman 2003; Wyatt 2003), although in these bright
Small grains undergoing PR drag can become trapped diskscollisionsarefarmoreimportantthanPRdrag(Wyatt
inmeanmotionresonanceswithplanets(Gold1975;Gonczi, 2005),whichmaywashoutsuchpatterns(Kuchner&Stark
Froeschle&Froeschle1982).Thisproducesclumpycircum- 2010). Possible detections of clumpy structures which may
stellar rings of grains, which are known to exist around the be associated with planetary resonances have already been
Sun due to trapping by the Earth (Dermott et al. 1994; reported (Holland et al. 1998; Greaves et al. 2005; Corder
Reachetal.1995;Reach2010),andVenus(Leinert&Moster etal.2009).Intheinnerparts(<<10au)ofsystems,bright
2007; Jones, Bewsher & Brown 2013). The dust detector dust disks are much rarer (Fujiwara et al. 2013), and tran-
aboard New Horizons may soon detect dust in resonance sient in nature (Wyatt et al. 2007). However, bright outer
©0000RAS
2 Shannon et al.
disks may leak dust inwards due to PR drag (Mennesson gravitational force. If the star is also orbited by a planet of
etal.inpress),whereresonantcapturemayoccur.PRdrag massm ,andorbitalparametersa ,e ,i ,(cid:36) ,Ω ,M ,res-
p p p p p p p
canalsoberelativelymoreimportantindimmerdisks.The onant interactions between the planet and dust grain may
KeckNullingInterferometer(Millan-Gabetetal.2011;Men- alsobeimportanttoitsorbitalevolution.Inthissectionwe
nesson et al. 2013) and the Large Binocular Telescope In- consider these effects analytically. These expectations will
terferometer(Hinz2013;Kennedyetal.2014)canallowfor then be compared to simulations in §3.3.
the detection of these dimmer disks. These interferometers
detectthefluxleakfromadustdiskthroughaninterference
2.1 Poynting-Robertson Drag
pattern. With repeated measurements of the same system,
rotating structures, such as resonantly trapped dust, may Assumingdustgrainstobespherical,theyexperienceaforce
be able to be inferred. from the stellar radiation of the form
Tomodelresonantlytrappeddust,variousgroupshave
used N-body simulations of dust grains. Dermott et al. (cid:12) (cid:12)(cid:20)(cid:18) r˙(cid:19) (cid:18)rθ˙(cid:19) (cid:21)
(1994) simulated the inspiral of 912 dust grains with diam- F(cid:126)radiation =βPR(cid:12)(cid:12)F(cid:126)gravity(cid:12)(cid:12) 1−2c rˆ− c θˆ , (1)
eters of 12 µm from the asteroid belt, and showed that the
pattern produced by the capture and persistence of those where βPR is the ratio of the force of radiation pressure to
grainsinresonanceiscompatiblewiththatobservedbythe gravity acting on the grain, which will depend on its size,
Infrared Astronomical Satellite. Deller & Maddison (2005) shapeandcomposition.Forexample,asphericalblackbody
integrated 500 dust grains under PR drag for 120 differ- grain of size s and density ρ has
ent systems, generating a catalogue of images for different 3L
β = ∗ , (2)
planetmass,planeteccentricity,dustsize,dustparentbody PR 8πcρGm s
∗
eccentricity, and dust parent body semimajor axis. Stark
where L is the stellar luminosity.
& Kuchner (2008) required a 420 core cluster to perform ∗
Averaged over an orbit, this causes dust particles to
simulations of 5000 particles spiralling past a planet under
evolveinsemimajoraxisandeccentricityas(Wyatt&Whip-
PRdragfor120differentparametercombinations.Similarly,
ple 1950)
Rodigas, Malhotra & Hinz (2014) used a suite of 160 sim-
uInlatthioenssetcoasperso,dthueceloangfitctoinmgpfuotramtiuolnatfoimredsisokfaNp-pbeoadryanscimes-. ddat(cid:12)(cid:12)(cid:12)(cid:12)PR = −Gma ∗βPcR((cid:0)12−+e32e)232(cid:1)
ulations makes their use to generate potential disk images
on a system-by-system case somewhat impractical, hence = −0.624auβ m∗ 1au(cid:0)2+3e2(cid:1), (3)
theaforementionedgroupshavegeneratedcataloguesofdisk kyr PRm(cid:12) a (1−e2)32
images and fitting formulae to be matched with disks once
discoTvehreesde.high computation times are set by the need to ddet(cid:12)(cid:12)(cid:12)(cid:12)PR = −25Gam2∗βPcR(1−ee2)12
resolve individual orbital periods of the planet(s) and dust
grains. However, the evolution under PR drag and the evo- = −1.56kyr−1β m∗ (cid:18)1au(cid:19)2 e . (4)
lution in resonance occur on much longer time-scales than PRm(cid:12) a (1−e2)12
theorbitalperiod.Thus,inthispaperwedevelopamodelof
Integratingequation3fromaninitialsemimajoraxisa and
resonant capture and evolution, with which dust grain tra- 0
eccentricitye =0,adustparticlewillreachthestarintime
jectoriescanbecalculatedonthePRdragtime-scale,rather 0
thantheorbitaltime-scale,resultinginordersofmagnitude a2c
τ = 0
decreases in the computation time. PR 4Gm β
∗ PR
In §2, we outline the physics of a single dust grain un- (cid:16) a (cid:17)2 m 1
= 0.40 0 (cid:12) kyrs. (5)
dergoing PR drag around a star where it may be captured 1au m β
∗ PR
intomeanmotionresonancewithaplanet.In§3weperform
Particleswithhighere reachthestarmorequickly,butthe
a suite of N-body simulations, and use them to validate 0
problem does not lend itself to a compact analytic form.
and calibrate our model. In §4, we describe how to use the
During the orbital evolution, there is a constant of in-
equationsofevolutiontoproducediskimages,andcompare
tegration which can be obtained from equations 3 and 4
them to the N-body simulations images. In §5, we compare
(Wyatt & Whipple 1950)
the output of the model to measurements of the structure
of Earth’s resonantly trapped ring. In §6, we discuss other K =a(cid:0)1−e2(cid:1)e−54. (6)
possible applications of the model.
ThevalueofK isconstantovermanyorbits,butsmall
variations occur during an orbit. Differentiating equation
6 with respect to time, substituting in a˙ and e˙ from Burns,
2 SINGLE PARTICLE PHYSICS Lamy&Soter(1979),andtakingtheexpressiontofirstorder
in e, we find
Considerastarofmassm ,orbitedbyadustgrain,theor-
∗
1 dK 8β v
bit of which is described by six parameters: semimajor axis ≈ PR kepler cosf, (7)
a, eccentricity e, inclination i, longitude of pericentre ((cid:36)), K df 5 e c
ascending node (Ω), and mean anomaly (M). The orbital and thus a libration in the constant K with a magnitude:
evolution of the dust grain will be dictated by the param-
∆K 8β v
eter β , the ratio of the force of radiation pressure to the ≈ PR kepler (8)
PR K 5 e c
©0000RAS,MNRAS000,000–000
Resonant capture of dust 3
over the course of an orbit. Asnoted,forafixeda thereisasmallrangeofpossible
p
aforthedustgrainsforwhichitisstillpossibletoconstruct
a resonant angle that can librate. Defining C as equation
r
2.2 Mean motion resonance
13, and otherwise following the derivation from Murray &
Adustgraincanbeinak-thordermeanmotionresonance Dermott (1999), the maximum libration width of a particle
with a planet if their mean motions are close to the ratio in a k=1 resonance is
j :j+k,wherejandkarepositiveintegers.Fromthemean (cid:18) (cid:19)1 (cid:18) (cid:19)1
longitudes λ=M +(cid:36) we construct a resonant angle ∆amax =± 16|Cr|e 2 1+ 1 |Cr| 2 − 2 |Cr|
a 3 n 27j2e3 n 9je n
j:j+1
ϕ=jλ −(j+k)λ+k(cid:36), (9) (19)
p
and for a k=2 resonance the width is
whereλandλ arethemeanlongitudeofthedustgrainand
p
(cid:18) (cid:19)1
the planet respectively. For our purposes, (cid:36) evolves slowly ∆amax =± 16|Cr|e2 2 . (20)
compared to λ and λp, so the j :j+k resonance occurs at aj:j+2 3 n
1 (cid:18)j+k(cid:19)32 We do not investigate resonances with higher k, as our N-
aj:j+k ≈(1−βPR)3 j ap, (10) body simulations produce very few captures at k>1.
which differs from the usual expression due to the effects
of radiation pressure. Here the match is only approximate 2.3 Resonance capture
as the k(cid:36) term allows for a small offset in a, and there
As particles migrate past a planet’s resonances due to PR
is a small libration of a about the nominal resonance lo-
drag, the question of whether they become caught in res-
cation corresponding to the libration in ϕ. While the dust
onance is deterministic. However, as it can depend on the
grain is near resonance, the perturbations from the planet
relative orbital phases of the planet and dust grain as res-
can be approximated by taking the terms in the disturbing
onance is approached, it can be treated probabilistically
function that are first order in eccentricity, and using La-
(Henrard1982;Quillen2006).Mustill&Wyatt(2011)usea
grange’s equations for the evolution of the orbital elements
simpleHamiltonianmodelwhichusesthelowestorderterm
to give (Murray & Dermott 1999)
of the disturbing function (as we did in §2.2) to calculate
(cid:12)
ddat(cid:12)(cid:12)(cid:12) =−2jCraeksinϕ, (11) tphaertcicalpestutrheartaetenscoaunndteirnifitirasltlaibnrdatsieocnonwdidotrhdserformmeaanssmleoss-
j:j+k
tionresonancesduringmigration,aswellastheeccentricity
ddet(cid:12)(cid:12)(cid:12)(cid:12) =kCrek−1sinϕ, (12) ktuicrkes. Wduhriinlegtrheesomnaondtelcrisosdsienvgeslotpheadtidnotnhoetcroenstuelxttinofctahpe-
j:j+k
capture of planetesimals by a migrating planet, the calcu-
where lations depend only on the rate of change of the bodies’
C = − Gmp (αf (α)+f (α)) (13) orbital separation, and so can be applied to the capture of
r na2a d i migrating dust particles by a fixed planet. However, that
p
√
model assumes β = 0, which requires a slight correction
− Gm (αf (α)+f (α)) PR
= p d i , (14) for our context. This is because when β > 0, the reduc-
√m∗(cid:16)j+jk(cid:17)12 (1−βPR)23 a2p tion in effective stellar mass seen by thePdRust grains due to
radiationpressuremeansthatresonancesoccurclosertothe
(cid:112)
where n = Gm∗(1−βPR)a−3, α = ap/a ≈ planet than when particles do not feel radiation forces, and
(1−β )−1/3((j+k)/j)−2/3, f (α) is the direct term of therefore the strength of a given resonance is different for
PR d
the disturbing function, and fi(α) is the indirect term of particles of different βPR.
the disturbing function. WepresentdetailsoftherequiredchangestotheHamil-
NeglectingPRdrag,andusingequation10forthedust tonian model, including changes to the resonant strengths
grain’ssemimajoraxis,takingthesecondtimederivativeof andeccentricitydampingeffects,intheAppendix.Basedon
ϕ, neglecting the contribution from λ¨ and (cid:36)¨ terms, one thesechanges,wenowpredictcaptureprobabilities,aswell
p
finds ϕ evolves as asinitiallibrationwidthsandeccentricitykicksduringnon-
capture resonance crossings, for real scenarios. We adopt a
ϕ¨ = 3j2C neksinϕ, (15)
r stellar mass of 1M , planet mass of 1m , and planet semi-
(cid:12) ⊕
= 3 j3 (cid:112)Gm (1−β )a−1C eksinϕ. (16) majoraxisof1au.Weintegratetrajectoriesofdustparticles
j+k ∗ PR p r with β = 0.005,0.01,0.02,0.04,0.08,0.16 and 0.32, with
PR
As C < 0, this is the equation for a pendulum with fre- initial eccentricities of 0.01. Encounters with resonances
r
quency fromthe2:1uptothe19:18areintegrated,unlesssaidreso-
nanceliesinsidetheplanet’sorbit.Infigure1,weshowthe
ω2 = 3j2|Cr|nek (17) captureprobabilitiesforfirst-orderresonancesasafunction
= 3 j3 (cid:112)Gm (1−β )a−1|C |ek. (18) ofdimensionlessparametersJ(cid:0)∝e2(cid:1)anddB/dt(∝da/dt).1
j+k ∗ PR p r In this figure we overplot the parameters corresponding to
capture in the first-order resonances of a 1M planet for
Thus, ϕ can librate or circulate, with libration in the case ⊕
that the dust grain is in resonance with the planet. This is
similartothecasewithoutradiationpressure,but|Cr|con- 1 We use B for the parameter denoted β in MW11, to avoid
tainsdependenceonβPR,andα,whichalsodependsonβPR. confusion.
©0000RAS,MNRAS000,000–000
4 Shannon et al.
azimuthalvelocity(v ),thechangeinspecificenergyofthe
θ
5:4 dust grain is
δE ≈ 2 Gm∗ (cid:18)mp(cid:19)2(cid:16)1−(cid:112)1−β(cid:16)1− (cid:15)(cid:17)(cid:17)−3. (23)
(cid:15)2 a m 2
p ∗
Foraninitial(cid:15) andafinal(cid:15) ,thespecificorbitalenergy
1 2
change is also
10:9
Gm (1−β)
δE ≈ ∗ ((cid:15) −(cid:15) ), (24)
2a 1 2
p
11:10 and equating the two energies yields
19:18 19:18 (cid:15)1−(cid:15)2 = 1−4β (cid:18)mmp(cid:19)2 (cid:15)12 (cid:16)1−(cid:112)1−βPR(cid:16)1− 2(cid:15)(cid:17)(cid:17)−3,
PR ∗
β =0.005 (25)
βPR=0.01 19:18
ββPPRR==00..0042 which assumes (cid:15)1−(cid:15)2 <<(cid:15).
βββPPPRRR===000...310268 19:18 aparAt finterloonngeitucodnej,umncetainonin,gththeaptlatnheetnaenxdt dcounstjugnrcatiinonmwovilel
PR
happen at
(cid:18) (cid:112) (cid:18) 3 (cid:19)(cid:19)−1
λ =2π 1− 1−β 1− (cid:15) . (26)
c PR 2
Figure 1. Capture probabilities calculated by integrating the If we assume that interactions between the dust grain
Hamiltonianmodel(white=100%,black=0%),includingeccen- and the planet will be dephased if the change in the longi-
tricity damping at k =5×10−5. The x-axis measures a dimen- tudeofconjunctionδλcisgreaterthansomeamount∆,this
sionless “eccentricity” (cid:0)J ∝e2(cid:1), and the y-axis a dimensionless yields a constraint that
migrationrate(dB/dt∝da/dt).Overplottedarethedimension-
lessmomentumandmigrationratecorrespondingtothecompar- √ 12π (cid:18)mp(cid:19)2 1 (cid:16)1−(cid:112)1−β (cid:16)1− (cid:15)(cid:17)(cid:17)−3
isonintegrationsdescribedinthetext. 1−βPR m∗ (cid:15)2 PR 2
×(cid:18)1−(cid:112)1−β (cid:18)1− 3(cid:15)(cid:19)(cid:19)−2 (cid:62)∆, (27)
PR 2
dust particles with β = 0.005 up to β = 0.32. We can
PR PR
see that many of these resonances lie in the regime of pos-
sible but not certain capture, meaning that a population of which reduces to the scaling of Wisdom (1980) in the limit
drifting dust grains will populate a range of resonances. βPR →0 of
(cid:18) (cid:19)1 (cid:18) (cid:19)2
2.4 Maximum j of capture (cid:15)(cid:54) 128π 7 mp 7 . (28)
3∆ m
∗
Wisdom (1980) showed that first order resonances over-
A reasonable expectation is that ∆ = 2π, which results in
lap when they are at semimajor axes less than ∼
a leading coefficient of ≈1.5, slightly larger than the value
1.3(m /m )2/7a greater than that of the planet. Bodies
p ∗ p obtained by Wisdom (1980), though numerically, Duncan,
movechaoticallybetweentheoverlappingresonances,which
Quinn&Tremaine(1989)foundthiscoefficienttobe≈1.5,
makes them unstable. Consequently, we would expect such
and Chiang et al. (2009) preferred a coefficient of ≈ 2.0 to
resonances to not capture migrating dust grains. Here we
fittheirsimulationsofFomalhaut’sdisk.Recentresultshave
employ a different approach to that of Wisdom (1980) to
also shown some dependence on system age and collision
calculate the width of this chaotic zone, but which recovers
time(Nesvold&Kuchner2014;Morrison&Malhotra2014).
their result, at the same time accounting for the effects of
Forourownsimulations,wewillfitthiscoefficientin§3.3.2.
radiation pressure.
Intheinertialframeofaplanetwithsemimajoraxisa ,
p
a dust particle with radiation pressure co-efficient βPR on 2.5 Resonant evolution
a circular orbit with semimajor axis a (1+(cid:15)) approaches
p
In the absence of drag a dust grain resonating with an in-
from infinity with velocity
terior planetary perturber librates around ϕ=π (Equation
(cid:115)
v ≈ Gm∗ (cid:16)1−(cid:112)1−β (cid:16)1− (cid:15)(cid:17)(cid:17), (21) 15). With the addition of PR drag, the centre of libration
∞ a PR 2 changes. Momentum and energy loss to PR drag must be
p
offsetbyagainofmomentumandenergyfromtheresonant
assuming (cid:15) << 1, an assumption we make throughout this interaction. Equating equations 3 and 11, we can estimate
derivation. The particle’s path is deflected by an angle θ, that ϕ librates around ϕ , given by
0
given by (cid:32)n v β (cid:0)2+3e2(cid:1)(cid:33)
θ ≈ 2mmp1(cid:15) (cid:16)1−(cid:112)1−β(cid:16)1− 2(cid:15)(cid:17)(cid:17)−2 (22) ϕ0 = sin−1 jC:jr+k keplecr,j:j+k2jPeRk (1−e2)32 , (29)
∗
assuming the angle is small. Converting back to the frame = sin−1(cid:32)(1−βPR)31 βPR 2+3e2 Gm∗ j13 (cid:33).
ofthestar,andconsideringthefinalradialvelocity(vr)and 2ekCr c (1−e2)32 a2p (j+k)43
©0000RAS,MNRAS000,000–000
Resonant capture of dust 5
Note that this expression does not hold for the 2 : 1 motion resonance with the planet. For a particle not cap-
resonance, or other resonances which exhibit asymmetric turedinlowerj resonances,themethodof§2.3calculatesa
libration for which higher order terms in the disturbing probabilityofcaptureintothe6:5resonanceof∼0.3.Allow-
function are important (Message 1958). For those reso- ing for capture into lower j resonances, the overall capture
nancestherearetwopossiblecentresoflibrationwithdiffer- probabilityintothe6:5meanmotionresonanceforparticles
ing capture probabilities (Murray-Clay & Chiang 2005). In withsimilarinitialconditionsis∼0.24.Thus,weexpectthis
this case, a numerically calibrated prescription from Wy- is a typical outcome. After capture, the semimajor axis of
att (2003) may be employed, which we do in §4. They the dust particle oscillates, but remains within the possible
found that the centre of libration is given by cosϕ2:1 = values for a resonant dust grain as outlined in equation 19
0
0.39−0.061/e. Using the PR-drag migration rate of equa- (top left panel of figure 2). After capture, the evolution of
tion 3, their results would predict that the relative chance the eccentricity is well described by equation 31 (top right
of capture into the lower 2 : 1 resonance is 1/2 − panel of figure 2). The expected centre of libration of the
0.85(cid:112)Gm∗βPR/c(mpa)−0.25(cid:0)2+3e2(cid:1)0.5(cid:0)1−e2(cid:1)−0.75. resonance angle (ϕ0) changes quickly at early times, as ex-
While the dust particle is in resonance, it undergoes pected from equation 29, while the centre of libration of
only small variations in a. The eccentricity, however, may theresonanceangleinthesimulationcatchesuponamuch
evolvesubstantially.Takingthevalueofϕ fromequation29 slowertimescale(bottomleftpaneloffigure2).Thelibration
0
andsubstitutingitintoequation12givesanaveragee˙ from widthgrowsexponentiallywithtime(bottomrightpanelof
resonant interaction. Combined with the e˙ from PR drag figure2).At3.6×105 years,theparticleescapesfromreso-
(equation 4), the eccentricity evolves according to nance. Later on we will quantify when and how dust grains
escape from resonance.
de = 1 Gm∗ βPR (cid:32) k(cid:0)2+3e2(cid:1) − √ 5e (cid:33). After escaping from resonance, the particle continues
dt 2a2j:j+k c (j+k)e(1−e2)23 1−e2 to evolve under PR drag decreasing its semimajor axis as
(30) equation3anditseccentricityasequation4.Atabout3.9×
Expandingtosecondorderineandsolvingyields(Weiden- 105 years, the particle is removed from the simulation once
schilling & Jackson 1993) a∼0.05,atwhichpointtheintegrationtimestepistoolarge
to resolve the orbit.
(cid:115)
e(t)= 5(j2+k k)(cid:16)1−e−τte(cid:17). (31) derstFarnodmintghiosfetxhaemcpalpe,tuwreecparnocseeses.thHaotwweevehra,vaesswomeehuavne-
noted, some elements are not well understood (libration
Note that the expression in Weidenschilling & Jackson widthgrowth,resonanceescape),andothers(captureproba-
(1993) has a typographical error: it is missing the square bility,resonancewidthsatcapture)wouldbenefitfromaddi-
root. Here the growth time of the particle’s eccentricity is tionalvalidation.Tothisend,weperformasuiteofN-body
simulations varying system parameters of interest.
τ =0.2a2 c/(Gm β ). (32)
e j:j+k ∗ PR
3.2 Suite of simulations: Setup
3 N-BODY INTEGRATIONS
The central star is characterised by one property, its mass
3.1 Single particle example: Evolution under PR
(m ). In all simulations, m =m .
∗ ∗ (cid:12)
drag
Weincludeoneplanetofmassm ,withsemimajoraxis
p
Toenhancequalitativeunderstanding,webeginbypresent- ap,witheccentricityep,inclinationip,longitudeofpericen-
ing the evolution of a single dust grain (Figure 2). We per- tre((cid:36)p),ascendingnode(Ωp),andmeananomaly(Mp).In
form this simulation, and all other simulations, using the all cases, ep =0, ip =0. In all cases (cid:36)p, Ωp, and the initial
hybridintegratorofMERCURYsuiteofN-bodyintegrators Mp are chosen randomly and uniformly between 0 and 2π.
(Chambers1999),whichwehavemodifiedtoincludetheef- Unless otherwise stated mp =m⊕ and ap =1au.
fectofradiationofdustgrains,asgiveninequation1.Totest The dust grains are assumed to have zero mass, and
thatitwasimplementedproperly,weconfirmedthatourN- the orbits are characterised by the same properties as the
bodycodeproducesaK thatisconstantonlongtimescales planet: semimajor axis a, with eccentricity e, inclination i,
with short timescale variations of size ∆K (equation 8). longitude of pericentre ((cid:36)), ascending node (Ω), and mean
Wepresenthereaparticlefromasimulationwithasolar anomaly (M). As with the planet, (cid:36), Ω, and the initial
mass star, an earth mass planet on a circular orbit at 1 au, M are chosen randomly and uniformly between 0 and 2π.
and a dust grain with β = 0.01, e = 0.01, a = 2.22au, Unless otherwise stated, the dust grains begin with eccen-
PR 0 0
i = 0.0628. This particle is from the simulation labelled tricity2 of 0.01, inclination3 of 2π ×10−2, and semimajor
0
B in table 1, where we describe the total range of system axis chosen randomly and uniformly between 2.20 ap and
parameters we consider in section 3.2. The particle evolves
initially under only PR drag, decreasing its semimajor axis
as per equation 3 and its eccentricity as per equation 4. At
2 Thisistheinitialeccentricityofthedustgrain,andneglectsthe
∼ 105 years, it encounters the 2:1 mean motion resonance,
which does not capture the particle, but provides an eccen- cthonansidtehreairtiopnartehnattbhoidgihesβ(P§R4.1a)r.eWcreeaintecdoraptorhaitgehetrhaetccoennltyriwcihtieens
tricity kick. Similar jumps are discernible when it crosses
consideringafulldiskmodel(section4)
the 3:2, 4:3, and 5:4 mean motion resonances. 3 Comparisons to runs at i = 0.0314, and i = 0.0157 show no
At 1.5×105 yrs, theparticle iscaptured in a6:5 mean significantdifferences.
©0000RAS,MNRAS000,000–000
6 Shannon et al.
33 11
111...111555 Dust grain
)) Weidenschilling & Jackson 1993
UU 22..55
AA
axis (axis ( 22 111...111 111 000555 222xxx111000555 333xxx111000555 ricityricity 00..11
r r 11..55 ntnt
oo ee
ajaj cc
mm 11 EcEc 00..0011
ii
mm
SeSe 00..55 Dust grain
a
6:5
00 00..000011
110055 22xx110055 33xx110055 110055 22xx110055 33xx110055
TTiimmee ((yyrrss)) TTiimmee ((yyrrss))

222666000 90
222444000 ) 80
222222000 D j( 70
eee 222000000 dth 40
glglgl 111888000 wi 50
nnn
AAA n
111666000 o
i 40
t
111444000 a
r
RReessoonnaannccee aannggllee ((jj )) b
111222000 EExxppeecctteedd CCeennttrree ooff LLiibbrraattiioonn ((jj )) Li 30
Simulated Centre of Libration (j 00)
111000000 0

000 111000555 222xxx111000555 333xxx111000555 0 105 2x105 3x105
TTTiiimmmeee (((yyyrrrsss))) Time (yrs)
Figure 2.Evolutionofadustgraininsemimajoraxis(topleft),eccentricity(topright),resonanceangle(bottomleft),andresonance
width(bottomright).ThedustisaffectedbyPRdragandradiationpressure,andperturbedbyaplanetwithap=1auandmp=m⊕.
At105yrs,thedustgraincrossesthe2:1resonanceandisperturbedtoahighere.Similarjumpsarediscerniblewhenitcrossesthe3:2,
4:3, and 5:4 mean motion resonances. At 1.5×105 yrs, the particle is captured in a 6:5 mean motion resonance with the planet. The
eccentricitygrowsquicklyinwhileinresonance,flatteningoutataconstantvalue,closelyfollowingequation31.Thedustgrainescapes
theresonanceat3.6×105yrs.AfterescapePRdragbringstheparticleintothestarat3.8×105yrs.
2.25 a . We choose this initial semimajor axis as we are in- cases spacing trial points by factors of two in the relevant
p
terested in first and second order mean motion resonances, parameter (table 1). Thus we perform a total of 25 simu-
themostdistantofwhichisthe3:1,whichforadustparticle lations, each with 104 particles. Simulations are run with
with β = 0 is located at 2.08 a , and closer for β > 0 an8(0.01/β )(a/1au)1.5 daytimestep,withdataoutputs
PR p PR PR
(equation 10). Thus, a = 2.20−2.25 a is effectively an every104(0.01/β )(a/1au)1.5 days.Particlesareremoved
p PR
arbitrarily large distance for all choices of β .4 fromthesimulationwhentheyaremorethan100au,orless
PR
We aim to be able to simulate disks over a wide pa- than 0.05 au, from the star.
rameter space in βPR, mp, ap, and e0. A complete sur- As all results are for m∗ =M(cid:12), any fits to results may
vey of the parameter space is computationally prohibitive. have unknown stellar mass dependence. Because our simu-
We choose a canonical simulation βPR = 0.01, mp = m⊕, lations generate a negligible amount of k (cid:62) 2 captures, we
ap = 1au, e0 = 0.01 (simulation B), and explore these di- only analyse the j+1:j captures in our simulations.
mensionsbyholdingthreeoftheparametersconstantwhile Simulations continue until all particles have been re-
varying the fourth one5. We consider βPR = 0.005...0.32, moved.
m =1m ...256m ,a =1au...16au,e =0.01...0.64,inall
p ⊕ ⊕ p 0
4 Because PR drag decreases e and a for particles with larger
initiala,thisisn’tstrictlytrue.Wecanuseequations3and4to
evolve dust grains at higher a to an e at 2.2au, so this choice is
effectiveforourpurposes.
5 Non-zeroinitialeandiaremorephysicalthanexactlyzero,and
avoid possible degeneracies/singularities that might occur when components (λp,λ,(cid:36)) are not well defined for zero eccentricity
they’repreciselyzero;inparticular,theresonantangleϕandits andinclination
©0000RAS,MNRAS000,000–000
Resonant capture of dust 7
Simulation e0 βPR mp ap y β
ABCD 0000....00001111 0.0000...00005142 mmmm⊕⊕⊕⊕ 1111aaaauuuu probabilit 11 ..241 βββββ======0.000000.....000001512486
E 0.01 0.08 m⊕ 1au re
F 0.01 0.16 m⊕ 1au ptu 0.8
G 0.01 0.32 m⊕ 1au ca 0.6
H 0.02 0.01 m⊕ 1au ve
I 0.04 0.01 m⊕ 1au ati 0.4
J 0.08 0.01 m⊕ 1au mul 0.2
K 0.16 0.01 m⊕ 1au u
L 0.32 0.01 m⊕ 1au C 0
M 0.64 0.01 m⊕ 1au 5:4 10:9 15:14
Resonance
N 0.01 0.01 2m⊕ 1au
O 0.01 0.01 4m⊕ 1au
P 0.01 0.01 8m⊕ 1au Figure 3. Cumulative capture probability for particles with
Q 0.01 0.01 16m⊕ 1au e0 = 0.01, i0 = 2π×10−2, ap = 1au, mp = m⊕, for various
R 0.01 0.01 32m⊕ 1au βPR (simulations A-G). Comparing the outcomes of our simula-
S 0.01 0.01 64m⊕ 1au tions (red solid lines) with the calculations in §2.3 (blue dotted
T 0.01 0.01 128m⊕ 1au lines).Thefitworkswelloverall,althoughthecriticalj cutoffis
U 0.01 0.01 256m⊕ 1au perhapssmootherintheN-bodysimulationsthaninourHamil-
V 0.01 0.01 1m⊕ 2au tonianapproach.
W 0.01 0.01 1m⊕ 4au
X 0.01 0.01 1m⊕ 8au
TtioYanblBe f1o.rmAstaabc0le.e0n1torfe,alwl0i.tt0hh1eostihmeru1lsamitmi⊕ounlastwi1o6ensapuvearfroyrimngedb.ySfiamctuolars- eeeee00000=====00000.....1000061248’ 11
m33of..a332s.s1a,looSrnTugprilataaepxnpeeoistnfosgeSfmpiinmrimoitubiaaaljlbaoiertlciiatcoyxennisst.r:icRitey,supaltrtsicle size (βPR), planet Cumulative capture probability 000 111
0
Toautomatedlyidentifyresonancecapture,welookfordust
grains that are continuously between a +2∆a and
j:j+k max
0
aj:j+k−2∆amax foratleasttwicethetimepredictedbydi- 5:4 10:9 Re son an c e15:14 20:19
viding the crossing distance 4∆a by the instantaneous
max
migration speed when the grain first encounters the reso- Figure 4. Cumulative capture probability for particles with
nance,givenbyequation3.Weemploythefactorof2inthe βPR = 0.01, i0 = 2π×10−2, ap = 1au, mp = m⊕, for var-
expected width to allow for deviations from our expression ious e0 (simulations B, H-M), comparing the outcomes of our
for the width owing to higher order terms. We also require simulations(redsolidlines)withthecalculationsin§2.3fordif-
thedustgraintobeatthatsemimajoraxisforatleasttwice ferentinitialeccentricities(dottedbluelines).Eachinitialeccen-
as long to avoid falsely identifying dust grains as captured tricity is offset vertically by 0.75 for clarity. The e0 = 0.32 and
either when they are slowed by the resonant perturbation e0 =0.64 case both have <1% capture in both N-body simula-
tionsandthepredictionsof§2.3,andarenotplottedhere.Atlow
butnotcaptured,orwhenoscillations(asequation8)cause
eccentricitytheagreementisverygood,anditremainsreasonably
variations in a that result in dust grains staying slightly
goodatalleccentricities.
longer at the resonance than equation 3 predicts, or when
closeencounterswiththeplanetcauseatoevolveotherthan
aspredictedbyequation3,usingtheeccentricityofthedust of β (figure 3 - simulations A-G) and for different val-
PR
grain when it first encounters the resonance. A spot check ues of e (figure 4 - simulations B, H-M). This comparison
0
of 80 particles found that this approach correctly identified considers only the first capture of particles in the N-body
57 of 58 times when a particle crossing a j +1 : j reso- simulation;anysubsequentcapturesareneglectedaseccen-
nance was captured, and 359 of 363 times when a particle tricity growth while in resonance (§2.5) makes such com-
crossingaj+1:j resonancewasnotcaptured,considering parisons unsuitable for evaluating particles with the initial
onlycrossingsuptothefirstreportedcapture.Thetrapping conditions of the experiment. This comparison accounts for
lifetime (or resonance escape time) is recorded as the total the eccentricity evolution as the semimajor axis decreases,
time the dust grain spends between a +2∆a and andtheeccentricitykicksduetoresonancecrossings(figure
j:j+k max
a −2∆a lessthetimepredictedbyintegratingequa- 2);thequotedeccentricitiesaretheinitialvalues.Thecalcu-
j:j+k max
tion 3 across that change in semimajor axis. lated capture probabilities match the model well, although
Wecompareourobservedcaptureprobabilitiestothose slight differences exist.
calculated in §2.3, and plot a comparison of the cumulative The 2:1 resonance presents a particular challenge as
capture probability as a particle crosses for different values asymmetric libration results in capture into both the 2:1(l)
©0000RAS,MNRAS000,000–000
8 Shannon et al.
1
14 0 1
n
o
12 cti 0 1
ra 0 1
10 f
al j ve 0 1
Critic 8 ulati 0 1
6 m
0 1
4 N−Nb−odbyo deys tliimmaittes Cu 0 4567::::3456
Hamiltonian model 8:7
2 ap = ad 190::89
0 20 40 60 80 100
Initial resonance width (δφ )
0.01 0.1 0
b
Figure5.Thecriticaljabovewhichcapturecannotoccur.Lower Figure6.Cumulativedistributionoftheresonantanglelibration
bounds come from simulations A, B, and C, where the N-body width at capture (δϕ0). Each resonance is offset by 1 vertically
simulationscapture∼100%ofparticles,socapturemustbeable from one another for clarity. Width is measured using the max-
to occur until the j where the cumulative capture probability is imum and minimum libration angles (ϕ) while 0.04<e<0.06.
∼100%. The upper bound is from simulation F, where the lack Thick lines are simulations, thin lines are predictions from sec-
of any captures constrains the maximum j to be less than the tion 2.3 (for the most part, the two lines overlap too much to
lowestjforwhichcaptureispredictedintheHamiltonianmodel. be distinguished). These are all taken from the default case of
The N-body estimates are from simulations D and E, with the βPR=0.01,e0=0.01,mp=m⊕,ap=1au(SimulationB).The
critical j taken from where the Hamiltonian model and N-body overallagreementisexcellent,exceptinthe4:3case,wherethe
simulationspredictionsforcaptureprobabilitydiverge.Thegreen simulationhadonly55captures.
points show the j given by equation 27, with ∆≡2.3 chosen to
bringtheHamiltonianmodelintoagreementwiththeN-bodyre-
sults. The solid pink shows where the dust grain would orbit at simulationresultsfromsimulationB,thecanonicalcase,us-
thesamesemi-majoraxisastheplanet;heretheimpactparame- ing the maximum and minimum libration angles (ϕ) while
tergoestozero,andhencethecalculatedkickbecomesinfinitely 0.04<e<0.06. This is chosen to be confident the particle
large. isinresonance,buthasnotexperiencedsignificantlibration
width growth. The agreement is generally excellent, except
inthelowestjcase,wherethepredictedlibrationwidthsare
and 2:1(u) resonance, with differing probabilities (Message
slightly too large. We will see that it is important this dis-
1958;Chiang&Jordan2002;Murray-Clay&Chiang2005).
tribution is well characterised, as the initial libration width
WeemulatetheapproachofWyatt(2003)andcalibratethe
sets the length of time the particle remains in resonance
relative capture probabilities against our simulations in the
after capture.
form
P (l)=0.5−aθbµc, (33)
2:1
(cid:112) 3.3.4 Eccentricity evolution
where θ = (a˙/1au/Myr)/ (a/1au)(M /m ). We per-
(cid:12) ∗
formed limited extra simulations of capture where capture Insection2.5,wepresentedanexpressionfortheeccentricity
into both the (u) and (l) resonance occurred, with 1000 evolution of a particle caught in resonance (equation 31).
particles; sampling planet mass every 25m from 150 − Here,weevaluatehowcloselythatexpressionmatcheswhat
⊕
400m (with β = 0.01) and migration rate with factors of we find in N-body simulations.
⊕
2 in β from β =0.0025 to β =0.02 with m =256m . Our In figure 7, we compare the evolution of eccentricity
p ⊕
best fit parameters were a = 0.01, b = 0.25, and c = −0.4, while in resonance for particles from simulation B to equa-
whichweemployin§4topredicttherelativeratesofcapture tion31,plottingallparticlescaughtin5:4,6:5,7:6,8:7,9:8,
into the 2:1(l) and 2:1(u) resonances. and10:9meanmotionresonancesforatleast5×104 years.
The evolution timescale is in excellent agreement, and the
magnitudeisinreasonableagreement,withonlya∼4%de-
3.3.2 Maximum j cutoff
viation of the simulation from the analytic solution. The
Our derivation of the highest j at which capture can occur analytic solution assumes particles are captured at zero ec-
hasanumericalcoefficient∆thatwasnotspecifiedin§2.4. centricity, but as the growth is an exponential decay away
Comparison of our Hamiltonian capture model with our N- from 0, this difference is not significant.
body simulations favours a coefficient of ∆ ≈ 2.3 to bring
the two into agreement (figure 5). This is slightly higher
3.3.5 Centre of Libration
than values previously reported (see section 2.4).
Asseeninfigure2,aftercapture,theeccentricityrises,and
the predicted centre of libration moves quickly towards a
3.3.3 Libration width at capture
value set by the final eccentricity. The centre of libration
The model discussed in section 2.3 makes a prediction for found in simulation lags the predicted value, however, as
what the distribution of libration widths should be at cap- seen in the bottom left panel of figure 2. The same general
ture (δϕ ). In figure 6 we compare that prediction to our behaviourofaninitialoffsetwithalinearrelaxationtowards
0
©0000RAS,MNRAS000,000–000
Resonant capture of dust 9
40 4:3
5:4
0.3 6:5
35 78::67
9:8
0.25
30
y
entricit 0 .01.52 um offset 2205
c m
Ec 0.1 axi 15
M
0.05 10
5
0
0 50000 100000 150000 200000 250000 300000
Time in resonance (yrs) 0 0 20 40 60 80 100 120
Initial libration width
Figure 7. A comparison between the expected evolution of the Figure 8. Maximum offset between the predicted centre of li-
eccentricity e predicted from equation 31 (solid lines), bration and that found in simulations, compared to the initial
predicted
andthatfoundinthesimulationsesim (dots).Allparticlesplot- libration width, for particles from simulation B. The solid lines
ted here are from the canonical case of βPR = 0.01, e0 = 0.01, areafittoallsimulations,givenbyequations35,36,and37.
mp = m⊕, and ap = 1au (simulation B). Plotted are 1000 ran-
dom particles caught in 5:4 (black), 6:5 (red), 7:6 (green), 8:7
(blue), and 9:8 (magenta) mean motion resonances for at least 40
4:3
5×104 years. 56::45
35 7:6
8:7
9:8
30
theequilibriumvalueisalsofoundfortheotherparticlesin et
the simulations. Thus, we set out to characterise this be- offs 25
haviour with a phenomenological model of the form um 20
m
ϕsim =ϕcal+∆ϕ (cid:18)1− t (cid:19). (34) Maxi 15
0 0 o τ
relax
10
Herethemaximumoffsetinthecentreoflibration,∆ϕ and
0 5
the relaxation time, τ, are parameters that we will fit phe-
nomenologically, and t is the time after resonance capture. 0
0 100000 200000 300000 400000 500000 600000
We begin with a fit of the maximum offset (∆ϕ ) be- Relaxation time (years)
0
tween the equilibrium and actual centre of libration, from
Figure9.Thetimeittakesforthelibrationcentretorelaxtothe
simulationB,plottedagainsttheinitiallibrationwidth,for
equilibriumvalueofequation29,comparedtotheinitiallibration
several resonances (figure 8). Here, the maximum offset is
width,forparticlesfromsimulationB.Thesolidlinesareajoint
found to depend only on two parameters, the j of the reso- fittoallthesimulationsoftheform∆ϕ0=C3(cid:0)1−e−τrelax/C4(cid:1).
nance,andtheinitiallibrationwidthofthecapturedparti-
cle, δϕ (§3.3.3). The best fit is of the form
0
with
∆ϕ =C cos(δϕ /C ). (35)
0 1 0 2
C =7605j−1.03β0.9m−0.94a−0.45 (39)
Fitting the same for all simulations, we find 3 PR p p
and
(cid:18)m (cid:19)−0.864(cid:16) a (cid:17)−0.423
C1 =4475βP0.R847j−0.81 m⊕p 1apu (36) C4 =13949j−1.54βP−R0.79m−0.385a1p.73. (40)
and An example of this fit for simulation B is plotted in figure
(cid:18)m (cid:19) (cid:16) a (cid:17) 9. Here, of course, the combination of equations 35 and 38
C =163.9−1.76j−1.4 p +0.73 . (37) allows one to express the relaxation time as a function of
2 m 1au
⊕ the initial libration width in the form
We neglect the β dependence of C as the best fit coeffi- (cid:18) (cid:18) (cid:19)(cid:19)
PR 2 C δϕ
cient of −1.93±1.76 is consistent with zero. τrelax =C4ln 1− C1 cos C0 , (41)
3 2
Plotting the maximum offset against the relaxation
time,wefindthetwohaveanapproximatelyfunctionalcor- and thus for particles with the same system parameters
respondence (figure 9). We fit the maximum offset between caught in the same resonance, their evolution is dictated
the predicted centre of libration and that found in simula- entirely by their initial libration width.
tions against time it takes for the centre of libration to to
relaxtheequilibriumvalueofequation29.Thebestfitvalue
is 3.3.6 Libration width evolution
(cid:18) ∆ϕ (cid:19) The resonance width (δϕ) was found to grow exponentially
τ =C ln 1− 0 (38)
relax 4 C with time by Sidlichovsky & Nesvorny (1994), but they did
3
©0000RAS,MNRAS000,000–000
10 Shannon et al.
τResonace width e-folding time / orbital period (/P) φp 11 ..246802xxxxxx111111000000555566 2 4 6 8 WAHNOQUVXYBCRLTPJI 10 sest Planet-Dust Grain approach (AU) 000000......123456 ’box.dat’ using 2:1:3 11000
o
j Cl 0 1
0 0.2 0.4 0.6 0.8 1
Figure 10. The e-folding time of the resonance width for Time (resonance escape times)
dust grains in all simulations. We fit the resonance width
as exponential growth between 160(0.01/βPR)(ap/1au)2 kyrs
and 320(0.01/βPR)(ap/1au)2 kyrs for those particles which are Figure 11. Closest approach between the planet and the dust
trapped in resonance for at least 480 (0.01/βPR)(ap/1au)2 kyrs particle in simulation B. Dust grains are binned with all grains
trappedinthesameresonanceinresonanceescapetimebinswith
(all simulations not plotted have zero particles meeting this cri-
width 104 yrs. Closest approach is sampled across 1000 years
terion).PointsarevaluesfromtheN-bodysimulations,thesolid
for each bin. Once the particles come within ∼ 0.04 au, they
lines are equation 43. Here the error bars cover the smallest to
are ejected from the resonance. We measure this escape radius
largestvalues,withthepointcentredonthemedianvalue.Rare
acrossallsimulationsusingthesecondlasttimebin,where10%
catastrophic failures of the automatic fitting produce extremely
ofmeasurementsarebelowtheescaperadius.Onlythoseparticles
largevalues,soweemploythemedianvalue,ratherthanthemean
caughtformorethan104 yearsareincluded.
value.
one another as
not quantify this behaviour. We also observed exponential
growth with time (Figure 2). This growth is of the form (cid:18) (cid:19)2
δϕ=δϕ et/τϕ (42) τϕ ≈3 j+j 1 3 (1−βPR)−23 τPR (44)
0
andherethej andβ termswillbeclosetoone,andthus
where δϕ is characterised in §3.3.3. To fit the reso- PR
0 the e-folding time will always be larger than the PR drag
nance width of a particle in resonance, we bin all out-
timebyafactorofafew.Becauseanyresonantobjectmust
puts of ϕ within a 4000 (0.01/β ) (a/1 au)1.5 year pe-
PR have δϕ < π, unless dust grains are captures at very small
riod (as in §3.3.5), taking the largest difference between
libration widths (which is not found by Mustill & Wyatt
two values as the resonance width at that time. For those
(2011) for any parameter space they investigate), they can
particles that persist in resonance for more than 4.8 ×
only undergo a few e-folds of their libration width before
105(0.01/β )(a/1 au)1.5years,wefitanexponentialcurve
PR they escape the resonance. Therefore, resonantly captured
tothegrowthbetween1.6×105(0.01/β )(a/1 au)1.5 and
PR ringsshouldneverexceedthebackgrounddiskinbrightness
3.2×105(0.01/β )(a/1 au)1.5 yearstoassignaτ tothat
PR ϕ by more than a factor of O(10).
particle’s evolution. The τ of different particles with the
ϕ
sameβ caughtinthesameresonanceclustertightly(Fig-
PR
ure 10).
3.3.7 Escape from resonance
Byfittingallofthevaluesderivedthiswayfordifference
resonances and simulations, we find At early times, dust grains undergo close encounters with
the planet as they are pushed into the resonance. The res-
(cid:18)0.01(cid:19)(cid:18)j+1(cid:19)2(cid:16) a (cid:17)2 onance then protects the dust grain from close encounters,
τϕ ≈1.14×105 β j 1apu years, (43) butastheeccentricityandlibrationwidthincrease,encoun-
PR tersbecomeincreasinglycloser,untiltheparticleislostfrom
afitofwhichcanbeseeninplot10.Fixingtheconstants theresonance.Toanalyseifandwhencloseencounterscause
of equation 43, varying one while keeping the rest fixed, a a dust grain to be lost from resonance, we bin together all
non-linear least-squares Marquardt-Levenberge fit, weight- dust grains from the same simulation which are caught in
√
ing each point by n where n is the number of particle- the same resonance and which are trapped in that reso-
histories represented by a point, gives the normalisation as nance for the same amount of time (rounded to the near-
1.14±0.01×105, the exponent of the (j+1)/j term as est 104(a/1 au)2(m /m )0.5 years). We define the closest
p ⊕
2.02±0.01,theexponentoftheβ termis−0.94±0.05,the approachforeachbinofdustgrainsbythesmallestplanet-
PR
exponent of the planet mass term would be −0.006±0.002 dust grain separation output by the simulation over each
if we included a power law term dependent on the planet’s 103(a/1 au)2(m /m )0.5yearinterval(Figure11).Because
p ⊕
mass,andtheexponentofthesemimajoraxisis1.99±0.02. theplanet-dustgrainclosestapproachesdecreaseduringthe
Note that the scalings in equation 43 are very similar evolution, and escape occurs for all particles at a fixed en-
to the PR drag lifetime in equation 5. They are related to counterdistance,weconcludeescapefromresonanceoccurs
©0000RAS,MNRAS000,000–000