Shadows, gaps, and ring-like structures in protoplanetary disks

0.82 MB·English

by  Ralf Siebenmorgen

Astronomy & Astrophysics manuscript no. aa201118493˙final c ESO 2012
(cid:13)
January 18, 2012
Shadows, gaps, and ring-like structures in protoplanetary disks
R. Siebenmorgen1 and F.Heymann1,2
1 European Southern Observatory,Karl-Schwarzschild-Str. 2, D-85748 Garching b.Mu¨nchen, Germany
2 Department of Physics and Astronomy,Universityof Kentucky,Lexington, KY 40506-0055, USA
Received November21, 2011 / Accepted December30, 2011
ABSTRACT
2 Westudythestructureofpassively heateddisksaround T Tauriand HerbigAestars, and present avectorized Monte
1 Carlodustradiativetransfermodelofprotoplanetarydisks.Thevectorizationprovidesaspeedupfactorof∼100when
0 compared to ascalar version of thecode. Disksare composed of eitherfluffycarbon and silicate grains of varioussizes
2 or dust of the diffuse ISM. The IR emission and the midplane temperature derived by the MC method differ from
modelswheretheradiativetransferissolvedinslabgeometryofsmallringsegments.IntheMCtreatment,dustyhalos
n
abovethedisksareconsidered.HalosleadtoanenhancedIRemissionandwarmermidplanetemperaturethandopure
a
J disks. Under the assumption of hydrostatic equilibrium we find that the disk in the inner rim puffs up, followed by a
shadowedregion.Theshadowreducesthetemperatureofthemidplaneanddecreasestheheightoftheextinctionlayer
7
of the disk. It can be seen as a gap in the disk unless the surface is again exposed to direct stellar radiation. There
1
the disk puffs up a second time, a third time and so forth. Therefore several gaps and ring–like structures are present

in the disk surface and appear in emission images. They result from shadows in thedisks and are present without the
]
R need to postulate the existence of any companion or planet. As compared to Herbig Ae stars, such gaps and ring–like
structures are more pronounced in regions of terrestrial planets around T Tauri stars.
S
. Key words.radiative transfer – diffusion – dust, extinction – planetary systems: protoplanetary disks – infrared: stars
h
p
-
o 1. Introduction and gaps, spiral arms, and rings. Such disk structures are
r
inferredfromfittingthespectralenergydistribution(SED)
t
s Low–mass stars are surrounded by gaseous dust disks; by applying models considering pure disks (Calvet et al
a prominent examples are T Tauri and Herbig Ae/Be stars
2002; Espaillat et al. 2010; Mulders et al. 2010) or disks
[
(Waters & Waelkens, 1998). The existence of disks is de- with halos (Schegerer et al. 2008; Verhoeff et al. 2011).

1 rived from observations of the submillimeter continuum Recently, several disks with rings have been observed
v (Beckwithetal.,1990),theIRexcessoverthephotospheric
(van Boekel et al. 2005; Fukagawa et al. 2006; Ratzka et
7 emission (Kenyon & Hartmann, 1995), and direct imaging
al. 2007; Brown et al. 2007, 2008; Mayama et al. 2010;
7 (McCaughrean & O’dell, 1996). Protoplanetary disks with
Thalmann et al. 2010; Olofsson et al. 2011). In the gaps,
5 weak mid–IR and strong far–IR emission are called transi-
between the rings, companion candidates are detected
3
tionaldisks (Najita etal.2007;Furlanetal.2009;Ciezaet
(Huelamo et al. 2011; Kalas et al. 2008; Lagrange et al.
.
1 al. 2008).They are of interest because the missing mid–IR 2010).Hydro-modelsaccountfortheobservedradialdistri-
0 emission could be caused by grain growth or planet for-
butionofplanetsandlifetimesofthedisks,whichareafew
2 mation (Testi et al., 2003). Their optical extinction is flat-
Myr (Strom et al. 1989;Bertout et al. 2007).The strength
1 ter than the standard ISM extinction curve (Fitzpatrick &
of the disk–planet interactions are influenced by the gas
:
v Massa, 2007) and shows time variability (Schisano et al. pressuresothatthestructureofprotoplanetarydisksplays
i 2009). This indicates that large grains or clumpy material a key role in understanding the planet formation process.
X
orbits and partially obscures the star.
There is a rich literature on radiative transfer models
ar The evolution of the disk is followed theoretically by of protoplanetary disks (e.g. Wolf et al. 2003; Whitney et
studying the time dependence of the surface density of gas al. 2003; Pinte et al. 2006; Robitaille et al. 2006; Pinte
and dust. Hydro–dynamical models of viscous protoplane- et al. 2008). The structure of the disks is discussed by
tarydisksincludingplanetshavebeendeveloped(Lubowet Nomura (2002), Walker et al. (2004), and Dullemond &
al. 1999; Ozernoy et al. 2000; Klahr & Bodenheimer 2003; Dominik (2004a), and self–shadowing ripples in the disks
Edgar & Quillen 2008; Alexander & Armitage 2009). The havealreadybeenmentionedbyD’Alessioetal.(1999)and
modelspredictthatthediskisclearedbyanorbitingplanet Dullemond (2000); for a recent review, see Dullemond &
and that the surface density shows structures with holes Monnier(2010).Inthispaperthethermalstructureofdisks
isdiscussedbysolvingthehydrostaticandradiativebalance
Send offprint requests to: [email protected] equations. We are interested in deriving the structure and
1
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
IRappearanceofastarlight–heated,so–calledpassivedisk. results using dust opacities as derived for protoplanetary
Dust models used in this paper are described in Sect. 2, disks and the diffuse ISM. For the dust in the ISM, we use
and a vectorized three–dimensional Monte Carlo (MC) ra- a power–law size distribution: n(a) a−3.5 with particle
∝
diative transfer code in Sect. 3. We compare results of the radii between a− a a+ of bare spherical particles of
≤ ≤
MCmodelagainstamethodthatsamplesthediskinsmall silicates (a− = 320˚A, a+ = 2600˚A) with optical constants
ring segments at radius r from the star and in which the by Draine (2003) and carbon (a− = 160˚A, a+ = 1300˚A),
radiative transfer of the segments is solved in a slab ge- withopticalconstantsbyZubkoetal.(1996),andbothwith
ometry(Sect.4.3). Theinfluence ofdifferentdustopacities bulkdensityof2.5g/cm3.Weusedustabundances(ppm)of
on the SED of the disks is discussed in Sect. 4.4. We ap- 31 [Si]/[H] and 200 [C]/[H], which are in reasonable agree-
ply the models to T Tauri and Herbig Ae stars. The effect ment with cosmic abundance constraints (Asplund et al.,
of a dust halo on the emission and temperature structure 2009).This gives a gas–to–dustmass ratio of 130.For pro-
ispresented(Sect.4.5).Theverticaltemperaturestructure toplanetary dust we consider fluffy agglomeratesof silicate
introducesvariationsoftheheightofthedisksurfacewhich andcarbonassubparticles.Otherparametersareasforthe
causes shadows. They emerge as gaps and ring–like struc- ISM dust except that grains grow to 100 times larger radii
tures at distances from the star, which are thought to be ofa =33µm.Asrelativevolumefractionofthecomposite
+
the birthplaces of terrestrial planets (Sect. 5). Strikingly grain we use 34% silicates, 16% carbon, and 50% vacuum,
the structures are caused only by the coupling of the hy- whichtranslatestoarelativemassfractionof68%silicates
drostaticandradiativetransferequationswithouttheneed and 32% aC. Absorption and scattering cross-sections and
to postulate a planet. thescatteringasymmetryfactoriscomputedwithMiethe-
ory for the ISM grains, and we apply the Bruggeman rule
forthefluffy composites.Thisgivesatotalmassextinction
2. Dust models
cross section for the large grains in the optical (0.55µm)
of Kext = 30600 of the ISM dust and 4000[cm2/g–dust]
In protoplanetary disks the properties of dust are thought V
of the fluffy composites. The wavelength dependence Kext
to evolve, and the mineralogy, composition, porosity and λ
of both models is similar to those displayed in Kru¨gel &
sizedistributionofthegrainsdiffersfromthatofthediffuse
Siebenmorgen (1994, Fig.12).
ISM. Thanks to an enhanced particle collision rate, fluffy
grains are produced that grow in size over time (Natta et
al.2004;Ackeetal.2004;Natta etal.2007;Lommenetal.
2007; Ricci et al. 2010). Disk particles settle towards the
midplane under gravity, a process counteracted by grain 3. Vectorized MC radiative transfer
fragmentation and disk turbulence (Mizuno et al. 1988;
Weidenschilling et al. 1993; Sterzik et al. 1994; Dullemond TheradiativetransferproblemissolvedbyaMCtechnique
& Dominik 2004b, 2005). It is found that grains of radius byfollowingtheflightpathofmanyrandomphotons.Basic
>10µmsettletowardsthemidplane,anddependingontur- ideasaregivenbyWitt(1977),Lucy(1999),andBjorkman
∼bulentmotions,willbeblownuptothesurface,wherethey & Wood (2001), and the original version of our code was
remainforlongtimeperiods,andwellcoupledwiththegas developed by Kru¨gel (2008). Our model space is a three–
(Alexander & Armitage 2007; Charnoz et al. 2011). dimensionalCartesiangrid(x,y,z) that is partitioned into
Dust gets annealed in high–temperature regions of the cubes and, where a finer grid is required, further divided
disk,andcrystallinestructuresarebuilt.Theyaredetected into subcubes (cells). The star of luminosity L emits a to-
and the observed profiles of the 10µm silicate band differ tal of N = n Nν photon packages per second and in each
·
from that of the interstellar medium (Malfait et al. 1998; of the Nν frequency bins n photon packages are emitted.
Bouwmanetal.2008;Schegereretal.2006;Kessler–Silacci Each package has a constant energy ε = L/N. A package
et al. 2006;Furlan et al. 2011;Oliveira et al. 2010;Watson entering a cell on its flight path may be absorbed there or
etal.2009).About4%to7.5%ofthedustmassoftheSolar scattered. The probability of such an interaction is given
System is in crystalline silicates whereas crystallization is by the extinction optical depth along the path within the
2%inthediffuseISM(Kemperetal.,2005).Icecoatings cell, ∆τ, andoccurs when ∆τ log(ξ)using unified ran-
∼are built in frosty regions (Siebenmorgen & Gredel 1997; dom number ξ. The package is≥sc−atteredif the dust albedo
Pontoppidan et al. 2005; Visser et al. 2009). The location A>ξ′, using randomnumber ξ′; otherwise, it is absorbed.
wherewatericecanexistisimportantforterrestrialplanet Whenthepackageisscattered,itonlychangesdirectionde-
formation (Ida & Lin 2008; Min et al. 2011). All of these termined in a probabilistic manner by the phase function.
processes are combined by mixing and transport mecha- When it is absorbed, a new package of the same energy,
nisms and are valid for the neutral part of the disks where but usually different frequency ν′ is emitted from the spot
grain charging can be ignored. Ionization becomes impor- of absorption. The emission is isotropic. Each absorption
tant in the disk surface;therefore,the true situation of the event raises the energy of the cell by ε, and accordinglyits
mineralogyinprotoplanetarydisksisacomplicatedmatter temperature.
with a rich chemical network and evolution mechanisms in The computational speed of the code is increased con-
place (Gail, 2004). siderablybycalculatingtheflightpathsofphotonssimulta-
TheSEDandthethermalstructureofthedisksdepend neously.Thistypeofvectorizationisrealizedintwoflavors
on the applied dust model and for comparison we present usingeitherconventionalcomputerprocessingunits(CPU)
2
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
within FORTRAN 90 and OMP environment1 or graphics packagesthroughthemodelspaceisvectorizedinso–called
processing units (GPU) with NIVIDA cards and CUDA2 threads.Trajectoriesthathitthecellsofhighopticaldepth
(Heymann,2010).Thespeed–upscalesalmostlinearlywith have many more interactions than others, most likely the
thenumberofprocessingcoresavailable.Forparallelization majorityoftrajectories.Thisresultsinaratherunbalanced
particular attention needs to be paid to the random num- workloadover all threads, so that the advantage of vector-
ber generator, for which we choose the Mersene Twister ization is lost. In such cases the run time increases by the
algorithm (Matsumoto & Nishimura, 1998). number of CPU or GPU cores available, which in our case
Verification of the vectorized MC code against bench- are 8 and 480, respectively. This factor of efficiency loss is
marks (Ivezic et al. 1997;Pascucciet al. 2004)was done in ontopofthespeedupgainof5–10oftheMRWprocedure
anaccompanyingpaper(Heymann&Siebenmorgen,2012). achieved in scalar MC treatments (Pinte et al., 2009).
Thecodeis alsotestedagainstthesphericsymmetricalray
tracing code by Kru¨gel (2008). The SED is computed by
4. Protoplanetary disk models
counting the photon packets leaving the cloud towards the
observer using a large beam. The observed intensity can
Thetimescalefordiskstoattainthermalequilibriumisset
also be computed by following the radiative transfer of the
bythebalancebetweenheatingandcooling.Thistimescale
line of sight from the observer or a pixel of the detector
ismuchshorterthantheevolutionofthediskortheheating
plane, through the model cloud using the MC computed
source. The number density of the gas, even in the upper
dusttemperatures andscatteringevents(Heymann, 2010). layers of the disk, exceeds 105cm−3, so that gas and dust
Theraytracerisalsovectorizedandallowsustoderivedust
arethermallycoupled.Matterwillspiralinthedisktowards
scattering and emission images with high signal to noise,
the star and dissipate gravitational energy. For a thin disk
which is not possible by photon counting procedures.
that extends to the stellar surface, an accretion luminosity
The number of interactions of photon packets with the
of about a quarter of the stellar luminosity is converted in
dust increases exponentially with the optical depth of the
this way into heat. For classical T Tau stars, observed ac-
cpehlol.toInncpeallcskaogfeesxtarreemterahpigphedopatnicdalindteeprathct, τwVit>∼h1t0h0e0,dtuhset cretionratesare∼10−7...−9M⊙/yr(Muzerolleetal.,1998).
Atearlyepochswhentheaccretionisstrongest,dissipation
over and over again before they have a chance to escape.
dominates, but later heating by stellar radiation becomes
In this way MC treatments are slowed down considerably.
more important. As fiducial model of the T Tau star, we
Unfortunately, this situation appears in cells close to the
take a mass of 1M⊙, a luminosity of 2L⊙, and a photo-
midplaneofprotoplanetarydisks.Amodifiedrandomwalk
spheric temperature of 4000K. For comparison we treat a
(MRW) procedure for improving the computational effi-
disk heated by a Herbig Ae star with 2.5M⊙, 50L⊙ and
ciency of MC methods has been presented by Fleck & 104K (van den Ancker et al., 1997).
Canfield (1984). The MRW enlarges the mean free path
¿From near–IR interferometry (Millan-Gabet et al.,
length of packets by a diffusion approximation whenever
2006) of T Tauri and Herbig Ae stars, it has been found
necessary, and has been tested by Min et al. (2009) and
that the inner disk scales withstellar luminosity r √L.
Robitaille (2010). Let r be the distance of the interaction in ∝
Such a dependency is expected assuming that r is set
point to the nearest site wall of a cell. The MRW is con- in
by evaporation of grains at temperatures of 1000 1500K
sidered in a cell when the photon package has a distance −
(Dullemond & Monnier, 2010). For Herbig Be stars, larger
r γ/ρK , where γ > 5 is a user specified constant, ρ
≥ R deviations of the simple relation are measured at L >
the dust density, and K the Rosseland mean extinction
coefficient. The photon trRavel distance d is derived from 103L⊙. Spatially resolvedspectroscopy,at milliarcsec reso-
lution,allowsdetectionofahottergaseousemissioncompo-
1 r 2 nentclosertothestar(Eisneretal.2009;Najitaetal.2009;
d= ln δ (1) Pontoppidanetal.,2011).Evaporationtemperaturesofsil-
−D (cid:16)π(cid:17)
icateandcarbonparticlesarearound1500K.Thegrainsin
with diffusion constant D = 1/3ρK and a pre–computed the disk are assumed to be fluffy and to have been formed
R
sample of δ given by bycoagulationofmuchsmallercompactinterstellargrains.
The fluffy agglomerates are held together by weak van der
ξ =2 ∞ ( 1)i+1 δi 2 (2) Waals forces (binding energies <0.1eV), and they are ex-
− pected to fall apartinto their re∼fractoryconstituents when
Xi=1 (cid:0) (cid:1) thetemperatureofthecompositegrainexceeds1000K.The
constituent particles will then be hotter than the average
with unified random number ξ. The number of absorption
porous grains because they are much smaller. We assume
eventswithinacellisusedtocomputethedusttemperature
thatthediskextendsinwardsuptothepointwhereporous
and, within r of the MRW, is estimated by n =d ρ K ,
abs P
grainsreach1000K,orequivalentlysmallinterstellargrains
where K is the Planck mean mass extinction coefficient
P
would be at about 1500K.
of the dust. Ignorance of these acceleration methods has a
Instar–formingregionsdustisheatedtoabout30K.At
tremendouseffectontheruntimerequirementsofavector-
a large distance from the star grain heating by the ambi-
izedMCtreatment.Inourschemethetrajectoryofphoton
ent radiation field of nearby low or massive stars becomes
1 www.openmp.org important. We do not make additional assumptions on the
2 www.nivida.com outer radiation field and consider only the primary central
3
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
Table 1. Parameters of the fiducial T Tauri and Herbig Ae disks.
Parameter T Tauri Herbig Ae
Stellar luminosity L∗ [L⊙] 2 50
Stellar mass M∗ [M⊙] 1 2.5
Photospheric temperature T∗ [K] 4,000 10,000
Column density Σ(r)= τ⊥(1AU)( r )γ [g–dust/cm2] r<1AU: γ =0.5
KV AU
r≥1AU: γ =−1
Vertical optical depth τ⊥(1AU) 10,000
Dust density in halo ρ [g-dust/cm3] 0 or 1.5×10−18
halo
Innerdisk radius r evaporation
in
Outerdisk radius r [AU] 22.5 40
out
heating source. Models are computed to an outer radius of The disk is divided into small ring segments of width
r = 22.5AU for T Tau and 40AU for Herbig Ae disks ∆r at radius r from the star. For each ring the radiative
out
where the midplane temperature drops below 30K. transfer is solved under all inclination angles θ, for incom-
Initially weassumethatthediskisisothermalinz,and ing I− and outgoing I+ radiation, assuming a slab geome-
the density is given by try. In the vertical direction the opacities are so high, e.g.
τ⊥(1AU)=104,that eachring segmentis split into a com-
ρ(r,z)= 2 Σ(r) e−z2/2H2 (3) pletelyopaquemidlayerandamuchthinnertoplayerofop-
rπ H(r) tical thickness τtop. The midlayeris assumedto be isother-
mal at temperature T and the top layer extends so far
mid
with scale height H2 = kTmidr3/GM∗m, surface density downthat the temperature atits bottom approachesTmid.
Σ(r), molecular mass m = 2.3m and midplane tempera- For computing the spectral energy distribution it is suffi-
p
ture, Tmid, for which we use a power law as initial guess. cient to choose τtop ∼ 20. The ray tracer is an iteration
The height of the disk is set to z = 4.5 H. In models procedureandyieldsthetemperaturestructureT(r,z)and
with pure disks, the density ρ 0above z > z is 0 and theintensitiesoftheupwardI+ anddownwardI− directed
halo o
constant when a halo is considered (Table 1). The surface radiation. By observing the disk at viewing angle θobs, the
density is adjusted to reach a given optical depth in the received flux is computed from the emission of the star,
verticaldirectionτ⊥ =Σ(r)KV,withvisualextinctionKV. the intensity I+ at the disk surface at θ = θobs as well as
We use τ⊥ =τ1AU (r/AU)γ with γ = 1 for r 1AU and by summing up contributions from all rings. The radiative
γ = 0.5 otherwise (Min et al. 2011). F−or the p≥orous grain transfer is solved for each ring segment at radius r sepa-
model with a = 33µm, we take a vertical optical depth rately. However, the propagation of the radiation in radial
+
from the surface to the midplane at 1AU of τ⊥ = 10,000. direction from one ring into the next is ignored. The pro-
This choice holds computational time requirements of the cedure is often called the 1+1D disk model.
MCcodewithinreasonablelimits.Ittranslatestoasurface
density ofΣ(1AU)=5g-dust/cm2 oragassurfacedensity,
4.2. Disks with axial symmetry
which is half the value estimated for the minimum mass
of the early solar nebulae (Hayashi, 1981). Higher surface In disks that are in hydrostatic equilibrium in a vertical
densities can be obtained considering larger or ice–coated direction,thegravitationalforceisbalancedbythepressure
grains because K is reduced in such models (Kru¨gel & gradient3:
V
Siebenmorgen, 1994).
zGM∗ 1dP
= (4)
− r r2 ρ dz
4.1. Disks with slab geometry
with pressure P = ρkT(z)/m. For each height z and ra-
A solution of the radiative transfer of a hydrostatic and
dius r2 = x2 +y2 from the star, the temperature T(r,z)
geometrically thin disk is developed by Kru¨gel (2008,
is computed using the MC code of Sect. 3. The star is set
Sect.11.3).The disk is symmetric withrespectto the mid-
atorigin.Azimuthaldependencies inthe diskstructureare
planeatz =0.Thedensitystructureisgivenincylindrical
not considered, so that the problem can be solved in axial
coordinates ρ(r,z). Light from the star falls on the disk
symmetry.Itallowsustosolvetheradiativetransferinone
under a grazing angle α so that the star is visible from
gr octant of a cube and reduces the computational time by a
everywhereonthedisksurface(Armitage,2007).Flatdisks
factor 8 when compared to a full three dimensional treat-
with a constant grazing angle of α =2o and flared disks,
gr ment. For each height z 0 and radial distance at x 0
similar to Chiang & Goldreich (1997), are considered with ≥ ≥
and y 0, we compute the azimuthal average of absorbed
3o αgr r2/7 7o. The flat disk has a smaller grazing ≥
≤ ∝ ≤
angle and intercepts less stellar radiation than the flared 3 The protoplanetary disk masses are well below the stellar
disk, and even more so at large distances from the star. mass Mdisk≤0.01M∗.
4
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
photon packets. This average is used to derive the dust Theextinctionlayeris sampledwithcubes havingτ <
V
temperature T(r,z), and with Eq. 4 a new estimate of the 0.6. This is achieved by sampling cubes located close t∼o
density structure of the disk ρ(r,z). Obviously there is an the extinction layer in up to 20 subcubes. In total the disk
iterative scheme: starting with an initial guess of the mid- has 107 cells, which provides a sampling with high spatial
plane temperature and a first set up of the disk structure resolution.IneachMCrunanumberof108photonpackets
ρ(r,z) using Eq. 3, we compute the temperatures T(r,z) areemittedfromthestar.Thefairlyhighnumberofpackets
after a MC run. Once the temperatures are inserted into provideahigh–energyresolutionandisrequiredtoarriveat
Eq. 4, we arrive at a new density structure ρ(r,z), which a stable estimate of the temperatures T(r,z), in particular
is used to update T(r,z) with a new MC run. Typically near the midplane.
after less than a dozen iterations, the program converges
to a stable disk configuration. The problem we are solving
4.3. Comparison: 1+1D versus 2D
is intrinsically two dimensional, but our MC code is based
oncubiccellsthatarethreedimensional,thereforewelabel
TheSEDandthe midplane temperaturesofaflatandflar-
such computations 2D disk models hereafter.
ing1+1Ddiskarecomparedwiththoseofa2Ddiskmodel.
Particular attention has to be paid to a good set–up As an example, we take a T Tauri star with parameters as
of the grid used in the MC code. In the positive (x,y,z) specifiedinTable1.Foreaseofcomparison,thedistribution
directions, we use a basic grid of (350,350,50) cubes. At ofthe dustcolumndensityΣ(r) ischosentobe identicalin
locations where the radiation field varies strongly, cubes allthreemodels.Wetaketheonecomputedbythe2Ddisk
aresampledintosubcubes.Largegradientsoftheradiation as evaporation radius. The SED is compared in the upper
fieldareexpectedclosetothesurface,withintheextinction panel of Fig. 1. In a 1+1Ddisk the heating is proportional
layerorattheinnerwallofthedisk.Estimatesoftheheight to the grazingangle. Since flaring disks have a largergraz-
oftheextinctionlayerz0 anditsthicknessℓ H/2isgiven ing angle than flat disks they show higher dust re-emission
∼
by Siebenmorgen & Kru¨gel (2010). luminosities and appear warmer.The assumption of 1+1D
disk models is that stellar photons impinge upon the disk
at a given radius with constant grazing angle. Such a sim-
plepicturefailsclosetothedustevaporationzone.Inearly
1+1D models the inner disk region was treated as an op-
2 tically thick vertical wall at constant temperature (Natta
3 et al. 2001; Dullemond et al. 2001). Then the dependency
y) 10 1 of gas pressure on the dust evaporation temperature is in-
x (J cluded, resulting in a curved evaporation zone (Isella &
u Natta 2005; Kama et al. 2009). In the MC scheme the ra-
fl 1D
diationtransportin radialdirectionis included, which,be-
cause of the fairly difficult geometry, is otherwise not easy
to implement in ray–tracing techniques as used in 1+1D
1
disks. The SED of a 2D disk is very distinct with a much
10 100
wavelength (µm) stronger mid–IR emission and overall flatter appearance.
We also show in Fig. 1 the convergence of the SED and
the midplane temperature of the MC model as computed
during the first nine iterations of the density structure fol-
1000 2D 12 lowingEq.4.Themidplanetemperaturedistributionsofthe
3
4 1+1D disks are described well by a power law, whereas in
5
K] 67 the 2D diskstriking upanddownswith radiusarederived.
T [mid 1D 89 Wtemeaetxitceenffsievcetlsyttheasttecdoutlhdabtethinetrsotrduuccteudrebyisthneotndoiuseectoaussyesd-
100 by a pure photon statistics or inadequate sampling of the
MC cells. Following Eq. 3 the midplane temperature is di-
rectly linked to the density, and one expects a hilly path,
as derived in the lower panel of Fig. 1, to be reflected in
1 10
radius (AU) the overallobserved disk structure (Sect. 5).
4.4. SED dependency on dust opacities
Fig.1. Top: IR emission of a T Tau disk viewed at 45◦,
∼
whichiseithercomputedusingthe1+1Dmodel(label’1D’) A realisticdescriptionofdust inprotoplanetarydisks is an
in a flat (green full line), a flaring (dashed green line) disk open and controversial issue, as are the applied dust ab-
configuration, or the 2D disk. Iterations 1 to 9 of Eq. (4) sorption and scattering efficiencies. Because these parame-
of the 2D disk are shown. Bottom: Midplane temperatures ters directly enter into the radiative transfer equation and
of the 1+1D flat disk and the various iterations of the 2D theMCmethod,wewishtostudytheirimpactontheSED.
disk. As anexample we take a T Tauristarwith the parameters
5
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
25
fluffy a=33 µm
fluffy a++=0.3µm 10 T Tauri disk + halo
ISM
20
)
y
J y)
x ( 15 (J
u x
fl 10 flu
1 pure disk
5
10 100
wavelength (µm)
Herbig
Fig.2. IR emission of 2D disks viewed nearly face–on ( 100
20◦) computed by applying different dust opacities. ∼
in Table 1. Three different dust models are distinguished:
a) ISM like grains,b) fluffy grainswith particle radiiup to 10
a = 0.3µm, and c) fluffy grains with a = 33µm. Other
+ +
parameters remain as given in Sect. 2. The SEDs of the
disks are compared using identical distributions of the col-
umn density Σ(r). Since the dust models have different vi-
1 10 100
sual extinction coefficients KV, the vertical optical depth wavelength (µm)
at 1AU of 104 for our standard dust (model c) needs to be
increased by a factor 7.7 for ISM dust and by a factor 8.5
for fluffy grains with size distribution, as for ISM grains. Fig.3. IR emission of 2D disks viewed nearly edge–on (
SEDs of the 2D disks are displayed in Fig. 2, where one 20◦)with(dashed)andwithout(fullline)adustyhalo.Th∼e
notices that: disks are either heated by a T Tauristar (top) or a Herbig
i) The flux in the far–IR peaks at longer wavelength Ae star (bottom).
when applying fluffy grains.
ii) The silicate emissionprofile depends ongrainporos-
ity, particles size, and dust temperature (Voshchinikov &
above the disk (Vinkovic et al., 2006). Dust in a halo de-
Henning, 2008). In the band the cross section peaks at
cays through radiation pressure and gravity and may be
9.5µmforISMdustandat9.8µmforthe fluffycomposites.
replenished by dynamical processes, such as outflows and
Thefeature–to–continuumratiooftheextinctioncrosssec-
winds from the disk surface, or by planets in highly ec-
tion in the 10µm band, when averaged over the dust size
centric orbits (Krijt & Dominik, 2011). Miroshnichenko et
distribution, is 25% higher for ISM than for fluffy grains
al. (1999) added the emission of an optical thin halo onto
of same size. Nevertheless, because of the detailed temper-
the SEDs computed for optically thick disks. They found
ature structure, model b), which is made of fluffy grains
that the halo dominates the IR emission and provides ad-
with sizes typical of the ISM, shows the most striking sili-
ditional disk heating on top of the direct stellar radiation.
cate emission feature.
The MC scheme allows a self–consistent treatment of the
iii) The disk with fluffy grains up to a = 0.3µm have
+ radiative transfer in a configuration of a disk plus halo. In
a similar submillimeter slope as the one with 100 times
theTTauridiskmodels(Table1)weaddahalosothatthe
larger, a = 33µm, particles (Fig. 2). This agrees with
+ minimum dust density of Eqs. 3 and 4 in each cell of the
observations of protoplanetarydisks where the silicate fea- MCgridexceedsρ(r,z) ρ =1.5 10−18(g-dust/cm3).
halo
tureisnotcorrelatedwiththemillimeterslopeoftheSEDs ≥ ×
This gives an optical depth for the halo of τ 0.4. The
V
(Lommen et al., 2010). ∼
resulting SED of a disk plus halo is shown in Fig. 3. The
SED appears warmer than the pure disk with ρ = 0.
halo
In the example, the disk plus halo model produces a fac-
4.5. Disks plus halos
tor 3 stronger IR peak emission and less near–IR flux
∼
Protoplanetary disks are thought to be formed from their for wavelengths below 6µm. The halo provides additional
initialspheroidalconfigurationoftheprotostellarenvelopes heating on top of the direct stellar light. Their midplane
(Watson et al., 2007). It is speculated that there is some temperature becomes a factor 2 warmer than derived in
∼
residual dust left or that it is replenished at high altitudes a disk without halo within 1 – 10AU of the star.
6
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
WerepeattheexerciseforthemoreluminousHerbigAe
stars with the parameters given in Table 1. With the same
0.3
minimum dust density, ρ , the optical depth of the halo TTS disk+halo
halo
is τV 0.6 as for the T Tauri star. The resulting SEDs
∼ gaps
are shown in Fig. 3. The halo dominates the IR emission
for wavelengths >3µm and produces warmer dust in the r 0.2
mthiedpHlaernbeigthAane sint∼aprsuraeredidskissp.lTayheedminidFpilga.ne6.teTmhpeemraidtuprleanoef z / bot pure disk
temperature in the pure disk shows a rapid decline in the
inner3–4AU,followedbya1/r0.6 dependencyupto12AU 0.1
1d
and a strong decline farther out. The disk plus halo model gaps
is smoother and fit by T 1/r in the 2–15AU region
mid
∼
(Fig. 6). 0.5 1 10
radius (AU)
5. Ring structure of protoplanetary disks
0.30
We notice that the midplane temperature in 2D disks of
T Tauri stars show a wavy structure in particular near the 0.25 0 25 50 75 100
zonewhereterrestrialplanetsaresupposedtoform(Fig.3).
Hfesetrseiwtseelifnivnetshtiegaatpepheaorwantcheisoftetmhepderisakt.uFreirbstehwaevsioturdmyatnhie- z / rbot0.2

height of the photosphere which we define as the height of 0.15
the bottom of the extinction layer, z , where the vertical
bot
optical depth τ⊥(V) = 1. In Fig.4 we show the height of 10
thediskphotosphereasafunctionofradius,expressedasa
unitlessratiozbot/r.WedothisfortheTTauridiskmodels y (AU)1
with and without a halo using parameters of Table 1. For 1 x (AU) 10
1+1D disks the midplane temperature is a smooth func-
tion declining with radius, which can be approximated by 0.30
a power law, and so is the scale height which smoothly in-
0.25 0 25 50 75 100
creases with distance (Fig.4). The geometrical thickness of
tthheenT1h+ite1s2DtDhdidciskiksnkeinnssceraderaetscheliesnewevsiatshpiomrraialdtaiiruonst.ozotnheeisdpecurffeeadseupin,atnhde z / rbot0.2
0.15
midplane temperature.The temperature reductionis to be
understand by shadowed region of the disk surface where

light from the star is extincted and where grain heating 10
becomes less efficient. The shadow is located close behind
y (AU)
from the puffed–up inner rim. The disk becomes thicker 1
with increasing distance, because the gravity decreases, so 1 x (AU) 10
that there is a point where the shadow becomes less effi-
cient, and the disk is exposed to direct stellar light. There
Fig.4. Top: The relative height of the disk photosphere,
again the disk is puffed up followed by a second shadow,
z /r,asafunctionofradiusofadiskheatedbyaTTauri
which is located at about 1AU in our example and a third bot
one near 4AU. Farther out (>10 AU) the midplane tem- star without (full line) and with (dashed) halo. The scale
height over radius, H/r, of a flat 1+1D disk (dashed) is
peraturedropstobelow 30K∼becausethe stellarheatingof
shownforcomparison.Gapsbetweenringstructuresarein-
the passive disk is inefficient, and another distortionin the
dicated.Ashadedrepresentationofz /r alongthe (x,y)–
disksurfaceisnotvisibleinmostcases.Thedetailedstruc- bot
planeisdisplayedforadiskplushaloconfiguration(middle)
ture of passively heated T Tauri disks are shown in Fig.4
anda pure disk (bottom). Colorbarsgive the relativegray
within the regionofterrestrialplanets.Ifone considersha-
scale (%) in units of z /r.
los,then gapsandring–likestructuresaredimmedandthe bot
disk appears thicker because their midplane temperature
is warmer. The disk surface is visualized in a three di-
there is a wide gap opening farther out and smaller gaps
mensionalrepresentationbyrotatingthefunctionz (r)/r
bot visible closer to the star. This structure is dimmed when a
alongthe midplane (Fig.4). The innermostpuffed–up rim,
halo is considered.
several gaps, and an overall wavy surface structure of the
disk are visible. In Fig. 5 we display images and surface
profiles for the mid–IR. At other wavelengths the emission
structure is similar. In the image of the disk without halo
7
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
4.7 4.0 3.3 2.5 1.8 1.1 0.4 1000 Herbig Ae
+225
+135 r-0.6
]
K
[
+45 mid
mas) T 100 r-1
α ( -45

-135
0.3
disk+halo
-225
-225 -135 -45 +45 +135 +225 gap
δ (mas)
5.4 4.3 3.2 2.2 1.1 0.0 -1.0 z / r bot 0.2 pure disk

gap
+225
0.1
+135

+45 3 10
as) radius (AU)
m
α ( -45
0.30
-135
0 25 50 75 100
0.25
-225
-225 -135 -45 +45 +135 +225
2 x (AδU (m) as)6 8 z / rbot0.2

1000 0.15
disk+halo
2y/arcsec) 100 pure disk
ness (J 20
ht
brig 10 y (AU)
3
gaps 3 x (AU) 20
00 45 90 135 180
δ (mas) Fig.6.Top:Themidplanetemperatureandrelativeheight
ofthe disk photosphere,z /r, as afunction ofradius ofa
bot
Fig.5.Mid–IRimagesofdisksheatedbyaTTauristarat
2Ddiskheatedby aHerbigAe star.Models areshownina
distance of 50pc and viewed at 30◦. The image of a pure
configuration of a pure disk (full line) and a disk plus halo
disk (top) and a disk plus halo is shown (middle). Color
(dashed). Power law fits (dotted) to the midplane temper-
barsareinlog(Jy/arcsec2).Bottom:Correspondingsurface
ature distribution are as labeled. Gaps between ring struc-
brightness distribution measured as cut along δ through tures are indicated. Bottom: A shaded representation of
origin of the image. Gaps between emission rings are indi-
z /ralongthe(x,y)–planeisshownforapurediskmodel.
bot
cated.
Color bar gives relative gray scale (%) in units of z /r.
bot
8
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
6. Conclusion
5.6 4.8 4.1 3.3 2.5 1.8 1.0 The infrared appearance of passive disks around low and
+200 medium mass stars has been discussed. Dust particles are
made either ofbare materialorfluffy composites ofsilicate
+120 and carbon. Grain sizes range from 160˚A to 3000˚A as de-
rived for the ISM and up to 33µm assuming grain growth
inthedisks.Thedisksareheatedbystellarlightwithlumi-
+40
as) nosities between 2L⊙ and 50L⊙. The disks have an inner
m
α ( -40 radius that is roughly determined by the dust evaporation
temperature of the grain material. The surface density of
the disks was set close to what is estimated for the mini-
-120 mum mass of the solar nebula. Disks are configured with
or without halos. Two hydrostatic and geometrically thin
-200 diskscenarioswereexaminedforwhichweassumethatgas
-200 -120 -40 +40 +120 +200 and dust are mixed and have the same temperature.
δ (mas)
i)In1+1Dthedisksweredividedinsmallringsegments
in which the radiation transfer was solved. For each ring,
2 4 x (AU ) 8 10 a slab geometry was applied, and it was assumed that in
105
Herbig Ae deeper layers the disk is isothermal. The transport of radi-
ation in the radial direction was ignored. From the surface
2c) 104 of 1+1D disks, the star is always visible. A flat or flaring
arcse disk+halo structureofthe disk wasconsidered.This type ofdiskis in
ess (Jy/ 103 pure disk widieis)pIrnead2Dusdeisinksththeelitheyrdartuosrtea.tic and radiative transfer
n
bright 102 eaqnuaatcicounrsatweerMeCsolmveedthinodanwiatesraptrievseensctehdemteh.aFtocratnhebleatatper-
plied to arbitrary dust geometries. We used an adaptive
101
three–dimensional Cartesian grid where cubes are divided
0 50 100 150 200
δ (mas) into subcubes whenever required, for example, close to the
disk surface. The MC method is vectorized and the paral-
lelization is realized to work on graphics cards or conven-
Fig.7. Top: Mid–IR image of a disk heated by a Herbig tional processing units and provides a speed–up roughly
Ae star at distance of 50pc and viewed at 30◦. Color bar proportional to the number of processor cores or parallel
is in log(Jy/arcsec2). Bottom: Mid–IR surface brightness threads, available. On our conventionalcomputer, the vec-
distribution measured as a cut along δ through the origin torizedcodeisafactor100fasterthanthescalarversion.In
of an emission image of a pure disk (full line) and that regionsof veryhigh opticaldepth, for example close to the
computed in a disk plus halo configuration (dashed). midplane of the disk, photon packets may get trapped. In
this casethe algorithmisfurther acceleratedbya modified
randomwalkprocedure.The 2Ddisksprovideamorereal-
isticdescriptionthanthe1+1Ddisks:IntheRTsolutionof
the1+1Ddisksitisassumedthattheradialtransportofra-
diationcanbe ignoredandthatlayerswith verticaloptical
ofτ >20areisothermal.Theseassumptionsbreakforex-
top
Disks of the more luminous Herbig Ae stars have an ample∼neartheinnerregionsofthedisk.Ontheotherhand,
evaporation zone farther out than around T Tauri stars. inthe2Ddiskstheradialtransportofradiationisincluded
At these distances the gravitationis reduced and the disks andnoassumptionsaremadethatregionsofthedisksneed
become thicker than for T Tauri disks. The height of the to be isothermal. Therefore a self consistent treatment of
Herbig Ae disks and a shadowed presentation is displayed the hydrostatic balance is only provided in the 2D disks.
inFig.6.WhencomparedtoTTauristars(Fig.4),the ra- Our main findings are that:
tio z (r)/r is indeed higher and has less structure. There
bot
is a pronouncedpuffed–up inner rimcausinga far reaching
shadow, so that a second one farther away than 30AU is – 1+1D disks have a smooth surface structure, and their
only envisaged. However, this far away the midplane tem- midplane temperature is approximated by a power law
perature has dropped to below 30K (Fig. 3) and other ef- well. They gravely underestimate the mid IR emission
fects than stellarheating maydominate the disk structure. when compared to 2D disks.
Wedisplaythemid–IRimageandsurfaceprofileofthedisk – Halos substantially increase the IR emission and the
heated by a Herbig Ae star in Fig. 7. At least one striking temperature of the midplane, which is smoother than
gap and a ring in the disk are identified. obtained in disks without halos.
9
R. Siebenmorgen and F.Heymann: Shadows, gaps, and ring-like structuresin protoplanetary disks
– In 2D disks the midplane temperature shows up and Espaillat C., D’Alessio P.D., Hernandez J., et al., 2010, ApJ, 717,
downswith radius(Fig.1).Sucha hilly structureis also 441
revealed in their surface structure. FitzpatrickE.L.,MassaD.,2007,ApJ,663,320
FleckJr.J.A.,CanfieldE.H.,1984,J.Comput.Phys.,54,508
– Emission images of 2D disks display gaps and ring–like
Fukagawa M.,TamuraM.,ItohY.,etal.,2006,ApJ,636,L153
structures inparticularinthe regionofterrestrialplan- FurlanE.,WatsonD.M.,McClureM.K.,2009,ApJ,703,1964
ets. Disk structures are caused by puffing up the disk FurlanE.,LuhmanK.L.,EspaillatC.,2011, ApJS,195,3
surface and followed up by their shadows behind. The GailH.P.,2004,A&A,413,571
HayashiC.,1981,PThPS,70,35
disk structure is dimmed when halos are considered.
HeymannF.,2010,PhD,UniversityofBochum.
Rings and gaps are more pronounced in T Tauri stars
HeymannF.,SiebenmorgenR.,2012, ApJ,submitted
than in Herbig Ae stars. HuelamoN.,LacourS.,Tuthill,P.,etal.,2011,A&A,528,L7
IdaS.,LinD.N.C.,2008,ApJ,685,584
IsellaA,NattaA.,2005,A&A,438,899
The detected gap and ring structures of the disks are
IvezicZ.,GroenewegenM.A.T.,MenshchikovA.,SzczerbaR.,1997,
only caused by the assumption of hydrostatic equilibrium MNRAS,291,121
andthederiveddetailedverticaltemperatureprofilesofthe KalasP.,GrahamJ.R.,ChiangE.,etal.,2008,Science, 322,1345
disk. On the other hand, it is easy to imagine, without the KamaM.,MinM.,DominikC.,2009,A&A,506,1199
KemperF.,VriendW.J.,TielensA.G.G.M.,2005, ApJ,633,534
need of such detailed computations, that behind a puffed–
KenyonS.J.,HartmannL.,1995,ApJS,101,117
up disk the heating is dimmed and therefore that the dust
Kessler–SilacciJ.,AugereauJ.-C.,DullemondC.,etal.,2006,ApJ,
temperaturesarelower.Thiscausesthatinthoseshadowed 639,275
regions the height of the disk is smaller unless it is again KlahrH.K.,Bodenheimer P.,2003,ApJ,582,869
exposed to direct stellar light. Therefore the detected gap KrijtS.,DominikC.,2011,A&A,531,A80
Kru¨gelE.,SiebenmorgenR.,1994,A&A,288,929
and ring–like structures are a plausible effect of the shad-
Kru¨gelE.,2008AnintroductiontothePhysicsofInterstellarDust,
ows, but they have not yet been reported by others. We IoP
feel that this is important because ring–like structures of Lagrange A.-M., Bonnefoy M., Chauvin G., et al., 2010, Science,
protoplanetarydisksareofteninterpretedas“fingerprints” 329,57
Lommen D., Wright C.M., Maddison S.T., et al., 2007, A&A, 462,
oftheplanetformationprocess,whereasinthemodelspre-
211
sented planets have not been considered.
LommenD.,vanDishoeckE.F.,WrightC.M.,MaddisonS.T.,etal.,
2010, A&A,515,77
Acknowledgements. WearegratefultoEndrikKru¨gelfordiscussions Lubow S.H,SeibertM.,ArtymowiczP.,1999,ApJ,526,1001
andhelpfulsuggestions. Lucy, L.B.,1999,A&A,344,282
MalfaitK.,Waelkens C.,WatersL.B.F.M.,1998, A&A,332,25
Matsumoto M., Nishimura T., 1998, in: Monte Carlo and Quasi-
MonteCarloMethods,Springer
References
MayamaS.,TamuraM.,HanawaT.,etal.,2010,Science, 327,306
McCaughreanM.J.,O’dellC.R.,1996, AJ,111,1977
Acke B., van den Ancker M.E.,C. P. Dullemond, 2004, A&A, 422,
Millan–GabetR.,MonnierJ.D.,BergerJ.P.,etal.,2006, ApJ,645,
621
L77
Alexander R.D.,ArmitageP.J.,2007,MNRAS,375,500
MinM.,DullemondC.P.,deKoterA.,etal.,2009, A&A,497,155
Alexander R.D.,ArmitageP.J.,2009,ApJ,704,989
MinM.,DullemondC.P.,KamaM.,DominikC.,2011,Icarus,212,
Armitage P.J., 2007, Astrophysics of planet formation, University
416497,155
press,Cambridge
Miroshnichenko A., Ivezic Z., Vinkovic D., Elitzur M., 1999, ApJ,
Asplund, M.,Grevesse, N.,Sauval, A.J.,Scott, P.,2009, ARA&A,
520,115
47,481
MizunoH.,MarkiewiczW.J.,VoelkH.J.,1988,A&A,195,183
BeckwithS.V.W.,SargentA.I.,ChiniR.S.,Gu¨stenR.,1990,AJ,99,
MuldersG.D.,DominikC.,MinM.,2010,A&A,512,A11
924
MuzerolleJ.,HartmannL.,CalvetN.,1998,AJ,116,2965
BertoutC.,SiessL.,CabritS.,2007, A&A,473,L21
NajitaJ.R.,StromS.E.,MuzerolleJ.,2007,MNRAS,378,369
Bjorkman,J.E.,Wood,K.,2001,ApJ,554,615
NajitaJ.R.,DoppmannG.W.,CarrJ.S.,etal.,2009,ApJ,691,738
BouwmanJ.,HenningTh.,HillenbrandL.A.,2008, ApJ,683,479
NattaA.,PrustiT.,NeriR.,etal.,2001,A&A,371,186
Brown J.M., Blake G.A., Dullemond C.P., et al., 2007, ApJ, 664,
NattaA.,TestiL.,NeriR.,etal.,2004,A&A,416,179
107
NattaA.,TestiL.,CalvetN.,etal.,2007,in:ProtostarsandPlanets
BrownJ.M.,BlakeG.A.,Qi,C.,etal.,2008,ApJ,675,1009
V,B.Reipurth,D.Jewitt,andK.Keil(eds.),UniversityofArizona
CalvetN.,D’Alessio,P.,HartmannL.,2002, ApJ,568,1008
Press,Tucson,951,767
CharnozS.FouchetL.,AleonJ.,MoreiraM.,2011,arXiv1105.3440C
NomuraH.,2002ApJ,567,587
ChiangE.I.,GoldreichP.,1997,ApJ,490,368
OliveiraI.,Pontoppidan K.M.,MernB.,etal.,2010,ApJ,714,778
CiezaL.A.,SwiftJ.J.,MathewsG.S.,WilliamsJ.P.,2008,ApJ,686,
OlofssonJ.,BenistyM.,AugereauJ.-C,etal.,2011,A&A,528,6
L115
Ozernoy L.M.,Gorkavyi N.N., Mather J.C., Taidakova T.A., 2000,
D’AlessioP.,CantoJ.,HartmannL.,etal.,1999, ApJ,511,896
ApJ,537,L151
DraineB.T.,2003,ApJ,598,1026
Pascucci I.,WolfS.,Steinacker J.,etal.2004,A&A,417,793
DraineB.T.,LiA.,2007,ApJ,657,810
PinteC.,MenardF.,DucheneG.,BastienP.,2006,A&A,459,797
DullemondC.P.,2000,A&A,361,L17
PinteC.,Padgett D.L.,MenardF.,etal.,2008,A&A,489,633
DullemondC.P.,DominikC.,NattaA.,2001,ApJ,560,957
PinteC.,HarriesT.J.,MinM.,etal.2009,A&A,498,967
DullemondC.P.,DominikC.,2004a, A&A,417,159
PontoppidanK.M.,DullemondC.P.,vanDishoeckE.F.,etal.,2005,
DullemondC.P.,DominikC.,2004b, A&A,421,1075
ApJ,622,463
DullemondC.P.,DominikC.,2005, A&A,434,971
Pontoppidan K.M.,BlakeG.A.,SmetteA.,2011,ApJ,733,84
DullemondC.P.,MonnierJ.D.,2010, ARA&A,48,205
Ratzka Th.,LeinertC.,HenningT.,etal.,2007,A&A,471,173
EdgarR.C.,QuillenA.C.,2008,MNRAS,387,387
RicciL.,TestiL.,NattaA.,etal.,2010,A&A,512,15
Eisner J.A., Graham J.R., Akeson R.L., Najita J., 2009, ApJ, 692,
309
10