Jumping the Gap: The Formation Conditions and Mass Function of Pebble-Pile Planetesimals

2 MB·

by  Philip F. Hopkins

Mon.Not.R.Astron.Soc.000,000–000(0000) Printed15December 2015 (MNLATEXstylefilev2.2)
Jumping the Gap: The Formation Conditions and Mass Function of
“Pebble-Pile” Planetesimals
Philip F. Hopkins1,2∗
1TAPIR,Mailcode350-17,CaliforniaInstituteofTechnology,Pasadena,CA91125,USA
SubmittedtoMNRAS,January,2014
ABSTRACT
Inaturbulentproto-planetarydisk,dustgrainsundergolargedensityfluctuationsandundertherightcircumstances,
grain overdensities can collapse under self-gravity (forming a “pebble pile” planetesimal). Using a simple model
for fluctuations predicted in simulations, we estimate the rate-of-formation and mass function of self-gravitating
5
1 planetesimal-massbodiesformedbythismechanism.Thisdependssensitivelyonthegrainsize,disksurfacedensity,
0 andturbulentMachnumbers.However,whenitoccurs,theresultingplanetesimalmassfunctionisbroadandquasi-
2 universal, with a slope dN/dM∝M−(1−2), spanning size/mass range ∼10−104km (∼10−9−5M ). Collapse to
⊕
planetesimalthroughsuper-Earthmassesispossible.Thekeyconditionisthatgraindensityfluctuationsreachlarge
c
e amplitudesonlargescales,wheregravitationalinstabilityproceedsmosteasily(collapseofsmallgrainsissuppressed
D byturbulence).Thisleadstoanewcriterionfor“pebble-pile”formation:τ (cid:38)0.05ln(Q1/2/Z )/ln(1+10α1/4)∼
s d
0.3ψ(Q,Z,α) where τs =tsΩ is the dimensionless particle stopping time. In a MMSN, this requires grains larger
3
thana=(50, 1, 0.1)cmatr=(1, 30, 100)au.Thismayeasilyoccurbeyondtheiceline,butatsmallradiiwould
1
dependontheexistenceoflargeboulders.Becausedensityfluctuationsdependstronglyonτ (inverselyproportional
s
to disk surface density), lower-density disks are more unstable. Conditions for pebble-pile formation also become
]
P morefavorablearoundlower-mass,coolerstars.
E
Keywords: planetsandsatellites:formation—protoplanetarydiscs—accretion,accretiondisks—hydrodynamics
.
h —instabilities—turbulence
p
- 1 INTRODUCTION ispossible,ifgrainsaresufficientlylarge,abundant,andthedisks
o obeyotherconditionsone.g.theirgasdensitiesandsoundspeeds.
r Itiswidelybelievedthatdustgrainsformthefundamentalbuilding
Theconcentrationphenomenon(includingso-called“vortextraps”;
st blocks of planetesimals. But famously, models which attempt to Barge&Sommeria1995;Zhu&Stone2014)canoccurviaself-
a formplanetesimalsviathegrowthofdustgrainsinproto-planetary
excitation of turbulent motions in the “streaming” instability (Jo-
[ disksfacethe“meterbarrier”or“gap”:planetesimalsofsizes(cid:38)km
hansen & Youdin 2007), or in externally driven turbulence, such
arerequiredbeforegravityallowsthemtogrowfurtherviaaccre-
3 asthatexcitedbythemagneto-rotationalinstability(MRI),global
v tion, but grains larger than ∼cm tend to shatter rather than stick gravitationalinstabilities,convection,orKelvin-Helmholtz/Rossby
8 whentheycollide,preventingfurthergraingrowth.
instabilities(Dittrichetal.2013;Jalali2013;Hendrix&Keppens
5 Therefore, alternative pathways to planetesimal formation
2014).Thedirectnumericalexperimentshaveshownthatthemag-
4 havereceivedconsiderableattention.Ifadustgrains–evensmall
nitude of these fluctuations depends on the parameter τ =t Ω,
2 ones–couldbestronglyenoughconcentratedat,say,thediskmid- s s
theratioofthegas“stopping”time(friction/dragtimescale)t to
1. plane,theirdensitywouldbesufficienttocausetheregiontocol- the orbital time Ω−1, with the most dramatic fluctuations arousnd
0 lapsedirectlyunderitsownself-gravity,and“jumpthegap”todi-
τ ∼1.Theseexperimentshavealsodemonstratedthatthemagni-
4 rectly form a km-sized planetesimal from much smaller grains – s
tudeofclusteringdependsonthevolume-averagedratioofsolids-
1 whatwecalla“pebblepileplanetesimal”(Goldreich&Ward1973;
to-gas (ρ˜≡ρ /ρ ), and basic properties of the turbulence (such
: forreviewsseeChiang&Youdin2010;Johansenetal.2014).In d g
v astheMachnumber).Thesehaveprovidedkeyinsightsandmoti-
i general,though,turbulenceinthedisksetsa“lowerlimit”tothede- vatedconsiderableworkstudyingtheseinstabilities;however,the
X greetowhichgrainscansettleintoarazor-thinsub-layer,andthis
parameter space spanned by direct simulations is always limited.
r hasgenerallybeenregardedasabarriertothe“pebblepile”sce-
a Moreover, it is not possible to simulate the full dynamic range
nariodescribedabove(butseeLyraetal.2009b;Leeetal.2010;
ofturbulenceinthesesystems:the“topscales”ofthesystemare
Chiang&Youdin2010,andreferencestherein).
λ ∼AU, while the viscous/dissipation scales λ of the turbu-
However,itisalsowell-establishedthatthenumberdensityof max η
lenceareλ ∼m-km(ReynoldsnumbersRe∼106−109).Clearly,
solid grains can fluctuate by multiple orders of magnitude when η
someanalyticmodelsforthesefluctuations–evensimplisticones
“stirred” by turbulence, even in media where the turbulence is
–areneededformanycalculations.
highlysub-sonicandthegasisnearlyincompressible–ithasthere-
fore been proposed that in some “lucky” regions, this turbulent Fortunately, the question of “preferential concentration” of
concentrationmightbesufficienttotrigger“pebblepile”formation aerodynamic particles is actually much more well-studied in the
(seee.g.Braccoetal.1999;Cuzzietal.2001;Johansen&Youdin terrestrialturbulenceliterature.Therebothlaboratoryexperiments
2007; Carballido et al. 2008a; Lyra et al. 2008, 2009a,b; Bai & (Squires&Eaton1991;Fessleretal.1994;Rouson&Eaton2001;
Stone2010b,a,c;Panetal.2011).Thesestudies–mostlydirectnu- Gualtieri et al. 2009; Monchaux et al. 2010) and numerical sim-
mericalsimulations–haveun-ambiguouslydemonstratedthatthis ulations (Cuzzi et al. 2001; Yoshimoto & Goto 2007; Hogan &
Cuzzi2007;Becetal.2009;Panetal.2011;Monchauxetal.2012)
have long observed that very small grains, with stokes numbers
∗ E-mail:[email protected] St ≡ts/te(λη)∼1 (ratio of stopping time to eddy turnover time
(cid:13)c 0000RAS
2 Hopkins
at the viscous scale) can experience order-of-magnitude density
fluctuations at small scales (at/below the viscous scale). Consid- 6 τs=0.001
τ =0.01
erableanalyticprogresshasbeenmadeunderstandingthisregime: τs=0.1
s
demonstrating, for example, that even incompressible gas turbu- 5 τs=1.0
lenceisunstabletothegrowthofinhomogeneitiesingraindensity ) τs=10
fl
(Elperinetal.1996;Elperinetal.1998),andpredictingthebehav- Z
iorofthesmall-scalegrain-graincorrelationfunctionusingsimple /d4
Z
modelsofgaussianrandom-fieldturbulence(Sigurgeirsson&Stu- ρ
δ
art2002;Becetal.2007). g( 3
Thesestudieshaverepeatedlyshownthatgraindensityfluc- lo
tuations are tightly coupled to the local vorticity field: grains are
“flungout”ofregionsofhighvorticitybycentrifugalforces,and 2
collect in the “interstices” (regions of high strain “between” vor-
tices).1Studiesofthecorrelationfunctionsandscalingbehaviorof Q=60, α=10−4
1
higher Stokes-number particles suggest that, in the inertial range 3 2 1 0 1
log(λ/h )
(ignoring gravity and shear), the same dynamics apply, but with 6 d
the scale-free replacement of a “local Stokes number” t/t , i.e.
s e Q=60, α=10−6
whatmattersforthedynamicsonagivenscalearethevorticesof Q=60, α=10−4
thatscale,andsimilarconcentrationeffectscanoccurwheneverthe 5 Q=60, α=10−2
eddyturnovertimeiscomparabletothestoppingtime(e.g.Yoshi- Q=10, α=10−4
moto&Goto2007;Becetal.2008;Wilkinsonetal.2010;Gus- )fl Q=150, α=10−4
tavssonetal.2012).Severalauthorshavepointedoutthatthiscrit- Z 4
/
icallylinksgraindensityfluctuationstothephenomenonofinter- Zd
mittencyanddiscrete,time-coherentstructures(vortices)onscales ρ
δ 3
largerthantheKolmogorovscaleinturbulence(seeBecetal.2009; g(
o
Olla2010,andreferencestherein).Inparticular,Cuzzietal.(2001) l
demonstrateconvincinglythatgraindensityfluctuationsbehavein
2
a multi-fractal manner: multi-fractal scaling is a key signature of λ=h
d
well-tested, simple geometric models for intermittency (e.g. She λ=λ(t =t)
e e s
&Leveque1994).Inthesemodels,thestatisticsofturbulenceare 1
3 2 1 0 1
approximatedbyregardingtheturbulentfieldasahierarchicalcol-
log(τ )
lection of “stretched” singular, coherent structures (e.g. vortices) s
ondifferentscales(Dubrulle1994;She&Waymire1995;Chainais Figure1.Criticalgrainoverdensityδρ≡ρd(λ)/(cid:104)ρd(cid:105)fordynamicalcol-
2006). lapseunderself-gravity(Eq.18).WeplotδρZd/Z(cid:12),sinceitisthecombina-
Suchstatisticalmodelshavebeenwell-testedasadescription tionδρZdthatmustexceedsomecriticalvalueasafunctionofscale(λ/hd),
ToomreQ,stoppingtimeτs=tsΩ,andturbulencestrengthα.Top:Criti-
ofthegasturbulencestatistics(includinggasdensityfluctuations;
caldensityvs.scale,fordifferentgrainsizes(τs)inadiskwith“standard”
seee.g.Burlaga1992;Sorriso-Valvoetal.1999;Budaev2008;She MMSNpropertiesat∼1au(Q=60,α=10−4).Onmostscales,larger
& Zhang 2009; Hopkins 2013b). More recently, first steps have grainsrequiresmallerfluctuationstocollapse,becausetheinitialdustdisk
been taken to link them to grain density fluctuations: for exam- isthinner(higher-density)andresistancefromgaspressureisweaker.Bot-
ple, in the cascade model of Hogan & Cuzzi (2007), and the hi- tom:Criticaldensityvs.τs,evaluatedeitheratthediskscaleheight(hd,near
erarchicalmodelofHopkins(2013c).Thesemodels“bridge”be- wherethecollapseoverdensityδρ(λ)isminimized;solid)orthecharacter-
tween the well-studied regime of small-scale turbulence and that isticscalewheredensityfluctuationsaremaximized(λe(te=ts);dotted).
oflarge,astrophysicalparticlesinshearing,gravitatingdisks.The Bothgenerallydecreasewithτs untilτs∼1.Forsmallgrains(τs(cid:28)1),
key concepts are based on the work above: we first assume that thecriticaloverdensitynearλe(te=ts)(cid:28)hdislargebecauseofturbulent
support.WevaryQandα;thecriticaloverdensitiesincreasewithQ,asex-
grain density fluctuations are driven by coherent velocity “struc-
pected,andwithα(sinceturbulentsupportvs.gravityislarger),thoughthe
tures,”forwhichwecansolveanalyticallytheperturbationowing
lattereffectisweak.
toasinglestructureofagivenscale.BuildingonCuzzietal.(2001)
andothers,wethenattachthiscalculationtoawell-tested,simple,
cascade model for the statistics of velocity structures. In Hogan underwhich“pebblepile”planetesimalformationmayoccur,and
et al. (1999); Cuzzi et al. (2001); Hogan & Cuzzi (2007); Teitler inthatcase,toestimatethemassfunctionofplanetesimalsformed.
etal.(2009)andHopkins(2013c),thesemodelsareshowntogive
agoodmatchtobothdirectnumericalsimulationsandlaboratory 2 THEMODEL
experiments.
2.1 OverviewandBasicParameters
Inthispaper,wecombinetheseanalyticapproximationswith
simplecriteriaforgravitationalcollapse,tocalculatetheconditions Consider a grain-gas mixture in a Keplerian disk, at some (mid-
plane)distancer fromthecentralstarofmassM .Assumethat
∗ ∗
thegrainsareinadiskwithsurfacedensityΣ andexponentialver-
d
tical scale-height h ((cid:104)ρ (z)(cid:105)∝exp(−|z|/h )), so the mid-plane
1 Itissometimessaidthatanti-cyclonic,large-scalevortices“collect”dust d d d
density (cid:104)ρ (cid:105) ≡(cid:104)ρ (z=0)(cid:105)=Σ /(2h ). This is embedded in a
grains.Butitismoreaccuratetosaythatgrainspreferentiallyavoidregions d 0 d d d
withhighmagnitudeofvorticity|ω|.Ifsufficientlylargevorticesareanti- gasdiskwithcorrespondingΣg,Hg,(cid:104)ρg(cid:105)0,andsoundspeedcsand
cyclonicandalignedwiththediskplane,theyrepresentalocalminimumin theglobalgrain-to-gasmassratioisdefinedbyZd≡Σd/Σg.Being
|ω|,sograinsconcentratebybeingdispersedoutofhigher-|ω|regions. Keplerian, the disk has orbital frequency Ω≡(GM∗/r∗3)1/2 and
(cid:13)c 0000RAS,MNRAS000,000–000
JumpingtheGap:Pebble-PilePlanetesimals 3
2001),massfractionZ .Themid-planestoppingtimeis
d
6 τ =0.001
τs=0.01 (cid:40)
τs=0.1 t = ρ¯dRd × 1 (Rd≤9λσ/4) (1)
) 5 τss=1.0 s (cid:104)ρg(cid:105)0cs (4Rd)/(9λσ) (Rd>9λσ/4)
fl τ =10
Z s whereλ =1/(n σ(H ))=µm /((cid:104)ρ (cid:105) σ(H ))isthemean-free
/ 4 σ g 2 p g 0 2
Zd pathinthegas(ngisthegasnumberdensity,mptheprotonmass,µ
δρ themeanmolecularweight,andσ(H2)thecrosssectionformolec-
log( 3 uthlaerdcuosltlissicoanles)h.Weigehcta,nhdth=en(cid:112)deαfi/n(eατs+≡τst)sΩ≈.T(cid:112)hαis/aτnsdHαg,daetgeermneirnael
resultthatholdsforbothlargeandsmallτ (Carballidoetal.2006).
2 s
Nowallowafluctuationρ (k)=δ (cid:104)ρ (cid:105)(δ (cid:54)=1)ofthemean
d ρ d ρ
β=0 graindensityaveragedonthescalekcenterednearthemidplane,
1 where k ≡ 1/λ is the wavenumber (λ the wavelength) of the
3 2 1 0 1
log(λ/h ) fluctuation. For the incompressible (Kolmogorov) turbulent cas-
d cade,weexpectanrmsturbulentvelocityoneachscale(cid:104)v2(k)(cid:105)≡
g
αc2f(λ/λ )with f ∼(λ/λ )2/3whereλ isthetop/driving
6 s t max t max max
scaleofthecascade(wetakeλ ≈H ).3 Wecandefinethecor-
max g
respondingeddyturnovertimet (k)=λ/(cid:104)v2(k)(cid:105)1/2.Notethatby
e g
5 assumingaKolmogorovspectrum,weareimplicitlyassumingthat
)
fl grainsdonotmodifythegasvelocitystructures.Thisistruewhen
Z
/ 4 thedensityfluctuationsareweak,butlessclearwhentheyarelarge
d
Z (so that the local grain density becomes large compared to the
δρ gasdensity).However,preliminarysimulationresultssuggestthe
( 3
g powerspectraofthegasturbulencearenotmuchdifferentinthis
o
l limit(seee.g.Johansen&Youdin2007;Panetal.2011;Hendrix
2 &Keppens2014).
The grains will also have a scale-dependent velocity disper-
β=1
sion following the turbulent cascade, for which we can define
1
3 2 1 0 1 (cid:104)v2(k)(cid:105)≡αc2g(λ/λ ). However, grains are partially-coupled
d s t max
log(λ/hd) togas,sogt isingeneralanon-trivialfunctionwhichwederivein
AppendixA. Onlarge scales(for smallgrains)wheret (cid:28)t (k),
s e
Figure2.AsFig.1(top),butsimplyforcingβ=0(nogaspressure)or
thegrainsarewell-entrainedbythegas,soweexpectg ≈ f,but
β =1 (treating the gas as perfectly-coupled). For β =1, the predicted t t
on small scales wheret (cid:29)t (k), the grains are effectively colli-
thresholdsarenotdramaticallydifferentfromourfullcalculationexceptat s e
sionless,sohaveaconstant(scale-independent)minimumvelocity
intermediatescalesorforlargegrainsonsmallscales.Takingβ=0,how-
dispersion.
ever,wouldleadonetoinfermuchsmallercollapsedensities(byanorder
ofmagnitudeormore)onmanyscales.Onemustaccountforgaspressure Inordertocalculatethemassfunctionof“pebblepiles”–the
resistingthecollapseofgrainsonlargescales.Note,though,thatevenne- numberdensityorprobabilityofforming“interesting”fluctuations
glectinggaspressureentirely,collapseonsmallscales(cid:46)10−3hdrequires –werequirefourthings:
enormousdensityfluctuationsδρ(cid:38)104.
• (1)Amodelfortheproto-planetarydisk.Thisisgivenbythe
descriptionaboveand§2.2
• (2)Amodelforthestatisticsofgraindensityfluctuations.This
epicyclicfrequencyκ≈Ω(KepleriancircularvelocityV =ΩR).
K isoutlinedin§2.4,andisbasedonthedirectnumericalsimulations
In the regime of interest in this paper, the Mach numbers of gas
andexperimentsdescribedin§1.
turbulencewithinthemid-planedustlayeraresmall((cid:28)1;Voelk
etal.1980;Laughlin&Bodenheimer1994;Gammie2001;Hughes
et al. 2011), so the gas density fluctuations are much smaller
than the grain density fluctuations (Passot et al. 1988; Vazquez- 3 Moreaccurately,wecancorrectlyincludetheenergy-containingrange
Semadeni 1994; Scalo et al. 1998) and we can treat the gas as (λ > λmax) by taking the isotropic turbulent power spectrum E(k) ∝
approximately incompressible ρ ≈(cid:104)ρ (z)(cid:105). This also gives the k−5/3(1+|kλmax|−2)−(2−[−5/3])/2(Bowman1996).Thisgives
g g
gas scale height Hg =cs/Ω, and the usual Toomre Q parameter (cid:104)v2g(λ)(cid:105)=αc2s ft(λ/λmax) (2)
Q≡csκ/(πGΣg)=Ω2/(2πG(cid:104)ρg(cid:105)0).Wecandefinetheusualtur- 4Γ[11/6] (cid:90) x
bulentα≡(cid:104)v2(cid:105)/c2,where(cid:104)v2(cid:105)isthermsturbulentvelocityofthe ft(x)≡ √ t−1/3(1+t2)−11/6dt (3)
gasaveragedognthselargestscgalesofthesystem.2 πΓ[1/3] 0
≈1.189x2/3(1+1.831x7/3)−2/7 (4)
WewillfocusonamonolithicgrainpopulationwithsizeR =
d
R cm, internal density ρ¯ ≈2gcm−3 (Weingartner & Draine We use this (the approximate form is accurate to ∼1% at all k) in our
d,cm d
fullnumericalcalculations.Somecutoffisnecessaryatlargescalesorelse
apower-lawcascadecontainsadivergentkineticenergy(andwedonot
expectλmax(cid:29)Hg).However,wedonotincludeanexplicitmodelforthe
2 Westressthatα=(cid:104)v2g(cid:105)/c2s isherepurelytofunctionasausefulparam- dissipationrange(scalesbelowKolmogorovλη)–i.e.weassumeinfinite
eterdefiningtheturbulentvelocities.Wearenotspecificallyassuminga Reynoldsnumber–sinceallquantitiesinthispaperareconvergedalready
Shakura&Sunyaev(1973)-typeviscous“α-disk”noraGammie(2001)- onmuchlargerscales.Fortheconditionsofinterest,λη∼0.1km(Cuzzi
typegravito-turbulentdisk. etal.2008).
(cid:13)c 0000RAS,MNRAS000,000–000
4 Hopkins
• (3)Acriterionforan“interesting”fluctuation.Wetakethisto andre-radiatedemission(inwhichcasetheexternalradiationpro-
beafluctuationwhichissufficientlylargethatitcanundergody- ducesahotsurfacedustlayerwhichre-radiates∼1/2theabsorbed
namicalcollapseunderself-gravity(overcomingresistancetocol- lightbackintothedisk,maintainingT4 ≈T4 /2;seeChiang
mid,∗ eff,∗
lapse from gas drag and pressure, turbulent kinetic energy, shear &Goldreich1997).5
andangularmomentum).Wederivethiscriterionin§2.3andAp-
2.3 CriteriaforDynamicalGravitationalCollapse
pendicesB-C.
• (4) A mathematical method to “count” the interesting fluc- Now,todefinethemassfunctionof“interesting”graindensityfluc-
tuations, given the assumptions above. This is provided by the tuations. Here, we will define “interesting” as those fluctuations
excursion-setformalism,assummarizedin§2.5 whichexceedsomecriticaldensityρcrit,abovewhichtheycancol-
lapseunderself-gravityonadynamicaltimescale.Thisisnotthe
2.2 PhysicalDiskModels onlychannelbywhichdustoverdensitiescanformplanetesimals!
Therearesecularinstabilities(e.g.Youdin2011;Shariff&Cuzzi
In order to attach physical values to the dimensionless quantities
2011), and grain overdensities could promote grain growth; but
above,werequireadiskmodel.Wewilladoptthefollowing,mo-
theserequiredifferentconsiderations(see§6),andareoutsidethe
tivatedbytheminimum-masssolarnebula(MMSN),foradiskof
scopeofourcalculationhere.
arbitrarysurfacedensityaroundasolar-typestar:4
For an grain overdensity or mode with size/wavelength λ=
Ω=(cid:114)GrM∗3∗ ≈6.3ra−u3/2yr−1 (5) 1ρ/crkit,=thρecrcitr(iλtic=alkd−e1n)s.itIytiρsccriotnwvielnlibeentatofudnecfitinoenthoeftlhoactalwgaavse-tloen-dguthst,
massratioaveragedonascaleλaroundapointx(e.g.averagedin
Σ =Σ 1000r−3/2gcm−2 (6)
g 0 au asphereofradiusλaboutthepointx)
Teff,∗=(cid:16)(0.054rra2u∗2/7)R2∗(cid:17)1/4T∗≈140ra−u3/7K (7) ρ˜=ρ˜(x,k)≡1+ρd(ρxg,k) (15)
whereTeff,∗istheeffectivetemperatureofthedisk(Chiang&Gol- Ifweconsidergrainswhicharepurelycollisionless(nograin-
dreich1997);R∗ andT∗ aretheeffectivesizeandtemperatureof gas interaction), then a Toomre analysis gives the following cri-
thestar. terion for gravitational instability of a mid-plane perturbation of
FollowingChiang&Youdin(2010),Eqs.3-15,thesechoices wavenumberk:
determinetheparameters
|kh |
(cid:115)k T 0>ω2=κ2+(cid:104)v2d(k)(cid:105)k2−4πGρgρ˜1+|kdhd| (16)
c = B mid ≈0.64r−3/14kms−1 (8)
s µm au Theω hereisthefrequencyoftheassumed(linear)perturbations
p
(∝exp(−ıωt);seeAppendixB);forω2<0,themodeisunsta-
H c
g = s ≈0.022r2/7 (9) ble.Notethatthisisidenticaltothecriterionforastellargalactic
r V au
∗ K disk (Binney & Tremaine 1987). We show in Appendix B that a
Σ
(cid:104)ρ (cid:105) = g ≈1.5×10−9Σ r−39/14gcm−3 (10) systematicdustsettling/driftvelocitydoesnotchangethiscriterion
g 0 2H 0 au
g significantly,solongasthedriftvelocityv ∼τ ηV (cid:28)V .Here
drift s K K
Q= csΩ ≈61Σ−1r−3/14 (11) theκtermrepresentsthecontributionofangularmomentumresist-
πGΣg 0 au ingcollapse,and(cid:104)v2d(k)(cid:105)isthermsturbulentvelocityofgrainson
1 ∂((cid:104)ρ (cid:105) c2) the scale k; for a derivation of the turbulent term here see Chan-
Π= 2(cid:104)ρ (cid:105) V c ∂gln0r s ≈0.035ra2u/7 (12) drasekhar(1951);Chavanis(2000).6ThenegativeterminGrepre-
g 0 K s
sentsself-gravity,andde-stabilizestheperturbationatsufficiently
1 1.2Σ−1r39/14
λ = ≈ 0 au cm (13)
σ n σ(H ) 1+(r /3.2)3/7
g 2 au
(cid:40) 5 Moreaccurately,wecantaketheeffectivetemperaturefromillumina-
τs≈MAX 00..0000414ΣR−0d12,crma3u/r2a−uR9/d,7c(m1+(rau/3.2)3/7) (14) tainogn&tobGeo:lTder4fef,i∗ch=19T9∗47α).TSRin2∗c/e(4wre∗2)alwloiwthnαoTn-≈ze0ro.0α05,rtha−uis1+im0p.l0ie5sra2au/n7e(fCfehci--
wanidthwµe≈tak2e.3th(eapmporolepcruialtaerfcorrosas-ssoeclatiromniσx(tuHre)o≈f m2×ole1c0u−la1r5(g1a+s) t1i9v7e3v)i,swcohsiictyhapnroddaucccersetainonefrfaetcetMiv˙e≈tem3πpαeract2suΣregTΩe4f−f,1ac(cS≈ha3kuM˙raΩ&2/S(u8nπyσaBev)
(T/70K)−1)(Chapman&Cowling1970).He2reTmid,∗ isthedisk ((cid:112)σBkBiTsmtihd/e(µBomltpz)m.Aanmnocroenasctacnutr)a.teNeostteimtahtiesodfeTpmeinddisstohnenthgeivetenrmbycs2solv=-
mgaisd-cpilracnuelatremvpeeloractiutyre,ananddtΠhedeKfienpelesrtihaenocffirsceutlbaertwveeleoncitthyemVea−n ingtheimplicitequationTm4id=(3/4)[τV+4/3+2/(3τV)]Te4ff,acc+[1+
(cid:104)V (cid:105)≡ηV ,whereΠ≡ηV /c;thisisrelatedtothegraindriftKve- τV−1]Te4ff,∗,whereτV =τV(Tmid)=κR(Tmid)Σg/2istheverticaloptical
locgiatsyv ≡K 2ηV τ [(1+ρK˜)2+s τ2/4]1/2[τ2+(1+ρ˜)2]−1(Nak- depthfromthemidplane.HereκR istheRosselandmeanopacity,which
drift K s s s wecantakefromthetabulatedvaluesinSemenovetal.(2003)(crudely,
agawaetal.1986). κR∼5cm2g−1atTmid>160KandκR∼2.4×10−4T2cm2g−1K−2at
The expression we use for Tmid,∗ is the approximate expres- lowerTmid.Weusethismoredetailedestimateforourfullnumericalcal-
sion for the case of a passive flared disk irradiated by a central culation,howeveritmakesalmostnodifferencefortheparameterspacewe
solar-typestar,assumingthediskisopticallythicktotheincident consider,comparedtothesimplescalingsabove.
6 More exactly, for grains on small scales – where they are locally col-
lisionless–weshouldcombinetheturbulentvelocityanddensityterms,
4 WeassumeastellarmassM∗≈M(cid:12);usingasizeR∗≈1.8R(cid:12)andeffec- takinginsteadρd→ρdF(ω/κ,k2(cid:104)v2d(k)(cid:105)/κ2)whereF isthereduction
tivetemperatureT∗≈4500Kforayoungstar,orR∗≈R(cid:12),T∗≈6000K factordeterminedbyintegrationoverthephase-spacedistribution.How-
foramorematurestargiveidenticalresultsinourcalculations.Variations ever,therelevantstabilitythresholdcomesfromevaluatingF nearω≈0;
intheassumedstellarageleadtopercent-levelcorrectionstothemodel, inthisregimewecanTaylorexpandF (assumingaMaxwellianvelocity
muchsmallerthanourotheruncertainties. distribution),andtoleadingorderwerecoverthesolutioninEq.18.The
(cid:13)c 0000RAS,MNRAS000,000–000
JumpingtheGap:Pebble-PilePlanetesimals 5
largeρ˜.Thetermsin|kh |ontherightaretheexactsolutionforan paperwemakethesimpleapproximation7
d
exponentialverticaldiskandsimplyinterpolatebetweenthetwo- (cid:114)
dcaimseenosniosncaalle(tshi(cid:46)n-dhisk()seceasEelmonegscreaelens1(cid:38)98h7d;aKndimthreete-adl.im2e0n0s2i,onfoarl β≈ ts+tgrtagvrav =(1+ts/tgrav)−1=(cid:16)1+4πτs 3ρ˜Q(cid:17)−1 (19)
d
derivations). where we assume t ≈(3π/32Gρ)1/2, for regions which meet
grav
In the opposite, perfectly-coupled (t →0) limit, we have a
s our dynamical collapse criterion (i.e. this does not apply to
singlefluid,sothedispersionrelationisidenticaltothatofapure,
secularly-sedimentingregions,orregionswhereself-gravityisnot
singlecollisionalfluid,inwhichwesimplydefine“dust”and“gas”
strongerthanallotherforcesincludingthesupportfromgasdrag).
sub-components
Alternativederivationsofthesescalingsfromthelinearequa-
tions for coupled gas-dust fluid, and including the non-linear
1
0>ω2=κ2+ρ˜(c2s+(cid:104)v2g(k)(cid:105))k2 (17) stochastic effects of turbulence, are presented in Appendices B-
C, respectively. If anything, we have chosen to err on the side of
ρ˜−1 |kh |
+ (cid:104)v2(k)(cid:105)k2−4πGρ ρ˜ d cautionanddefineastrictcriterionforcollapse–almostallhigher-
ρ˜ d g 1+|kh |
d order effects make collapse slightly easier, not harder. This crite-
=κ2+c2sk2+(cid:104)v2(k)(cid:105)k2−4πGρ ρ˜ |khd| rionissufficienttoensurethat(atleastintheinitialcollapsephase)
ρ˜ g 1+|kh | apebblepileisgravitationallybound(includingthethermalpres-
d
sure of gas being dragged with dust and turbulent kinetic energy
Here cs and (cid:104)v2g(k)(cid:105) represent gas pressure and turbulent support ofgasanddust),andgravitationalcollapseissufficientlystrongto
(and we used (cid:104)v2g(k)(cid:105)=(cid:104)v2d(k)(cid:105) for the perfectly-coupled case). overcomegaspressureforces,tidalforces/angularmomentum/non-
Notethatthetermsdescribingthegaspressure/kineticenergyden- linearshearingoftheoverdensity,turbulentvorticityand“pump-
sity have a pre-factor 1/ρ˜ =ρg/(ρg+ρd), since what we need ing”oftheenergyandmomentumintheregion,andram-pressure
forthemixed-grain-gasperturbationistheenergydensityperunit forcesfromthe“headwind”owingtoradialdrift.Similarly,when
mass in the perturbation. Likewise, the grain kinetic energy den- this criterion is met, the gravitational collapse timescale is faster
sity term has a pre-factor (ρ˜−1)/ρ˜=ρd/(ρg+ρd). Since both thantheorbitaltime,thegraindrifttimescale,theeffectivesound-
sitinthesameexternalpotentialandself-gravitateidentically,the crossingtimeoftheclump,andtheeddyturnovertime.
κ and G terms need no pre-factor. In this limit we can think of Usingthedefinitionsin§2.1,Eq.18canbere-written:
the 1/ρ˜ factor as simply an enhanced “mean molecular weight”
fromtheperfectly-draggedgasgrains(sotheeffectivesoundspeed 0>1+(cid:16)Hg(cid:17)2 1 (cid:104)β(cid:16)1+αf(cid:2) λ (cid:3)(cid:17) (20)
of the gas ceff →c/(cid:112)1+ρ /ρ =c/ρ˜1/2). In other words, in hd λ˜2ρ˜ t λmax
s s d g s
theperfectly-coupledlimit,thesystembehavesasgaswhichmust (cid:2) λ (cid:3)(cid:105) 2Q−1
+α(ρ˜−1)g − ρ˜
“carry”someextraweightindust(seealsoYoudin2011;Shariff& t λmax 1+λ˜
Cuzzi2011,2014). τ (cid:104) (cid:16) (cid:17) (cid:105) 2Q−1
Wecaninterpolatebetweenthesecasesbywriting =1+λ˜2sρ˜ β α−1+ft +(ρ˜−1)gt − 1+λ˜ ρ˜ (21)
0>ω2≡κ2+βρ˜(c2s+(cid:104)v2g(k)(cid:105))k2 whereλ˜ ≡λ/hd andweabbreviate ft = ft(λ/λmax).Thishasthe
solution
ρ˜−1 |kh |
+ ρ˜ (cid:104)v2d(k)(cid:105)k2−4πGρgρ˜1+|kdh | (18) ρ˜>ρ˜crit(λ)≡ψ0(1+(cid:112)1+ψ1/ψ0) (22)
d
Q (cid:104) (cid:105)
Theonlyimportantambiguityintheaboveisthetermβ,whichwe ψ0≡ 4 (1+λ˜) 1+τsλ˜−2gt (23)
introduce(inaheuristicandadmittedlyad-hocmanner)torepre- 2τ (cid:104) (cid:105)(cid:104) (cid:105)−1
sentthestrengthofcouplingbetweengrainsandgas(β=0isun- ψ1≡ λ˜2s β(α−1+ft)−gt 1+τsλ˜−2gt (24)
coupled/collisionless;β=1isperfectly-coupled).Ingeneral,βis
(Note, if β itself is a function of ρ˜, then this is an implicit equa-
someunknown,presumablycomplicatedfunctionofalltheparam-
tion for ρ˜ which must be solved numerically). Recall the di-
etersabove,whichcanonlybeapproximatedinthefullynon-linear crit
mensionless grain density fluctuation δ =ρ /(cid:104)ρ (cid:105) , so ρ˜=1+
casebynumericalsimulations(andtheexactcriterionforinterme- ρ d d 0
(cid:112)
δ ((cid:104)ρ (cid:105) /(cid:104)ρ (cid:105) )=1+δ (Σ /Σ )(H /h )=1+δ Z τ/α.So
diate cases between these limits may require terms beyond those ρ d 0 g 0 ρ d g g d ρ d s
intermsofδ ,thecriterionbecomes
whichcanbeapproximatedbytheβhere).However,thelimitsare ρ
straightforward:ifaperturbationcollapsesonafree-falltimetgrav, (cid:114)α(cid:16) (cid:17)
andtgrav(cid:28)ts,weexpectβ→0(sincethereisnotimeforgasto δρZd> τ ρ˜crit(λ)−1 (25)
s
decelerategrains).Converselyift (cid:29)t,β→1.Thereforeinthis
grav s
2.4 ASimpleRepresentationofGrainDensityFluctuations
inIncompressibleGas
Hopkins(2013c)presentsomesimple,analyticexpressionsforthe
exact solution can be determined for the purely collisional limit (identi-
caltoEq.18)orthepurelycollisionlesslimit(identicaltoastellardisk, statisticsofgraindensityfluctuationsinaturbulentproto-planetary
wheretheminimumdensityforcollapseρ˜issmallerbyafactor=0.935.
Giventheotheruncertaintiesinourcalculation,thisdifferenceisnegligible.
Moreover,inAppendixC,weshowthattheformoftheturbulenttermsin 7 ThisapproximationismotivatedbytheMaxey(1987)linearexpansion
Eq.16accountsforthenon-linear,time-dependentandstochasticbehavior ofthesolutionforthede-celerationofdustgrainsbymolecularcollisions
oftheturbulence(i.e.accountsforfluctuationsinthevelocitydispersions, fortimest(cid:28)tsinasymmetric,homogeneoussphere.Theratioofthedust
andrepresentsthecriterionforaregionwheretheprobabilityofsuccessful de-celerationtermtothetermfromgaspressureifthespherewerepuregas
collapseislarge). (ρ−1∇P)inthislimitisthesameasβinEq.19.
(cid:13)c 0000RAS,MNRAS000,000–000
6 Hopkins
disk.Forourpurposeshere,whatisimportantisthattheseexpres- 2.5 Counting“Interesting”DensityFluctuations
sionsprovideareasonable“fittingfunction”totheresultsofdirect
For any model for the statistics of grain density fluctuations, and
numerical simulations and laboratory experiments studying grain
any threshold criterion for an “interesting” fluctuation (both as a
density fluctuations resulting from a variety of underlying mech-
functionofscale),thereisawell-definedmathematicalframework
anisms (e.g. driven turbulence, zonal flows, streaming-instability
forcalculatingthepredictedmassfunction,sizedistribution,cor-
andKelvin-Helmholtzinstabilities;Johansen&Youdin2007;Bai
relationfunction,andrelatedstatisticsoftheobjects/regionswhich
&Stone2010c;Dittrichetal.2013;Jalali2013;Hendrix&Kep-
exceedthethreshold.Thisisthe“excursionsetformalism,”well-
pens2014).Weshouldnotethatthisdirectlybuildsonseveralpre-
known in cosmology as the “extended Press Schechter” method
viousanalyticmodelsforgrainclusteringaroundtheKolmogorov
bywhichdarkmatterhalomassfunctions,clustering,andmerger
(turbulentdissipation)scaleandintheinertialrange(e.g.Elperin
histories can be analytically calculated (there, the statistics are
et al. 1996; Elperin et al. 1998; Hogan & Cuzzi 2007; Bec et al.
givenbytheinitialGaussianrandomfieldandcosmologicalpower
2008;Wilkinsonetal.2010;Gustavssonetal.2012).Wehaveex-
spectrum, and “interesting” regions are those which turn around
perimented with some of these models and find, for small grains
from the Hubble flow; see Bond et al. 1991). Recently, the same
wheretheclusteringoccursonsmallscales(wheret (cid:28)Ω−1,and
e frameworkhasbeenappliedtopredictthemassfunctionofstruc-
shear/rotationcanbeneglected),theygivequalitativelysimilarre-
tures (e.g. giant molecular clouds and voids) formed on galactic
sults;wediscussthisfurtherin§5.Howeverthesepreviousmodels
scalesbysuper-sonicinterstellarturbulence(Hopkins2012a);the
didnotconsiderthecaseofarotatingdiskwithlargegrainswhere
initial mass function (Hennebelle & Chabrier 2008, 2009; Hop-
t ∼Ω−1,inwhichcasethebehaviorcandifferdramatically(see
s kins2012b,2013d)andcorrelationfunctions/clustering(Hopkins
e.g.Lyraetal.2008).
2013e)ofcoresandyoungstarsinsidemolecularclouds;andthe
massspectrumofplanetswhichcanformviadirectcollapseintur-
bulent,low-Qdisks(Hopkins&Christiansen2013).Forreviews,
We briefly describe the model in Hopkins (2013c) here, but
seeZentner(2007);Hopkins(2013a);Offneretal.(2013).
interestedreadersshouldseethatpaper.Fundamentally,itfollows
Therearemanywaystoapplythismethodology.Probablythe
Hogan & Cuzzi (2007) in assuming grain density fluctuations on
simplest,andwhatweusehere,isaMonteCarloapproach.Con-
different scales can – like gas density fluctuations in supersonic
sidersomeannulusinthediskatradiusr .Selectsomearbitrarily
∗
turbulence (Hopkins 2013b) – be represented by a multiplicative
largenumberofMonteCarlo“samplingpoints”;eachoftheserep-
randomcascade.Atanyinstant,agroupofgrains“sees”agasvor-
resentsadifferentrandomlocationxwithintheannulus(i.e.they
ticityfieldwhichcanberepresentedasasuperpositionofcoherent
randomlysamplethevolume)–really,eachisadifferentrandom
velocitystructuresor“eddies”withawiderangeofcharacteristic
realizationofthefield,givenitsstatistics.Now,wecanaskwhat
spatial scales λ and timescales t . If we consider a single, ideal-
e thedensityofdustgrainsis,averagedinspheresofsizeλ,around
izededdyorvortex,wecananalyticallysolvefortheeffectithas
eachofthesepoints(ρ (x,λ)).Ifthis“initial”λ=λ issufficiently
d 0
onthegraindensitydistributioninandarounditself(assumingthe
large:
vortex survives for some finite timescale ∼t ). When t is much
e s
smallerthante,thegrainsaretightly-coupledtothegas,sothevor- λ=λ0(cid:29)Hg (26)
texhasnoeffectontheaveragegraindensitydistribution(sincethe ρ (x,λ )→(cid:104)ρ (x,λ )(cid:105) (27)
d 0 d 0
gasisincompressible);whent ismuchlargerthant ,thevortexis
s e forallx,i.e.allpointshavethesamemeandensityaroundthemon
unabletoperturbthegrains.Butwhent ∼t ,thevorteximprints
s e sufficientlylargescales(thisisjustthedefinitionofthemeanden-
large(order-unity)changesinthedensityfield.Thesechangesare
sity,afterall).Now,takeadifferentialstep“down”inscale–this
multiplicative,sototheextentthatthevorticityfieldcanberepre-
correspondstoshrinkingthesmoothingspherebysomeincrement
sentedbyhierarchicalcascademodels,thegraindensitydistribu-
∆λ,andaskwhatthemeandensityinsideeachsphereis.Forthe
tiononvariousscalesbehavesasamultiplicativerandomcascade.
spherearoundeachpointx,thisisgivenbytheappropriatecondi-
AssumingtheturbulenceobeysaKolmogorovpowerspectrum,and
tionalprobabilitydistributionfunctionP(δ ,λ −∆λ|δ [λ ],λ )
assumingsome“fillingfactor”ofstructureswhicheachbehaveas ρ 0 ρ 0 0
(essentially the probability of a given change in the mean den-
scaledversionsoftheidealvortex(constrainedtomatchthatpower
sity,betweentwovolumesseparatedbysomedifferentiallysmall
spectrum),witheachgrainencounteringarandomGaussianfield
smoothingsize),or:
ofvortexstructuresovertime,andadoptingasimpleheuristiccor-
rectionforthe“back-reaction”ofgrainsongas,thisleadstoapre- λ→λ −∆λ (28)
0
dictionforalognormal-like(randommultiplicative)distributionof δ (x,λ )→δ (x,λ −∆λ)=δ (x,λ )+∆δ (x) (29)
ρ 0 ρ 0 ρ 0 ρ
localgraindensities.
P(∆δ )=P[∆δ |δ (λ ),λ ,∆λ] (30)
ρ ρ ρ 0 0
The conditional probability distribution function
Quantitatively,foragivensetofglobal,dimensionlessprop- P[∆δρ|δρ(λ0),λ0,∆λ] is directly related to the power spec-
erties of the turbulent disk: τ, α, and Π, the model predicts the trumofdensityfluctuations(seeHopkins2013a,Eq.2-12),andis
s
distribution of grain density fluctuations (relative to the mean), determinedbythemodelforP(δρ,λ)describedin§2.4.Knowing
δρ(x,λ)=ρd(x,λ)/(cid:104)ρd(cid:105)0.Recall,thisisthegraindensityaveraged thatdistribution,wedrawarandomvalueof∆δρ foreachMonte
within a radius λ (as ρd(x,λ)=Md(|x(cid:48)−x|<λ)/(4πλ3/3)= Carlo point x, determining δρ(x,λ0−∆λ). We then repeat this
δ(λ)(cid:104)ρ (cid:105))aroundarandompointinspacex.Usingtheseparame- untilwereachλ→0(makingsuretotakesmallenoughsteps∆λ
d
ters,weobtainP(δ ,λ),theprobabilitythatanypointinspacelives sothatthestatisticsareconverged).
ρ
withinaregion,averagedonsizescaleλ,withgrainover-density Foreachrandompointx,wenowhavethevalueoftheden-
betweenδρ andδρ+dδρ.Thisisapproximatelylog-normal,with sity field smoothed on all scales, δρ(x,λ) – in the excursion-set
avariancethatdependsonscale(wheremostofthepower,orcon- language, this is referred to as its “trajectory.” We can now sim-
tribution to the variance, comes from scales where the “resonant plycomparethistothepredictedcollapsethresholdoneachscale
condition”t ∼t issatisfied).
s e
(cid:13)c 0000RAS,MNRAS000,000–000
JumpingtheGap:Pebble-PilePlanetesimals 7
ρ˜ (λ) (Eq. 25), to ask whether the region is “interesting” (ex- Reynolds-number flow, the fluctuations are approximately self-
crit
ceeds the critical density for dynamical collapse). To avoid the similar,becauseallgrains“see”alarge,scale-free(power-law)tur-
ambiguityof“double-counting”or“cloudsinclouds”(i.e.trajec- bulent cascade at both larger scales (t (cid:29)t) and smaller scales
e s
tories which exceed ρ˜ (λ ) but also have some λ >λ where (t (cid:28)t). As shown in many studies (Cuzzi et al. 2001; Yoshi-
crit 1 2 1 e s
theyexceedρ˜ (λ ),whichrepresentssmallerscalesthatareinde- moto&Goto2007;Hogan&Cuzzi2007;Panetal.2011;Hop-
crit 2
pendently self-gravitating/collapsing embedded in larger collaps- kins 2013c), the maximum local density fluctuations in this limit
ing regions), we specifically consider the “first crossing distribu- saturateatvaluesδmax∼300−1000.
ρ
tion” (see Bond et al. 1991; Hopkins 2012b,a). Namely, if a tra- Since t (λ < λ ) ∝ λ2/3, this “resonance” will occur at
e max
jectoryexceedsρcrit(λ)anywhere,weuniquelyidentifythelargest scales λ≈λmaxτs3/2α3/4(Hg/λmax)3/2 (cid:28)λmax (so λ˜ ∼α1/4τs2).
size/mass scale λ=λfirst on which ρ˜(x,λ)>ρ˜crit(λ) as the “to- Wecan,onthesescales,alsoapproximate ft≈gt≈(λ/λmax)2/3≈
tal” collapsing object. Since the trajectory δρ(x,λ) is continuous α1/2τs,anddrophigher-ordertermsinλ/λmaxorλ˜.Ifwetakeei-
in λ, ρ˜(x,λfirst)=ρ˜crit(λfirst), and there is actually a one-to-one therthetightlycoupled(β=1)orun-coupled(β=0)limits,we
mapping between the first-crossing scale and mass enclosed in a obtain
first-crossing,givenbytheintegralovervolumeinanexponential
disk(sincethatistheverticalprofileweassumed): (cid:16) Q (cid:17)1/2
M(λ )≡4πρ˜ (λ )(cid:104)ρ (cid:105) h3  2α3/2τs3 (β=1)
first crit first g 0 d ρ˜ ∼ (33)
crit
×(cid:104)λ2first+(cid:16)1+λfirst(cid:17)exp(cid:16)−λfirst(cid:17)−1(cid:105) (31) Q(cid:16)1+τ−2(cid:17) (β=0)
2h2 h h 2 s
d d d
(Hopkins2012a).Itiseasytoseethatonscalesλ <h ,thisis or
first d
justM=(4π/3)ρ λ3 ,onscalesλ >h ,justM=πΣ λ2
(whereΣ =2hcρrit fi)r.st first d crit first Z−1τ−2α−1/4(Q/2)1/2 (β=1)
crit d crit  d s
Finally, we can use our ensemble of trajectories to define δ (cid:38) (34)
ρ
thefunction f = f(λfirst),where f isthefractionofMonteCarlo Z−1τ−5/2α1/2(Q/2) (β=0)
“trajectories” that have a first-crossing on scales λ>λ . Since d s
first
M(λ ) is a function of λ , we can just as well write this
first first This requires extremely large density fluctuations: for
as a function of mass, f = f(M), where M ≡ M(λ ). Now,
first Z ∼ Z , and Q ∼ 60 (MMSN at r ∼ 1au), this gives
since each Monte Carlo trajectory represents the probability that d (cid:12) ∗
minimum δ of ∼ 3 × 105(τ/0.1)−2(α/10−4)−1/4 and
a random point in space – i.e. a random differential volume el- ρ s
∼5000(τ/0.1)−5/2(α/10−4)1/2,respectively.
ement – is embedded in such a region, the differential value s
Physically, even if we ignore gas pressure, and the density
|df(M)/dlnM|dlnM represents the differential volume fraction
fluctuationissmall-scale(soshearcanbeneglected),grainsmust
embeddedinsideofregionsofwithmassesM=M(λ )between
first still overcome their turbulent velocity dispersion in order to col-
lnMandlnM+dlnM.Sincethesefirst-crossingregionshavemean
lapse.AsimpleenergyargumentrequiresGM2(<λ)/λ(cid:38)M (<
internalmassdensity(bydefinition)ρ=ρ (λ ),thenumberof d d
crit first λ)(cid:104)v2(λ)(cid:105) (where M is the dust mass inside the region of size
independent“regions”or“objects”(perunitvolume)mustbe d d
λ); using M (< λ) ∼ ρ λ3 and (cid:104)v2(λ)(cid:105) ∼ (λ/td)2, this is just
d d d e
dn (M) ρ (λ [M])(cid:12)df(M)(cid:12) Gρ (cid:38)(td)−2.Inotherwords,thecollapsetimet ∼(Gρ )−1/2
first ≡ crit first (cid:12) (cid:12) (32) d e grav d
dlnM M (cid:12) dlnM (cid:12) mustbeshorterthantheeddyturnovertime(withinthegrains)ted
on the same scale. But recall, the clustering occurs characteristi-
Andthisisthedesiredmassfunctionofcollapsingobjects.Toturn
callyonascalewhereforthegas,t ∼t.Thus,thegrainsareat
thisintoanabsolutenumber(insteadofanumberdensity),wesim- e s
leastmarginallycoupled,andthegraintd ∼t ∼t –thesameed-
plyneedtointegrateoverthe“effective”volume(differentialvol- e e s
dies that induce strong grain clustering necessarily induce turbu-
umeinaradialannulusdRisjust(2h )(2πRdR));orwecandi-
d
lentgrainmotionswitheddyturnovertimeonthesamescale∼t
rectlyconvertfromvolumefractiontoabsolutenumberbasedon s
(see Bec et al. 2009). So collapse of even a “collisionless” grain
theargumentabove,foranassumeddisksize.
Forourcalculationshere,wetypicallyuse∼109MonteCarlo populationrequirestgrav(cid:46)ts.UsingQ∼Ω2/(Gρg)andρd∼ρgρ˜
(forρ˜(cid:29)1),weseethisisequivalenttotheβ=0criterionabove.
“trajectories”tosamplethestatistics,andsamplethosetrajectories
Since,inthislimit,t <t,takingβ=0isinfactagoodapproxi-
inlogarithmically-spacedsteps∆λ≈0.001λ.Mostoftheresults grav s
mation(andsincetheβ=1criterionrequiresastillhigherdensity,
hereconvergeatmuchcoarsersampling,butthemassfunctionat
sot (cid:28)t,itisnottherelevantcaselimithere).
thelowestmassesrequiresalargenumberoftrajectoriestobeprop- grav s
Thus even with no gas pressure effects
erlyrepresented.
(β = 0), collapse (δmax (cid:38) δcollapse) requires τ (cid:38)
ρ ρ s
0.2(α/10−4)1/5(δmaxZ /1000Z )−2/5(Q/60)2/5 – unless
3 APPROXIMATEEXPECTATIONS ρ d (cid:12)
the disks are extremely quiescent (α(cid:28)10−7), we are forced to
Wenowhaveeverythingneededtocalculatethedetailedstatistics considerlargegrains(whereτ (cid:28)1isnottrue).
s
ofcollapsingregions.Beforewedoso,however,wecangaincon- Beforegoingon,however,notethattheargumentswemake
siderableintuitionusingthesomesimpleapproximations. above apply only to dynamical collapse of small grains. Secular
collapse of small grains, through the slow, nearly incompressible
3.1 SmallGrains
“sedimentation” mode described in e.g. (Youdin 2011; Shariff &
Forsmallgrains(τ (cid:28)1),densityfluctuationsonscales∼h are Cuzzi2011)maystillbepossibleinthisregime.Asnotedabove,
s d
weak (since the grains are well-coupled to gas on these scales). thisrequiresadifferenttreatmententirelyandisoutsidethescope
Large density fluctuations are, however, still possible on small of this paper; however, it may present an alternative channel for
scales,wheret ∼t.Considerthislimit.Inthisregime,inalarge planetesimalformationifonlysmallgrainsarepresent.
e s
(cid:13)c 0000RAS,MNRAS000,000–000
8 Hopkins
Figure3.Predictedmassfunctionofcollapsing(self-gravitating)pebble-pileplanetesimalsformedbyturbulentgraindensityconcentrations.Weplotthe
cumulativenumberformedatvariousradialdistancesfromthestar(perunitorbitaldistance:dN/dlnrau),asafunctionofmass(inEarthmasses).Thediskis
ourstandardMMSNmodel(Z=Z(cid:12),Σ0=1;see§2.2),withα=10−4.Differentlinetypesassumethegrainmassisconcentratedingrainsofdifferentsizes
(aslabeled).Ifthegrainsarelarge(10cm),thenpebblepilescancollapsedirectlytomassesfrom∼10−8−1M⊕overarangeoforbitalradii∼0.1−20au.
Ifgrainsonlyreach1cm,thelowerτssuper-exponentiallysuppressesthisprocessatsmallerradii,anditcanonlyoccuratlargeradii(cid:38)20−30au,where
τs(cid:38)0.1(howevertherangeofmassesattheseradiiislarge,from∼10−4−10M⊕).Formaximumgrainsizes=1mm,thisispushedoutto(cid:38)100au.
√
3.2 LargeGrains thatthe“effective”soundspeedofthecoupledfluidis∼ceff/ ρ˜
s
(Safronov&Zvjagina1969;Marble1970;Sekiya1983),andthat
Forlargegrains,fluctuationsarepossibleonlargescales.Foraflat
τ/α=(H /h )2,c ∼ΩH ,andQ∼Ω2/Gρ .Thenweseethis
perturbationspectrum,themostunstablescaleisλ∼h (Goldreich s g d s g g
d criterion is equivalent to t (cid:46)t ≡h /ceff, i.e. that the col-
&Lynden-Bell1965;Toomre1977;Lau&Bertin1978;Laughlin grav cross d s
lapsetimeisshorterthantheeffectivesound-crossingtimeonthe
&Bodenheimer1994),sotakethislimitnow.Inthiscase f ≈g ≈
t t scaleh .Forτ (cid:46)1,thesegenerallydoallowt (cid:38)t,soβ∼1
(α/τ)1/3andλ˜ ≈1,giving d s grav s
s is the more relevant limit – but importantly, collapse of the two-
(cid:16)Qατs(cid:17)1/2 (β=1) oflvueidrdmenesdiituymcaenvefnoromntoinmessucfafilceisen(cid:29)tlytsliasrgaelloswcaelde,sp(rio.ev.idceodllaaplsaergies
ρ˜ ∼ (35)
crit “slow”comparedtothestoppingtime,but“fast”comparedtothe
Q(1+τ2/3α1/3) (β=0) dynamicalandeffectivesound-crossingtimes).Inotherwords,the
s
importantcriterionisthe“effective”Jeansnumberofthecoupled
or
dust-gasfluid,
Z−1Q1/2 (β=1)
δ (cid:38) d (36) J≡ ρG = (ρd+ρg)2G ∼ ρ˜2(cid:104)ρg(cid:105)0Gλ2 ∼ tcerfofss (cid:46)1 (37)
ρ c2 k2 ρ c2k2 c2 t
Zd−1Q(cid:112)α/τs (β=0) ρ˜s,ef(fλ)∼ gcss ∼Q1/2 css grav (38)
EvenatZd∼Z(cid:12) andQ∼60,thisgivesaminimumδρ of∼400 crit λ(cid:112)G(cid:104)ρg(cid:105)0 λΩ
a“neadsi∼er”10w0h(eτns/g0r.a1i)n−s1c/a2n(αin/d1u0c−e4fl)u1/c2t,uaretisopnesctoivnellayr.gCeoslclaalpesse. isfar ρ˜crit(λ∼hd)∼Q1/2(cid:112)α/cτsH Ω ∼(cid:18)Qατs(cid:19)1/2 (39)
Inthislimit,theβ=0criterionisjusttheaRochecriterion, s g
t (cid:46)Ω−1(theturbulenceissub-sonic,soitssupportisnotdom- (wherewehavedroppedtheorder-unityprefactors).Thiswasfirst
grav
inant on large scales). The β =1 criterion is more subtle: recall proposedasakeyrevisiontotheGoldreich&Ward(1973)mid-
(cid:13)c 0000RAS,MNRAS000,000–000
JumpingtheGap:Pebble-PilePlanetesimals 9
planedensitybySafronov&Zvjagina(1969)andMarble(1970),
and appreciated by Sekiya (1983) (see also Youdin 2011; Shariff
& Cuzzi 2011); recently, direct numerical simulations by Shariff
& Cuzzi (2014) following the full non-linear collapse of a dusty
sphere (no angular momentum or turbulent terms) in the perfect-
coupling(τ →0)limithaveexplictlyconfirmedthatforJ(cid:38)0.4,
s
collapsewilloccurandproceedsonthedynamical(free-fall)time,
whileforsmallerJonlythe“secular”sedimentationmodesurvives
(seealsoWahlbergJansson&Johansen2014,whoobtainconsis-
tentresultsintheperfectlyun-coupledlimit).Thisagreeswellwith
ourcriterion,ifwetaketheappropriatelimits(andkeeptherelevant
pre-factors).
In Hopkins (2013c), approximate expressions are given
for the maximum density fluctuations seen in simulations of
large grains on large scales (Table 2 therein), which pro-
vide a good fit to the results of numerical simulations
of MHD (zonal-flow), streaming-instability, and driven turbu-
lence (Hogan & Cuzzi 2007; Bai & Stone 2010b; Dittrich
et al. 2013; Hendrix & Keppens 2014). These are lnδmax ∼
ρ
6δ (1+δ )−1ln[1+b (1+δ )+b3/2((cid:112)1+δ2−1)]whereδ ≈
0 0 d 0 d 0 0
3.2τ/(1+τ2)andb ∼1dependsontheratioofdrifttoturbulent
s s d
velocities.Forτ (cid:46)1,thisbecomesδmax∼exp[20ln(1+b )τ];
s ρ d s
comparing this to the above (β = 1) criterion requires τ (cid:38)
s
0.05ln(Q1/2/Z )/ln(1+b ) ∼ 0.3. So sufficiently large grains
d d
canindeedachievethesefluctuations.
4 NUMERICALRESULTS
Nowweshowtheresultsofthefullnumericalmodeldescribedin
§2,forspecificchoicesofthediskparameters.
4.1 TheCollapseThreshold
Figure4.Predictedpebblepilemassfunction,asFig.3,forvariedα=
4.1.1 DependenceonSpatialScale:LargeScalesareFavored
10−6−10−2(topandbottom).Atlargemasses,αhaslittleeffectonthe
InFig.1weillustratehowthethresholdforself-gravityderivedin MF.Atlowmasses,increasingαmeanslargerturbulentsupportonsmall
scales,suppressinglow-masspebblepileformation.Forhigh-α,thisflattens
§2.3scalesasafunctionofvariousproperties.Recall,thecombi-
theMFslopeandeliminatespebblepileformationatsomesmallerradii
nationδ Z mustexceedsomevalue(Eq.25)whichisafunction
ρ d even for large grains. However, most of our conclusions about the radii
onlyof(τ,Q,α)inorderforanover-densitytocollapseonady-
s wherepebblepileformationcanoccur,asafunctionofgrainsize,arenot
namicaltimescale.Sothecollapsethresholdindimensionlessunits
changed.
of grain-density fluctuations (δ ) scales inversely with the dust-
ρ
to-gas mass ratio Z . We see that, as is generic for Jeans/Toomre
d
collapse and expected from the arguments in § 3, higher over- scalingδ ∝τ−5/2withweakresidualdependenceonα(andalso
densitiesarerequiredforcollapseonsmallscales,withaminimum ρ s
∝Q),asexpectedfromourderivationabove.
inδ aroundλ∼h .Onsmallscalesthethermalpressuretermin
ρ d
Eq.18(∝λ−2/ρ˜)dominatesthesupportvs.gravity(∝ρ˜),giving
4.1.3 ImportanceofGasPressure
ρ˜ ∝λ−1.Onlargescalesλ(cid:29)h angularmomentumdominates
crit d
Eq.18and,justasintheToomreproblem,ρ˜ ∝λ. InFig.2werepeatthisexercisebutsimplyforceβ=1orβ=0.
crit
Wecanseethateitherapproximationfailstomatchourusual“hy-
4.1.2 DependenceonGrainProperties
brid”interpolatedmodel,atsomerangeofscales,byaboutanorder
We also see that, generically, larger grains (larger τ) require ofmagnitude.Thisindicatesthatitisclearlynecessarytounder-
s
smallerδ forcollapse.Thisisbecause(withotherdiskproperties standthenon-linearbehaviorofcollapsingobjects,whenτ ∼t .
ρ s grav
fixed)theinitialdustdisksettlestoasmallerscaleheight(larger However, assuming β =1 does not much change the criteria for
density),andbecausetheresistancebygaspressureisweaker.The large-scalecollapse,andwhilethechangeatsmallscalesislarge
changeinthisbehaviorforlargegrainsτ (cid:38)1onsmallscalesowes werequiresuchlargevaluesofδ thatdynamicalcollapseisun-
s ρ
to the fact that the velocity dispersions of large grains de-couple likely in any any case. But assuming β =0 gives lower collapse
fromthegasandbecomescale-independent(donotdecreasewith thresholdsbyanorder-of-magnitudeormoreonlargescales(the
λ)onsmallscales. mostinterestingrangeforourcalculation).Inthisregimethecol-
If we focus on δ around scales λ∼h or λ∼λ (t =t), lapsethresholdsaresuchthatthecollapsetimeislongerthanthe
ρ d e e s
asin§3,weconfirmourapproximatescalingsabove.Near∼h , stoppingtime,soitisprobablynotagoodapproximationtoneglect
d
collapserequiresmodestover-densities∼100−1000,weaklyde- gaspressure.Thecollapsingdustwilldraggaswithit,increasing
pendentonτ orα(forsmallα(cid:46)10−3)and∝Q1/2,confirming thegaspressure(Youdin2011;Shariff&Cuzzi2011,2014),hence
s
our approximate scaling for β =1 (since in this limit, t (cid:38)t, our assumption of β =1, but the full non-linear behavior in this
grav s
β ∼1). Around λ∼λ (t =t), we see, as expected, a strong regimeremainspoorlyunderstood.
e e s
(cid:13)c 0000RAS,MNRAS000,000–000
10 Hopkins
Figure5.Predictedpebblepilemassfunction,asFig.3,forvariedprotoplanetarydiskproperties.Herewefixα=10−4,butvarythediskmass/surface
density(proportionaltoΣ0≡Σ/ΣMMSN)andmetallicityZd.Left:Verylow-density(buthigh-Zddisk);thiscorrespondstoaMMSNwhichhaslost∼99%
ofitsgasbutonly∼90%ofitslargegrains.Suchadiskisexpectedtoformcollapsingpebblepilesat∼1−10auevenwith∼mm-cmsizedgrains!Middle:
Lowdensitydisk(MMSNafterlosing∼90%ofitsgasand∼50%ofitslargegrains).IntermediategrainsRight:High-densitydisk(∼10×MMSN),with
solarabundances.Onlylargegrainscanformpebblepiles.Althoughthemassandmeandensityincrease,andToomreQdecreases,withincreasingΣ0,the
parameterτs∝Σ−01(atfixedgrainsize)decreases.Sincethemaximumamplitudeofgraindensityfluctuationsscalessuper-exponentiallywithτs(whilethe
thresholdforcollapseisonlylinearinQ),thismeanssmallergrainscanpreferentiallyformpebblepilesinlowerdensitydisks(whereτs∼1).
4.2 TheMassFunctionofResultingPebble-Pile Lin2013):alargein-planevortexinadiskperturbingarazor-thin
Planetesimals (two-dimensional)dustdistributioninthemidplane;forwhichthe
same scalings apply as the three-dimensional case. It simply be-
Givenourassumptions,wecannowestimatethemassfunctionof
comesdustsurfacedensityfluctuationsthataredrivenbythelarge
collapsingdustdensityfluctuations.Fig.3showstheresultsforour
vortices trapping or expelling dust, rather than three-dimensional
“default”MMSNmodel(Σ =1,Z =Z ,α=10−4),atvarious
0 d (cid:12) densityfluctuations.Sosurfacedensityfluctuationscan,inprinci-
radii,assumingdifferentgrainsizes.
ple, form over a wide range of scales; for a large structure with
4.2.1 DependenceonOrbitalDistanceandGrainSize λe(cid:38)hd,theenclosedgrainmassintheperturbationbecomesM≈
πΣ λ2=2π(cid:104)ρ (cid:105) ρ˜ h3λ˜2.Basedontheargumentsin§3,we
Asexpected,ifthegrainsaresufficientlylarge(τs∼1),themodel expecrcittt (cid:38)t (sogβ0∼cr1it)odnthesescales,soρ˜ ∼(Qτ/2αλ˜)1/2.
grav s crit s
predicts that self-gravitating pebble piles will form over a range
Ifweassumeτ ∼1,andthatthelargestpossiblestructuresreach
s
of orbital radii, with a wide range of self-gravitating masses. For
∼H ,thenweobtain
g
R ∼10cm, all radii r ∼0.1−10 have τ ∼1 and form peb-
d au s
blepiles.Atstillsmallerradii,largeQvaluesimplysoundspeeds Mmax (τ ∼1)∼√2πα1/4Q1/2(cid:104)ρ (cid:105) H3 (40)
sufficienttosuppresscollapse;atlargerradii,τs(cid:29)1,andsograin- collapse s τs1/4 g 0 g
densityfluctuationsareactuallysuppressedbecausethegrainsare (cid:16) α (cid:17)1/4
approximately collisionless (large density fluctuations cannot be ∼0.03 10−4 Σ10/2ra2u7/28M⊕ (41)
generatedbythegas,forτ (cid:38)3−5).Forsmallergrains,wemust
s So the maximum mass is only weakly-dependent on α, while it
gotolargerradiibeforeτ ∼1,andcollapsebecomespossible.For
s increases with disk surface density and is nearly proportional to
R ∼1cm, pebble pile formation at (cid:28)10au is completely sup-
d radius(becausethediskmassincreaseswithr ).9Notethatthisis
pressed – we stress that because the density fluctuations depend au
notthe“typical”objectorstructure–wedonotexpectthedriving
exponentiallyonτ,thepredictednumberdensityis<10−10here!
s scaleforturbulence,forexample,toreachH .Ratheritiswhere
Weseethisrapidthresholdbehaviorsetinbetweenτ ∼0.1−0.3, g
s we expect a very sleep cutoff in the largest possible object sizes
aparameterspaceweexplorefurtherbelow.
(inFig.3,notethatthiscorrespondstoanumberofsuchobjects
4.2.2 TheMinimumandMaximumMasses afactorof∼107 lowerthanthe“typical”objects).Ifwewantto
estimateamoretypicalmass,where,say,themassfunctionbegins
Wherepossible,thesecollapseeventsformobjectswitharangeof
toturnovermoresteeply,weshouldtakeasizescaleof∼h or
masses∼10−8−10M⊕. ofthevelocitystructureswithturnovertimes∼Ω−1(whichshoduld
The maximum mass is given by the behavior of the largest havesizes∼α1/2H ).Forτ ∼1,thesearethesame,andthisgives
g s
velocity structures. Recall, in this model, grains are essentially
usacharacteristicmass
passive, so if structures of non-zero vorticity exist with λe (cid:38)hd √
(with the appropriate te ∼ts), they still drive grain density fluc- Mcinotlelarmpseediate(τs∼1)∼ 2παQ1/2(cid:104)ρg(cid:105)0Hg3 (42)
tuations in the midplane dust layer (so long as the eddy inter- (cid:16) α (cid:17)
∼3×10−4 Σ1/2r27/28M
sects the midplane somewhere) on scales ∼λe, even if we take 10−4 0 au ⊕
thedustlayertobeinfinitelythin.8 Indeed,thisisjustoneofthe Thelowermasslimitisalsopredictedbecauseonsufficiently
toy-model cases considered in Hopkins (2013c) (see also Lyra & small scales,t =λ /v (cid:28)t (so grains do not have time to
cross e drift e
8 Notethat,evenifthegrainsthemselvesdrivetheturbulence,asinthe 9 Interestingly,ifwehadverylarge-scalefluctuations,λ(cid:38)Hg,thenthe
streaminginstabilitycase,suchlargevelocitystructuresformviatheshear- shear/angularmomentumtermwouldbethedominanttermresistingcol-
ingofsmallerstructures,albeitwithrelativelylimitedpowerintheirasso- lapseandwewouldobtainMcmoallxapse∼πQ(cid:104)ρg(cid:105)0Hg3,i.e.justthestandard
ciatedvelocityfluctuations(seee.g.Johansen&Youdin2007). Hillmass.
(cid:13)c 0000RAS,MNRAS000,000–000