Star Formation Newsletter#368 29-42(KASI)

Downloaded from https://academic.oup.com/mnras/article/52

nras/article/525/3/3998/7246898 by Dankook University / Library user on 10 October 2023

E. Hoemann, S. Heigl and A. Burkert.oup.com/mnras/article/525/3/3998/7246898 by Dankook University / Library user on 10 October 2023

4000

Downloaded from https://academic.oup.com/mnras/article/525/3/3998/7246898 by Dankook University / Library user on 10 Oct

f which|time-scalea √| 0.42+0theµ = km.s−1 )for(3)
et al. 2019; Bhadar
28Aends
and cs being the sound
speed
(0.19
mean molecular
Heigl& Burkert
2021
weainvestigated
on
2
=vcm
Myr
,
abo
e whichall.−20
filaments
wouldstart tocollapseTh(Clarke&
Whitworth2015;
Hoemannetal. 2023radially.
)
influ
l
conclusions are drawn in Sections 6 and 7.
2
c
−3
)
g cm
s
−1ρ̄/1.50×10
2weight of 2.36 and a temperature of 10 K), G the gravitationalThis
raisesthequ
2.1 Hydrostatic
filaments
2cµcrit = end-dominated
≈
16.×
4M
pc
,0.28
(4)ratio ofgrad
of a filament
are formed:
calculated
for 10 K. The criticality f is then thethe⊙Filament
density
gradients
suppress
−20
−3
s
1.69
10is
g√Acm
0.42
+
Ainfinity results in the
Integrating
thedensity
distribution
radiallyuntil
H 2 = fragmentation:
,
(2)
G
)
with
=
L
/2
R
being
the
aspect
ratio
and
ρ̄
the
mean
densitt=constant,andρ
thecentraldensityofthefilament.
The
filamentinmostcases.
Ho
w
c
col
=
Myr
,
(9)
line-mass% solution&Gρ̄ divided by the critical line-masscollapseasπGρcallow
Considering filaments
isothermal
cylinders,a hydrostatic
maximumline-mass
for
which
a
hydrostatic
solutionexists,
giv
en
by
f
ρ
1
2
κR
filament.
This
is
now
the
o
v
erall
collapsetime,
thetime-sisradially
constrainedbythepressureequilibrium
betweenthecdonevaluewith a sharp
above which all filaments would
to collapse radially. The
2 BA S I C P R I N C I P L E Sµ +0.28starttólkiewicz=external−1 −) 1for the radial profile
alreadyfoundby speedStod(s1963
and
in th
edge
0.42theA end cores move to the centre.
which
two
andwascs being
the
sound
(0.19
km
for
a
mean
molecularboundarypressure
andthepressure
P
=
P
.
Its
line-mass
f
=
.
√
ext
b
=
Myr
,
f
|
a
|
2 1,cminfluences
the acc
calculated
for
10
K.
The
criticality
f
is
then
the
ratio
of
the
actual1,2 ‹
1,3
2,3is with
Suppressing
theend-dominatedcollapse
3999
κ
=
1.66
and
a
being
the
accelerationofthecentre
of
mass
of
−20
−3
2
c
µ
Elena
Hoemann,
StefanHeigl
andAndreasBurkert
cm
ρ̄/
1
.
50
×10
g
cm
crit
Ostriker (1964): weight of is2.36
sthe gravitational−1
Although the edge effect is(4)predicted by theory, it is onl
givenandby itsoverall mass divided
by=itsGlength:a temperatureof
10critK),)µ≈ line-mass
16the.4 Germany,M, endbyregion1⊙ pc
2.1 Hydrostatic filaments
gradientthe
could infl
dividedthe criticalline-mass
1, D-81679 Munich,ät München, Universitäts-Sternwarte, Scheinerstr.
filaments
a
=
a
(
l
/2
−
κR
/2).
In
beginning,
cm
−20
−3
! 2 Ludwig-Maximilians-Universit$G
observed(Zernickel
et aspectal. 2013;theKainulainen
et beal.calculated
2016
; De
The
radialboundary
densityof ratiothefilamentcan−2 and
.
69
×
10
g
cm
"constant,
#
with
A
=
L
/2
R
being
the
and
ρ̄
the
mean
densi
ρ
the
centraldensityofthefilament.
The
filament
Max-Planck
Institutefor Extraterrestrial
Physics, Giessenbacherstr.
1, 1
D-85748
Garching,
Germanyc2
M
r Cluster ORIGINS,
2.3
allowing
perturbat
3 Excellence
=
Myr
,
(9)
to perturbation
theory(Larson
2.2
Edge
effect
two
cores
mo
v
e
in
free
fall
and
then
are
sloweddown
byrampressureµBoltzmannstrasse
2,D-85748
Garching,
Germany
µ
=
.
(3)
et
al.
2019
;
Bhadari
et
al.
2020
;
Yuan
et
al.
2020
;
Cheng
et al
Considering filaments as isothermalhydrostatic
solutionρ(r) = ρccylinders,, constrained(1)
1 + is aradially2ve whichfilaments
would filament.
startρbto=collapse
valuecollapse time,This
isfradially.
now
theTheoverall
the timeby the pressureabo
equilibriumthef ρc allfbetween
=
.
(5)
l
ρ
(1
−
)
,
c
in
the
next
sectionH
Miyamafor1992
Gehmanet wasal. already found by Stodólkiewicz (1963) and(Hoemann,
Heigl & This
Burkert2023).question
With thistwo-phaseapproach,raisestheratiowhywecentre.
do notobserve theInedgµcrit
the ;radial
profilecwhich
twotheend
cores
veactual
to theisdistribution
calculated
K.
The criticality
f isthethen
of mo
theof
pressure
andthepressure
Pext =until
Pbfor.infinity
Its10
line-mass
Since there is boundary
no hydrostatic
solution
forexternal
thedistribution
density
Integrating
thedensityradially
results
in
the
in
most
cases.
Ho
we
ver,
most
of
the
formertheoretical
studi
the
longitudinalcollapsetime-scaleofafilamentcanbedeterminedmostly Ostriker
not observ
ed,
sev
eral
with
κ
=
1.66
and
a
being
the
accelerationofthecentre
of
mass
of
and
the
radiusof
thisboundary
is
given
by
cm
Accepted
2023
August17.
Received2023
August16;
inoriginal
form
2023
July
20
(1964):
with
H the
Althoughtheedgeeffectisbe
predicted
by bytheory, it isfilamon
is given
its overall
mass divided
by itslength:dividedbythegivcriticalalong
the scalemain height:
axis
of a byfilament,
theline-massfilament
is expected
to collapseTheboundary
densityofthefilament
can2023
calculated
maximumforwhicha line-masshydrostatic
solution
eradialxists,
enby line-mass
done
with
a
sharp
density
cut-off
at
the
ends
of
the
filamen(Clarke&
Whitworth2015;
Hoemannetal.
)were detected (Tafalla
&
Hacar
the
filaments
end
region
a
=
a
(
l
/2
−
κR
/2).
In
the
beginning,
the
%
&
cm
1
/
2
!
2.3 Perturbation
observed(Zernickel
et al. 2013in; Kainulainen
et al. 2016line; De" r #2 $−2 longitudinally.
The M
acceleration along
a filament due to itsµselff acceleration2
2cs2
influencesthethe end region. A smoother
2
2
c
2
4000
E.
Hoemann,
S.
Heigl
and
A.
Burkert
two
cores
mo
v
e
in
free
fall
and
then
are
sloweddown
byrampressure
R
=
H
.
s
−1
ρ
=
ρ
(1
−
f
)
,
(6)
42√+0.28
A et al. 2019; Bhadari et al. 2020; Yuan(5) et al. 2020;(10)H= was,Aalready(2).
bt
c 0.(3)
ρ(r) = ρc 1 +, gravity
(1)
et a
B Sµ
T R=
AC
T . µcrit = by Burkert≈ 16.4 &M⊙Hartmann
pc , f (=(4)
investigated2004):− f influence and weaken
col =gradient 1could
theedge Cheng
effect sign
µincritGρ̄πGρc The onset ofl star formation
In
competitionto tHe effect itself cannot be stopped
G(Hoemann,
Heigl&
Burkert
2023
).
With
this
two-phaseapproach,
is
set
by
the
collapseoffilaments
theinterstellar
medium.
Fromatheoretical
point
of
view,
⎡
This0.09
raises thequestion of
why faster.
we doThis
notisobserve
the
edg
)
)%
f
(z)
is the perturbation ⎤strength which is observed to be typicallyallowing
perturbations
to
grow
what
we
will
inv
and
thedensity
radiusof
thisboundary
is
given
&
%
&al.Thefilament
can
also
4000
E.
Hoemann,
S.
Heigl
and
A.
Burkertan isolated
cylindricalfilament
forms
coresthevialongitudinaltheedgeeffect.
Duetotheself-gravity
ofavaluefilament,
thestrongincreasein by
2
2startgion of withthe filament.
Ho
we
ver,
collapsetime-scaleofafilamentcanbedetermined
−1
abo
v
e
which
all
filaments
wouldtocollapseradially.
The
Integrating
thedensity
distribution
radially
until
infinity
results
in
the
(Roy
et
2015
)
and
λ
is
the
dominant
fragmentation
length(seedomradialboundaryofthefilament
canbecalculated
by
0
.
42
+
0
.
28
A
H the zscale height: and cs being the sound speed
inismost
cases.
wever, most of the former theoretical
(0.19
km
s for
al meanmolecular
l
inMyr
the
next
sections.
=E), Ho(11) stud
2pile-up
2isacceleration
at bothisends
leads+toR
afor
of Thematterwhich
collapses
into
cores.
Ho√
we
ver,
this
ef−20
fect
Fischera&
Martin2012, Appendix⎣2maximum⎦then%&1g/2023−3 ) observed. Most
line-massa = −2
πG
ρ̄
z
−
+
z
+
−
z
+
R
,
calculated
10
K.
criticality
f
the
ratio
of
(Clarke&
Whitworth2015;
Hoemannetal.
y a density gradient,such that weight2 rarelyρ̄/the1by.50actual×10cm
line-mass
for
which
a
hydrostatic
solutionexists,
giv
en
done
with
a sharpdensity cut-off at the ends of thebyfilamenof 2.36
and amodels
temperatureof 10K), cut-off
G2 thegravitational! edge of the 2filament,f
2 3whereas 4a" smoother2
theoreticalconsiderasharp
densityatthetransitionismorerealistic2
ρ6line-mass
)f ,=
(6)
2cs is a natural! withλdom=observed.b25=− 6ρ.89
fto −
+be
9f
.18
− 3H
.44f 0.09
FWHMcyl(13)line-mass
dividedby
thecritical
c (1
is
the
perturbation
strengthwhich
istypicallyR.
(7)
the filament.
This
2
0of
.42
+0.28
A
also2decrease
the acceleration
attthe ends
thefilament.
We show withthat thecanbe significantly
Aedge
=1 −
Leffect/2fR being
the aspectratio(10)and ρ̄ inthethemean
influences
the slowed
acceleration
enddensity
region.ofAthesmootherconstant, andandρ cwould
the centralofThe
filamentH =,(2) the filament.
filam√
cs2 densitycol =2.3 Perturbations−1
using atheatver,
half-maximumof thecylindricalfilament
Gfull-widthρ̄ fragmentation
πGρcWe present
µ
f
(
z)
=
f
+
ϵc
1
(Roy
et
al.
2015
)
and
λ
is
the
dominant
length(see
ow an edge effect.
down
by
a
density
gradient,
althoughnotcompletelyv
oided.
Ho
we
this
allo
ws
perturbations
insidethefilament
to
gro
w
0
µ
=
≈
16
.
4
M
pc
,
(4)
dom
(8)
crit
⊙ equilibrium between
This(5)isgradientnow
overallcollapsetime,
the time-scaleon signinfluenceandweakenthe edge effectisradially constrainedbythethe radius
offilament.
thisboundary
given
bycouldfthe= pressure. to
FWHMcyland. Following
a perturbation
analysis(Nagasawa
1987), thetheistheperturbation
strengthwhich
isobservedbetypically
0.09
G
faster
than
the
edge.
We
determinethecriticaldensity
gradient
for
which
the
time-scales
are
equal
and
find
it
to
be
of
the
order
of
µcritpressure&
Martin2012, Appendix
E) .42+line-massunder whichedge should
In competition
to theperturbations
wi( be boundary?)athe0.28
A %time-scalecannow
be determinedon whichsuchperturbation
allowing
perturbations
toedgegroweffect,
faster.line-mass
This is what
we will inv
which
the
twoa endcores
move to thecentre.
pressure
andexternal
PalongPthe&
and cthe(0.19
kmmeans−1 Fischeraformeanmolecularext
b . 0Its√s being the sound speed
1/2 is an==start, The(11)
severaltimes
thefilament
radius.
Hence,
thedensity
gradient
at−20
the
ends
ofMyr
a filament
essential
parameter
for fragmentation
with
ρ̄
the
densityandzthepositionfilament.
The
abo
v
e
which
all
filaments
wouldtocollapseradially.
value
!
"
−3
would
grow
into
a
collapsing
core(compareHeigl,
Gritschneder&filament
can
also
grow.
Here
we
characterizeaperturbatio
with
f
the
unpertf(Roy canet 2.36al.
2015
λdombyisofitsthe
dominant
fragmentation
length(sees perturbations
grow
faster
ρ̄/of
1.50
×10
0
The
radialboundary
densityfilament
be calculatedby
2 length:
4 g cm canthe effectnext sections.
Althoughthe inedge
is predicted by theory,issteep
given
overall
dividedbyweightofanda) andtemperature
10
K),
Gmass
theof6end
gravitational
and
the
low
rate
of
the 9edge
effectbef
explained
gradients.
Rthen
=
Hthe
.actual(7) it is only rarelyt
λdom=.observed
25−
6.cases
89
fits
+
.18
f 3Burkert−could3inwards
.the
44
FWHM
(13)
2020
; naturally
Hoemann
et al.cyl2021)by shallow
increase
of
a
at
the
leads
to
a
large
velocityalongedge
line-massbyiscalculated
for
10
K.
The
criticality
f
is
ratio
of
the
1
2×1.66𝑅
1
ons.
1%the
− &faspect2
observed
al. 2013
al. 2016;filament.
Dewanganz = 0 c
withand
Aindynamics=theseL/2R#$
being
ratio(Zernickel
andmean
densityof the% et&$Fischera2012, Appendix
E)constant,and&ρ cMartinthecentraldensityoftheandfilament.
filamentρTheρac (1−thef )kinematics,of matter(6)ρ̄ the et
1 structure.
1!; KainulainenKeywords:
stars:
formation– ISM:–regions.
ISM:
b =
M
the
main
axis
as
a
resulttopileup𝑡=(−1)∝
line-mass
dividedbycriticalline-mass
t
!"#!
t
=
τ
log
−
1
(14)
.
s: First,isthe
basic principle
pert
domof the cylindricalusing the 𝑎full-widthat half-maximum
filament
2πz ; Cheng et al. 2021).
µthe=.
(3) theet2019; Bhadarietf (al.
eton
al. 2020
filament.
This is fwhen
now
ocol
val.
erall
collapsetime,
the2020ϵ𝑙 pressure𝑅radially
equilibrium
betweenthe𝑓 of by!constrained"z)=time-scale
f0 ; Yuan
1 + ϵ cos,The
resultingclumps
will
begin
to
collapseontothemselves!"#!
2.3
Perturbations
l
2
3
4
and
the
radiusof
thisboundary
is
given
by
µ
−+pressure
1.66)|
and perturbations
are6discussed
FWHM.
Followingaperturbation
analysis(Nagasawa
1987
),
the
λ
cyl
dom
This
raisesthequestion
of
why
we
do
not
observethe edge effectHereτdenotesthedominant
growthtime-scale,
which
can
be
which
the
two
end
cores
mo
v
e
to
the
centre.
boundary
external
P
=
P
.
Its
line-mass
λdom =pressure.25and− 6the.|𝑎(
89
f
9
.
18
f
−
3
.
44
f
FWHM
(13)
dom
b. mass. At the samef ext
= line(5)
2 reachthethedensity
2critical
they
time,
the
end
clumps% be determined&1/cyldeterminedwith
Fischera&
Martin(2012),
AppendixE
2
Integrating
distribution
radially
until
infinity
results
inthethe edge
effectcanneverbetime-scale
can
now
on
which
such
a
perturbationµ
suchpronouncedends
arerarelyobserved(Zernickel,Schilke&
Incompetition
to
the
edge
effect,line-mass
perturbations
cases. byHotheory,
wever,
ofunperturbed
therarely
former line-mass
theoreticalstudieswere
Althoughedge
effectinismost
predicted
itf0 istheonly
is given by its overall mass divided
byitswithmost
andz thecoordinate wal
1 length:Iards
N T ROD UcritCTIO
N centrefof the filament,
! the
" '
still
mo
v
e
inw
tow
ards
the
destroying
2the3
4 et al. 2016
Smith
2013
;
Kainulainen;
Dewangan
et
al.
2019
;
BhadariR=H.
(7)
maximumline-mass
for
which
a
hydrostatic
solutionexists,
by⁄/τdom= 4.08giv− 2en.Heigl,99f+ 1Gritschneder.46f done+ 0.40f
4
π
Gρ
.
(15)
derivationusingof the criticaldensity
𝑙
2
c
would
grow
into
a
collapsing
core(compare&filament
can
also
grow.
Here
we
characterizeaperturbatio
filament.
=
0 corresponds
centre whichof the perturb
witha sharpdensity
cut-offatthe ends of tothethefilament(Zernickel
2013;byKainulainenal.), which
2016;zDewangan
1density−f
the full-widthat +half-maximum
of1
cylindricalfilament
The
radialboundary
ofobserved
thefilament
canet2011
becalculatedObserv
ationssho; wToal
that
thethemolecularinterstellar
mediumisG
peral. et2020; Yuanet al.
2020; Chenget al.et2021
leads
to theMgrowá,ázquez-Semadeniómezfilament(Bastien
1983
V&;areal.
0.42
0.28
t
a
in the filament
faster
than
Since
most
of
the
cores
found
withinfilaments(Arzoumanian
pert
edge
Burkert2020;
Hoemann
et
al.
2021
)
line-mass
by
2
é
vaded
by
filaments(Schneider&
Elmegreen
1979
;
Andr
et
al.
question
of
why
we
do
not
observe
end-dominatedfilaments
more
µ=
.
(3)
influences
theal.acceleration
end region.
A smoother density
𝑅 ∝et al. 2019; Bhadari et al. 2020
; Yuan et
2020; Chenginettheal. 2021).
2al.
cs 2012
Pon;;ρClarke&
2015).2014This
iset1987
theso-called−1
FWHM. Following
a=et
perturbationönyvoften.
al.
;inKtheesedgeet al. 2015), this is expected
to be one of the𝑡$%& = µcrit
l thecylvalidationé et al.
et ,al.
Andr). They
are2011
key),%&(1Whitworth−
f 2013)2#$, ;(Nagasawa
(6)
!
%
&$
≈2010
16
.Arzoumanian
4b =M⊙ρanalysis(4)4, followed
by
cpc
could
weakenthesignificantly
raisesthe channels
question
ofgradient
why
wetheoretically.
donotinfluenceobserve andthe edge
effectmost
core
formation,
1This
𝑎log!"#!
2πzedge effectthe formationof starssince
pre-stellar
cores
areimportant
mainlyFor formagnetized,
disc-likeclouds the formation of an outer ring
G𝐺 𝜌understanding
MNRAS 525, 3998–4007
effect(Burkert&
Hartmann
2004).1suchInanearlierpaper(Hoemann,̅foundIntegrating
the density
distribution
radially
until
infinity
results
inthe, the results
are discussed
time-scalecanandnow
be determined
on
which
a
perturbation
t
=
τ
−
1
(14)
.
f
(
z)
=
f
+
ϵ
cos
1
,
pert
dom
0 studies
allowing
perturbations
grow
faster.
This is what we will investigate
withinthese
cold, dense,
filamentary
structures(Arzoumanian
iswe
thever,pendant
theeffect. Li theoretical(2001) toalreadyfound
thatwerecases.
mostof edgetheformerandtheradiusofstart
thisboundary
ismostgiven
byTheHo
ttoEring>tf).to
ϵ insev
λdom
Heigl& Burkert2021
we
ony spanwhich
time-scale
the
ends
pert
abo
vaehydrostatic
which
all
collapseradially.
value
et al.filaments
2011; K)önesinvestigatedet al.
eral
ordersofthecreationofGaedge
isF E
correlated
to the profile of the disc andline-mass for which
solution
eyvwould
xists,
giv2015
enbyTheand 7. maximum
3
T
H
E
I
N
E
V
I
TA
B
L
E
E
D
E
F
C
T
in
the
next
sections.
done
with
a
sharp
density
cut-off
at
the
ends
of
the
filament
which
would grow into a collapsing
core(compareHeigl,
Gritschneder&%et al. 2022&) frommagnitudeparsecs
in oflengththe sound speed. However, for filaments Seifried & Walch (2015)/2f hundredsa filamentfor
are10
formed:isof
calculatedK.(HacarThe
is dominant
then ofthe
ratio
the
actualwith
thetounperturbed
fregime 1(AndrHereτ domcriticality
denotesthetime-scale,
which
can
bepeakregion.
Ingrowth
formertheoreticalstudies
Clarke& Whitworth; Hoe2
influencestheaccelerationin theend
Af0 leadssmootherdensity line-mass and z the coordinate a
é
down
to
the
sub-parsec
et
al.
2010
;
Molinari
etal.
showedthat(e.g.
an initial
densityin 2015
thecentrea centrally)
2cs 2020; Hoemann
R
=
H
.
(7)
−1
%
&
Burkert
et
al.
2021
)
line-mass
dividedby
thecriticalline-mass
mann
et
al.
2021
),
filaments
areoften
considered
to
have
a
sharp
2010; determined2013; Goodmanetal. 2014
Matternet (al.
2018;), Appendix
dominated collapse
instead of anone.
In caseoffilament.
z =0 corresponds
to the centre of the perturbµcrit =
≈ 16.4 M⊙ pc ,(4)
with&; Martin2012
E weaken1−
f Fischera
gradient
could
influence
theend-dominatededge
effectsignificantly
1 Hacar etet2al.
κRedge, thus,
a hardcut-off
at theandends
the ratios
filament
which
makes
G
Schmiedeke
al.! 2021). Besides their obvious connection
to starshortfilaments
with ofaspectof 3:1the longitudinalcollapse was
2.3
Perturbations#$ tedge =
%
&
'
−
1
"
µ
predictions
possible(Burkert
Hartmann
2004is).
Ho
weallowing
grow
This
what
weshowed
will investigate
3analytic
oninvestigatedby&faster.
Keto& Burkert(2014). ver,Theythat forf τdomthere(5)4|.remain=|astill
08 The−many
2.value
99openf 2 questions+ 1.46f the+formation,
0.perturbations40
f 4 /also
4πto
Gρ
(15)
1f = µstart. to
1 formation,
cm
c.
above which all filaments would
collapseradially.
a more realistic setup
would
be a smoothtransitioncollapses
to the surrounding
evolution,andcollapseoffilaments.
lowline-masses
thefilament
longitudinallyintoacentralIn
competition
to
the
edge
effect,line-mass perturbations within thein thenext
sections.
t calculated= τdomlog− 1crit) then(14)
. the ratio
gas.collapse
Such a isdensitycore
gradient
to athroughdecreasea inacceleration
MNRAS 525, 3998–400
Particularly
the formationof coresand their consequent
that would
beginslead
to gocomplexpattern of oscillations.
is pertfor 10 K. Thefcriticality fϵisoftheactual
−20
−3
Since
of
the
cores
are
found
filaments(Arzoumanian
𝜖:
~0.09
can
alsohasgrow.
we characterize a perturbation of the
1interesting.69 × 10
gmost
cm
at
thewithinends of thefilament.
since thisisof
the thestagefilamentof early starformation.
For
isolated
Although suchfilament
a longitudinalcollapse
not yetHere
been confirmed
The
radialboundary
densitycanbecalculatedby=Myr(9)to dueline-mass
dividedby the critical line-massin accelerationtothe etend
filament oscillations of theylinders,
a hydrostatic
solutionfilaments,
supposed
via twocompetingby observations
(Roy
2015
et al.
; Kisönyves, toethappenal. 2015), Thethischange
is expected
beline-mass
one
ofal.profile
theby )ofthethepredictedffragmentation
ρ2011
c
2
is
depicted
in
Fig.
1
.
The
top
row
shows
three
simpledensity 2006; Aguti et al.
2.3
Perturbations
processes:growth
the edge effecttime-scale,(Bastien 1983; Burkertwhich& Hartmanncan
2004; beresultingµ τ domdenotesdominant(6) cores are (Redman, !Keto & Rawlings%&$ólkiewicz
d by StodHere(1963) and ρtheb = ρc (1 − f ) ,mostimportantchannels for coreformation,
theoretically.

Filament fragmentation: density gradients suppress end-dominatedcollapse
Figure 6. An exemplary simulation environment to validate the criticalgradient, in this case with R = 0.1 pc, f = 0.5, and dl/R = 6. The upper
part shows a section of a density slice, the second plot the criticality along
the z-axes, and the last one the gravitational acceleration along the z-axes.
Inside the filament perturbations with ϵ = 0.09 are set up followed by a flat
profile in the end region and a linear density gradient, which in this example is
dl/R = 6. Since there is a density cut to the ambient medium in the simulations
at the filament boundary density ρ b , the gradient is curved in the line-mass
projection.

1,2,3I O N S
5 VA1L,3I and
DATAndreas
ION BY SI M U L AT
Elena Hoemann ,1,2 ‹ Stefan Heigl
Burkert
1 Ludwig-Maximilians-Universität München, Universitäts-Sternwarte, Scheinerstr. 1, D-81679 Munich, Germany,

=t

AMR

RAMSES

R=0.1 pc, l=3 pc,
f=0.5, dl/R=6, 𝜖=0.09
àρlines present less dominant perturbations whereas the red curves necessary.
a =

−1

log

−1

.

(32)

(f<0.4)

FWHM~0.03 pc

f>>1

Downloaded from https://academic.oup.com/mnras/article/525/3/3998/7246898 by D

Downloadedfromhttps://academic.oup.com/mnras/article/525/3/3998/7246898
Dankook
University/ Lib
Downloadedfromhttps://academic.oup.com/mnras/article/525/3/3998/7246898
byby
Dankook
University/ LibraA=l/2Rµ/µcrit

dl/R

⾒t

µ/µ

~0.05 pc

Downloaded from https://academic.oup.com/mnras/article/525/3/3998/7246898 by Dankook University / Library user on 10 Oct

For validation, we performed isothermal hydrodynamic simulations
with the adaptive mesh refinement code RAMSES (Teyssier 2002).
3 Excellence Cluster ORIGINS, Boltzmannstrasse 2, D-85748 Garching, Germany
With a second-order Gudonov solver the Euler equations are solvedin their conservative form using a MUSCL (Monotonic Upstreamedge pert
Suppressing the end-dominated collapseCentered
4003Scheme for Conservation Laws, van Leer 1979), the HLLCAccepted 2023 August 17. Received 2023 AugustSolver
16; in originalform 2023 July
20 Contact, Toro, Spruce & Speares
4002 E. Hoemann, S. Heigl and A. Burkert(Harten-Lax-van
Leer1994), and the MC slope limiter (monotonized central-differenc,
Leer 1979).
4004 A BE.
S. Heigl and A. Burkert
We set up filaments with perturbations in the line-mass as disS THoemann,
RACT
playedin Fig.of6.filaments
In the upper
panel,
a densitymedium.
slice through
The
onset
of
star
formation
is
set
by
the
collapsein theinterstellar
Fromthea theoretical point of view,
Table 1. Results of the simulation ‘sim’ in comparison to the expected result
necessarytosuppress
the edgethe
effect.
However,the linear gradientfilament
is
given
and
in
the
lower
panelsaccordingcriticality
an isolated
cylindricalfilament
forms
effect.
Due to the self-gravity
of a filament,by the model‘exp’. ‘e’ means
simulations
with dominant
edgecores
effect orvia‘p’ the edge
is a goodfirst approximation,
especially for filamentswith f ≥ the
0.4. strong increase inandacceleration
along
the
z-axes.
Thebendingofthegradient
in
configurations
where perturbations
on thesame totime-scale
or faster.
acceleration
at bothgrow
ends
leads
a pile-upof Thematter which collapses into cores. However, this effect is rarely observed. Most
the end regionis due to the cut off at boundary density, which willedge effect and perturbation time-scale, as determined in the simulation,aretheoretical models consider a sharp density
cut-off atat thethe edgeof thefilament,
whereas werea smoothergiven in comparison to their expected values in brackets. An output
be time-step
discussed
section.
Perturbations
inserted transition is more realistic
6 DendI S Cof
U SthisSIO
N
and would
decrease
the acceleration
atequation
the
of
the
filament.
We
show
that
the
edge
in the simulation
takes !also
t = 0.1
Myr. The length
was always via
given
by
fiveends(12) going over into a flat region at the end,effectto notcan be significantly slowed
times thedown
dominant(equation13). although not completelyFor
the calculationof thethis
criticalseveral approximations
bywavelengtha density
gradient,
avoided.
Howe
ver,
allo
ws perturbations
disturb the edge effect.
At thefilamentend
the gradient,lineartransition
intoinside the filament to growhave been made. However, Fig. 4 already showed a comparison
thesurrounding
densityisdepicted,
in
this
case
with
a
gradient
faster
than
the
edge.
We
determinethecriticaldensity
gradient
for
which
the
time-scales
are
equal
and find it to be of the order off
R (pc) dl/R
tedge (Myr)tpert (Myr)exp
simbetween the numerical collapse times and the approximated collapse
dl/density
R = 6. gradientWe cut theoutermostpeaks
to notishave
an overlap
several times the filament radius. Hence,ofthe
at
the
ends
of
a
filamentanessential
parametertimes. For the tested cases, which show a reasonable spanover the for fragmentation
0.5
0.1
10
1.76 (2.02)
1.57 (1.55)
p between
p perturbations
andtheof filamentsedgeeffect.
Then,
we vofaried
thetime-scales
and
the
low
rate
of
observed
cases
of
the
edge
effectcould
benaturally
explainedbyshallow
gradients.
parameterspace
thecrossing
point
the
two
0.5
0.1
6
1.47 (1.71)
1.57 (1.55)
p
p
Figure 5. The grey line depicts the critical centre of mass acceleration
linear density
gradient
at both ends
transition between
is approximatedwell. toIn validate
addition, thethe determined
curve of the critical
0.1
2
1.18 (1.31)
1.57 (1.55)
ee(equation 32) depending on criticality. The different markers give numeric0.5
perturbationsandedge
effectasgiven
byequation(30calculated). We consider
Keywords:
stars:
formation–
ISM:
kinematics
anddynamics–
ISM:
structure.
crit
gradient
fits
well
with
the
numericallycurve (Fig. 3).
0.05
60.78 (0.84)
0.78 (0.78)ppresults of the centre of mass acceleration of different filament end profiles0.5
toAlthough
be dominant
if all corescollapsewithinfor
a timespanwith lower
0.05
101.27 (1.94)
1.67 (1.74)
p perturbationse
we observea systematicdeviation
filaments
given in the subplot. Above the grey line the edge effect is expected to be0.3
Figure 7. For simulations with f = 0.5, R = 0.1 pc the evolution of the
0.1
10
2.55 (3.99)
3.14 (3.49)
p of 10 eper cent of
dominant, whereas below, perturbations would grow faster. In comparison,0.3
aspectA, which
expected
to the firstused sinceapproximation,
eachratioother,or the isinnercoresduecollapse
criticality for the end cores against the criticality of the perturbations are
the violet crosses show the acceleration measured in the simulation by Heigl
the values
doconverge
towards
larger aspectratios.
However,shown
the for three different gradients. All simulations start at fpert = fedge = 0.5
thenno
pronounced
edges
can
be
detected.
In
contrast,theedgeeffectet al. (2022), where the edge effect was suppressed in agreement with oursuchpronounced
ends
are asrarelyobserved(Zernick
Schilke & f-values with time. Each marker represents one
deviationis limited
and
can cores
be interpreted
a transition
regionandmove toel,larger
maximalisconsideredtobedominating
if
the
end
collapse
beforethemodel. The used parameters were: ϵ = 0.09, R = 0.02 pc, l = 0.94 pc, dl/R =1 I N T RO D U C T I O N
betweenthetwoeffects.
Thus,
theapproximations
seem
applicable
Smith
2013
;
Kainulainen
et
al.
2016
;
Dewangan
et
al.
2019
;
Bhadarisimulation
time-step("t = 0.1 Myr). Above the solid black line, perturbations
9.
Figure 6. An exemplary simulation environment to validate
the criticalinnercores.
toisourpercase.
Observ
ationsthatthef =molecular
medium
et
al.
2020
;
Yuan
et
al.
2020
;
Cheng
et
al.
2021
),
which
leads
to
the
gradient,
in thiscase showith wR=
0.1 pc,0.5, and dl/Rinterstellar= 6. The
upper
grow
faster
in
f
than
the
edge. A clear edge effect can only be seen for dl/R =As a first test, we performed simulations for f = 0.5 and different
size of theof
endwhyregionrealisticis difficult
to2,say.
partvadedshows abysectionof a density
slice, the second
the criticality
along; André Whether
filaments(Schneider& plotElmegreen
1979
etFig.
al. 7.thequestion
wedldois notobserve
end-dominated
more
for the filaments
others perturbations
grow on similar time-scales or faster, which
5 VA L I DAT I O N B Y S I M U L AT I O N S
gradients,
given
in
Here
we
plot
the
criticality
of
the
cores
b
the z-axes, and the last one the gravitational acceleration along the z-axes.
Althoughthere often.
are some observations of the line-mass distribution
fits the prediction.
é
2010
;
Arzoumanian
et
al.
2013
;
Andr
et
al.
2014
).
They
are
key
in
Inside the filament perturbations with ϵ = 0.09 are set up followed
by a end
flat fend against the criticality of the cores formed fromat theFor validation, we performed isothermal hydrodynamic simulations
along filaments as Roy et al. (2015); Cox et al. (2016); Kainulainen
understanding
ofgradient,stars sincepre-stellar
aremainlyFor onemagnetized,
clouds in
the formation of an outer ringprofilein the end regiontheand formation
a linear density
which inthisexamplecoresisperturbations
fpert
Every
outputofdisc-likethe simulation
with the adaptive mesh refinement code RAMSES (Teyssier 2002).
et .al.
(2016marker); Yuan is
et al. (2020); Schmiedeke
et al. (2021), the end
d
l
/
R
=
6.
Since
there
is
a
densitycut
totheambient
mediuminthesimulations
found within these cold, dense, filamentary structures(Arzoumanianisthependant
totheedge
effect.
(2001) alreadyfound that
FigureTop panel: The
approximation
critical
gradient are(equation
With a3.second-order
Gudonov
solverofthetheEuler
equationssolvedtime-stepsof "regiont = 0.1
Myrwith
R = 0.1 pc,
l =enough0.3 pc, resolved
and f =Liagreesisoften
notdisplayed
or
not
to
make
a with the fact that in this case a dominant edge core would beatthefilament
boundary
densityρ,thegradientiscurvedintheline-massb30in) istheir
depicted
for atidifferent
criticalities
grey. The(Monotonic
different shades
depictet al. 2011; Könyves et al. 2015). They 0.5.
spanFor
sevdl
eral
ordersof dlstatement
thecreationof gradient.
a collapse
ring isThecorrelated
to theisexpected.
profileofThisthe disc
and the critical gradient given in Fig. 3 whichconservve form
using a inMUSCL
Upstreamprojection.
/
R
=
10
and
/
R
=
6,
both
cores
withinone
validates
reasonable
about
the
best
examplethesolutions for different filament aspect ratios, indicated by the given number.
Centered Scheme for Conservation Laws, van Leer 1979), the HLLCmagnitude(Hacar
et
al.
2022
)
from
hundredsofparsecs
inlengththesoundspeed.
Ho
we
ver,
for
filaments
Seifried
&
Walch
(
2015
)
The parameters (α, β) are fitted to the numerical solutions for a long filamentof Barnard
5 bySchmiedeke
al. the(2021). Here a clear
time-step, thus, observation
perturbations
and edge
effectcollapseet on
same
predicts the transition between
end core and perturbation domination
Solver (Harten-Lax-van Leer- Contact, Toro, Spruce & Speares
downtosametheassub-parsec
regime(Andrtoétheetapproximated
al.
2010
;
Molinari
etal.
showedthataninitialdensitypeakinthecentreleads
toacentrally(R = 0.1 pc, l = 200 pc), given by the black dots. Below the curve, a dominantlinear
trend
is
seen
at
the
end
of
the
filament.
The
decliningregionisFigure4. TheFig. 2, but adding
a comparison
time-scale and no clear edge effect can be detected. Only for dl/R =for these parameters to be at dl/R = 5. Further validations for different
1994
),
and
the
MC
slope
limiter(monotonizedcentral-differenc,edge effect is expected, whereas above the curve, perturbations grow faster,edge 2010
effect ;time(equation) with ;thedotted lines.etInal.
the2014
first row,
the
Hacar
et al.272013
Goodman; Mattern
et al.
2018; clearly
collapsean
one.canIn be
caseof in Table 1.
ofcollapses
the
order
ofdominated∼0.05
pc, thanwhereastheinstead
FWHMofiswhich∼ 0end-dominated.03 pc, thus,
Leer). edge effect is expected. Lower panel: The critical gradient is
2 the
edge
effectfasterthe perturbations
parameters
found
and
thus1979
no strongdependence on criticality is depicted, in the second row on radius and in thein tothisendfilaments
region iswith
indeedbiggerthethefilament
Schmiedeke
et
al.
2021
).
Besidestheirobviousconnection
starcase, theshortaspectratiosthanof 3:1longitudinal collapse wasup filaments
withperturbations
inlinetheshows
line-mass
as disgivenWefor setdifferent
perturbation
strengths
ϵ. The grey
the observed
last on the filament length. At the crossing points between tedge and tpert theradius, which also
is needed
for perturbations
be as fastas).theplayedofinϵ =Fig.
In thea density
slicethroughtheformation,
still formation
remain time-scale
many open
questions on the formation,
investigatedby Keto &toBurkert(2014
They showed that forstrength0.096.(Royet al.upper2015),panel,in comparisonin bluelowerstrengthapproximationof thethereedge effectis reasonable.
andin red stronger
ones.
filament
is given
and in the lower panels the according criticality
edge effect in thespace
investigated.
However, theMNRAS 525, 3998–4007 (2023)
evolution, and collapse of filaments.
lowparameterline-masses
thewefilament
collapses
longitudinally
into a central
cores. This can be calculated without any additional approximations
and acceleration along the z-axes. The bending of the gradient infilamentisis highly
supercriticaland tothus,
outsidethe arange
where pattern
ourParticularly
theformation
of
cores
and
their
consequentcollapsecore
that
beginsgothroughcomplexofoscillations.
Figure8.
Critical
gradient
for
a
cut
off
at
ρ
given
by
the
squaresinlightbusing the same ansatz as before and solving
for the acceleration:
the end region is due to the cut off at boundary density, which willmodel is applicable. Further observations of filament end regionsinteresting
thisis theany
stageof early
star formation.
Although such a longitudinal collapse has not yet been confirmed
the former since
solution
without
cut off
in comparisonas black For isolated
be discussed at the end of this section. Perturbations were insertedorange and!
are necessary to make a reasonable comparison and to test our
t
=
t
(31)
pert
line
of ϵ =(Roy etatal.the2015), which
For filaments,
line-masses
f ≥ 0.4 the deviations
are small. Hover, for via
lighterfragmentation
is supposed
to we
happentwo competingby observations (Roy et al. 2015) the predicted oscillations of theviaindicatesequationthe(12observed) going value
over into
a 0.09flat regionend,
to notdots. edge
prediction.
! estimation" is not sufficient!! and"much" smoother gradients are
isdisturb
also useddemonstrate
thelinearsolution.
The bluethethetoedge
effect. Atthe
theconvergencefilament endofthetransition
intofilaments,processes:
theresulting cores are (Redman, Keto & Rawlings 2006; Aguti et al.
1
2κRedge−2effect1 (Bastien
1 1983; Burkert & Hartmann 2004;
2 Max-Planck Institute for Extraterrestrial Physics, Giessenbacherstr. 1, D-85748 Garching, Germany

Barnard 5.oup.com/mnras/article/525/3/3998/7246898 by Dankook University / Library user on 10 October 2023

dominant, whereas below, perturbations would grow faster. In comparison,
the violet crosses show the acceleration measured in the simulation by Heigl
et al. (2022), where the edge effect was suppressed in agreement with ourmodel. The used parameters were: ϵ = 0.09, R = 0.02 pc, l = 0.94 pc, dl/R =9.

Triple Spiral Arms of a Triple Protostar System Imaged in Molecular Lines
Jeong-Eun Lee1,2 , Tomoaki Matsumoto3 , Hyun-Jeong Kim4 , Seokho Lee4 , Daniel Harsono5 , Jaehan Bae6 ,Neal J. Evans II7 , Shu-ichiro Inutsuka8 , Minho Choi4, Ken’ichi Tatematsu9 , Jae-Joon Lee4 , and Daniel Jaffe7

1

Department of Physics and Astronomy, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea; lee.jeongeun@snu.ac.kr
2
SNU Astronomy Research Center, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea
3
Faculty of Sustainability Studies, Hosei University, Fujimi, Chiyoda-ku, Tokyo 102-8160, Japan
4
Korea Astronomy and Space Science Institute, 776 Daedeok-daero, Yuseong-gu, Daejeon 34055, Republic of Korea
5
Institute of Astronomy and Department of Physics, National Tsing Hua University, Hsinchu 30013, Taiwan
6Department of Astronomy, University of Florida, Gainesville, FL 32611, USA
7
Department of Astronomy, The University of Texas at Austin, 515 Speedway, Stop C1400, Austin, TX 78712-1205, USA
8
Department of Physics, Nagoya University, Furo-cho, Chikusa-ku, Nagoya, Aichi 464-8602, Japan
9
The Astrophysical Journal, 953:82 (13pp), 2023 August 10
Lee et ofal.Natural Sciences,Nobeyama, Minamimaki, Minamisaku,
Nobeyama Radio Observatory, National Astronomical Observatory of Japan,National Institutes
Nagano 384-1305, Japan
Received 2022 September 15; revised 2023 April 30; accepted 2023 June 9; published 2023 August 4

IRAS04239+2436, ALMA, 857 µm, 0.1" (14 au), 42-au
Continuum(), 0.14-0.20 Mo.

SO 107 cm 3, 3 km/sAbstract
Most stars form in multiple-star systems. For a better understanding of their formation processes, it is important toresolve theàSO
individual protostellar components and the surrounding envelope and disk material at the earliest

SO 88-77

possible formation epoch, because the formation history can be lost in a few orbital timescales. Here we present
Atacama Large Millimeter/submillimeter Array observational results of a young multiple protostellar system,
IRAS 04239+2436, where three well-developed large spiral arms were detected in the shocked SO emission.
Along the most conspicuous arm, the accretion streamer was also detected in the SO2 emission. The observational
results are complemented by numerical magnetohydrodynamic simulations, where those large arms only appear inmagnetically weakened clouds. Numerical simulations also suggest that the large triple spiral arms are the result ofgravitational interactions between compact triple protostars and the turbulent infalling envelope.
Uniﬁed Astronomy Thesaurus concepts: Star formation (1569)
The Astrophysical Journal, 953:82 (13pp), 2023 August 10

Matsumoto+'15

1. Introduction
Lee et al.

masses of the individual central sources, may be important(Matsumoto et al. 2015b). In addition to turbulence, theFigure 1. The SO integrated intensity (a) and velocity dispersion(b) maps.binary
The whiteand black(b) contours
continuum
Althoughor (a)
multipleformationaredepictthe themainmodesemission distribution. The
ellipse in the lower left corner represents the beam size. The contour levels and the beam size of the continuum image are the same in all ﬁgures. magnetic
The SO integrated
ﬁeld must also play an important role in theof B,starformation(Offneralso
etincreases
al. 2022),
main
formation
intensity map shows a strong ﬂow-like structure toward Source
and thevelocity dispersion
towardthe
SourceB. Theindividual line proﬁles
can be found process itself, and in the structures insidefragmentation
mechanism
remains
theoreticalin Figure 3; the SO lines around Source B are very broad, while
the line intensity
itself elusive.
peaks at theObservational
ﬂow-like structure andbeneathSource B. The integrated intensity issince cores with weak magnetic ﬁelds tend to form 3
proportional to the column density multiplied by the excitation
temperaturein the moleculargas, circumstellar
and thus it increases
toward the centralfragments,
source. The peakstudiesthat zoom
into young
andrapidly
circumbinary
multiplesystems
intensity of the SO line is more sensitive to the excitation disks
temperatureofthemolecule,soitbettertracestheshocked
gas
along
the
spiralarms,
aspresentedin (Machida et al. 2005).
in the early evolutionary stage are critical for under"
Figure 12(a).
! To study the formation process of a close binary system

AMR

𝜌 𝑅 =

SFUMATO,'(,(1.56 pc) , T=10 K,, R =0.05 pc, 𝑃(𝑘) ∝ 𝑘 12 , <M>=1,

standing the formation process, before their history is erased
the environmental
0 of the natal cloud may still"when" " effectthrough dynamical interactions (Lee et al. 2017; Moe &)*+⁄+system,(,.are# ) we observed IRAS 04239+2436 with theremainin+theand its comparisons with numerical simulations (Section 3).
3.5,levels (σ) for the SO, SO2, HCN, and HCO+ images#
Kratter 2018). The interactions between the inner envelope
andWe discuss the effect of the magnetic ﬁeld on multiple
Large Millimeter/submillimeter
Array
4.45, 4.15, and 4.17 mJy beam−1, respectively.Atacama
The moment
5 ⁄613 (ALMA)
3
the central protostellar system, as well as the interactions 'cr and several molecularformation and the implications of close multiple protostars for$4maps of the molecular lines were generated3by theimmoments in thesimultaneously
857 μm continuum
among individual members, are important in multiple formaplanet formation in Section 4, and the ﬁnal conclusion is giventask using a threshold of 6σ. The continuum
image
lines,
with wasan angular resolution of ∼0 1 (∼14 au at 140 pc).
tion, since protostars grow their masses most aggressivelyinin Section 5.
produced by using line-freesink
channels with uniform
weighting.
IRAS
04239+2436 is an embedded close protobinary with athe early formation stage, when most of the mass still remains
The synthesized beam is 0 14 × 0 08, with a PAof 4°. 3,separation
andprojected
of 42 au at a distance of 140 pc (Reipurthin the envelope. Therefore, simultaneous imaging of the inner
the rms noise level is 0.15 mJy beam−1.
et al. 2000). IRAS 04239+2436 has a luminosity similar to thatenvelope and the central compact multiple system is also
2. ALMA Observations and Results
The ALMA observation (Figures 1, 2, and 3) of04239
of IRAS
the Sunand very rich near-infrared (NIR) spectrum (Greene
required. One prominent structure that shows the interactions+2436 reveals two compact continuum sources &corresponding
IRAS 04239+2436 was observed using ALMA during
Lada 1996). A high–spectral resolution NIR observation
between forming binary
circumbinary
is the any
spiralcontinuum emission
to andSources
A and material
B, without
Cycle 3 (2015.1.00397.S; PI: Jeong-Eun Lee) on 2016 August
kinematically revealed a sub-astronomical-unit-scale gaseousarm
structure(Takakuwa
et
al.
2017;
Alves
et
al.
2019).
16 10.UT.
Five associatedspectralwindows
in Band
7 wereset
to arecoverwiththe circumbinary material ordisk
beyondwith atheKeplerian rotation around a central mass of
FigureShock velocity
with thearms (left) and number
density distribution(right).
The data
rotated so that theyassociatedare oriented in a face-on
view. The𝑃 𝜌 = 𝑐 𝜌(1 + (𝜌⁄𝜌 )r =2.45 au.

, ρ =10SO

cm , ρsink=10 11 g cm 3,

Triple Spiral Arms of a Triple Protostar System Imaged in Molecular Lines
Jeong-Eun Lee1,2 , Tomoaki Matsumoto3 , Hyun-Jeong Kim4 , Seokho Lee4 , Daniel Harsono5 , Jaehan Bae6 ,Neal J. Evans II7 , Shu-ichiro Inutsuka8 , Minho Choi4, Ken’ichi Tatematsu9 , Jae-Joon Lee4 , and Daniel Jaffe7

1

Department of Physics and Astronomy, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea; lee.jeongeun@snu.ac.kr
2
SNU Astronomy Research Center, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea
3
Faculty of Sustainability Studies, Hosei University, Fujimi, Chiyoda-ku, Tokyo 102-8160, Japan
4
Korea Astronomy and Space Science Institute, 776 Daedeok-daero, Yuseong-gu, Daejeon 34055, Republic of Korea
5
Institute of Astronomy and Department of Physics, National Tsing Hua University, Hsinchu 30013, Taiwan
6Department of Astronomy, University of Florida, Gainesville, FL 32611, USA
7
Department of Astronomy, The University of Texas at Austin, 515 Speedway, Stop C1400, Austin, TX 78712-1205, USA
8
Department of Physics, Nagoya University, Furo-cho, Chikusa-ku, Nagoya, Aichi 464-8602, Japan
9
Nobeyama Radio Observatory, National Astronomical Observatory of Japan,National Institutes of Natural Sciences,Nobeyama, Minamimaki, Minamisaku,
Lee et al.
Nagano 384-1305, Japan
Received 2022 September 15; revised 2023 April 30; accepted 2023 June 9; published 2023 August 4

SO

HCN

deprojected, i=55 °, PA=140° (CSDs)

The Astrophysical Journal, 953:82 (13pp), 2023 August 10

Abstract
Most stars form in multiple-star systems. For a better understanding of their formation processes, it is important toresolve the individual protostellar components and the surrounding envelope and disk material at the earliest
possible formation epoch, because the formation history can be lost in a few orbital timescales. Here we present
Atacama Large Millimeter/submillimeter Array observational results of a young multiple protostellar system,
IRAS 04239+2436, where three well-developed large spiral arms were detected in the shocked SO emission.
Along the most conspicuous arm, the accretion streamer was also detected in the SO2 emission. The observational
results are complemented by numerical magnetohydrodynamic simulations, where those large arms only appear inmagnetically weakened clouds. Numerical simulations also suggest that the large triple spiral arms are the result ofgravitational interactions between compact triple protostars and the turbulent infalling envelope.
Figureconcepts:
12. ComparisonsStar
betweenobservations andsimulations. (a) The SO 8 –7 peak intensities in units of mJy beam and (b) the velocities of the intensity peaks
Uniﬁed Astronomy Thesaurus
formation(1569)

sim arm 1 <150 au
sim arm 2, 3 <100 au

8

7

−1

in units of km s−1 toward IRAS 04239+2436. The source velocity is VLSR = 6.5 km s−1. The SO peak intensity and velocity maps were generated using a threshold of
6σ (0.02
beam−1). The contours in (a) and (b) present the continuum emission at 857 μm. (c) The gas column density distribution and (d) the velocity distribution
Figure 13. Kinematics of the spiral arms of IRAS 04239+2436. (a) Peak velocities (relative to source velocity) along the spiral arms. The colored solidlinesJywithof the numericalsymbols are the velocities at the peak intensity positions of the spiral arm features (Figure 12) of IRAS 04239+2436. The ﬁlled symbols are the velocities
at the HD simulation. The arms that are clearly identiﬁable are named in (a), and the arms that may correspond to the simulation are marked in (c). The
white dots
positions used for the polynomial ﬁttings, as presented with the solid black curves in (b). The colored dashed lines are the velocities of the arm featuresfrom(c)theand crosses (d) mark the positions of the three protostars. The simulation was rotated and projected consistently with IRAS 04239+2436.
simulation, where the deprojected radii are scaled by a half to match the total mass of the companions. The black lines are the Keplerian rotation curves of 0.15 MeThe
AstrophysicalJournal,953:82(13pp),
2023
August10
Lee et al.
(solid), 0.30 Me (dotted), and 0.60 Me (dashed). (b) Peak intensity positions and velocities of the spiral arms (colors) on the deprojected peak intensity formation
map of SO 88of– multiples and the following development of large
considered according to previous simulations (Matsumoto 2011;
77 (gray image). The symbol “×” marks the center position deﬁned by the center of mass of the two continuum sources, assuming an equal mass. The black solid lines
spiral
arms,
we
also
ran
magnetohydrodynamic(MHD)
Matsumoto et al. 2017). The MHD scheme has second-orderindicate the polynomial ﬁtting results over the peak intensity positions that may belong to the arm structures with coherent kinematics. Gray circles are drawn fromsimulations, varying the initial magnetic ﬁeld strength, andaccuracy in space to maintain numerical stability, while the HD
0 5 to 2 5 in steps of 0 5.

MHD
(Matsumoto et al. 2015b). In addition to turbulence, theAlthough binary or multiple formation are the main modes
1. Introduction

masses of the individual central sources, may be importantmagnetic ﬁeld must also play an important role in theof star formation (Offner et al. 2022), the main formation
fragmentation process
and accuracy.
in the structures insidepresent the results
15.
scheme itself,
has third-ordermechanism remains elusive. Observational
andin Figure
theoretical
The initial condition and numerical models
are the samesinceasThewithMHDweakmodelsexhibit similar
the HDfragments,cores
magneticﬁeldsevolution
tend totoform
studies that zoom into young circumstellar
circumbinary
those of theandHD model,but the magnetic ﬁelds
are included.
model; the ﬁlamentary
cloud fragments into cloud cores, which
multiplesystems (Machida
et al. 2005).
disks in the early evolutionary stage
aremagnetic
criticalﬁelds
forperpendicularunder- to the ﬁlamentary cloud
Uniform
undergo gravitational collapse. A multiple- or single-star
To study the formation
processof aof close
binarysystem on
are imposedas the initialsystem forms
at the center
the cloudcore, depending
standing the formation process, before
their history
is condition,erased because recent observawhen
environmentalof ﬁeld
the natal
tions have suggested that dense ﬁlaments
tendtheto havethe initialeffectmagnetic
strength.cloud
Modelsmay
withstill
weakerthrough dynamical interactions (Leeet al. magnetic
2017; ﬁelds
Moe(Pattle& et al. 2022). Three models
perpendicularmagneticﬁelds tendIRAS
to form
multiple systems
withremain in the system,
we observed
04239+2436
withthemoreKratter 2018). The interactions between
theinitialinnerenvelope
and of 1, 5, and 10 μG arewith themagneticﬁeld strengths
companions (Figure 15) and to have more prominent and more
Large
Millimeter/submillimeter
Arrayﬁeld strength is Atacama
often describedinextended arms (Figure 16). The models
with 0, (ALMA)
1, 5, and 10 μG
the central protostellar system, asconsidered.
well asThethemagnetic
interactions
simultaneously
μm continuum
severalterms of the mass-to-ﬂux ratio, and the criticalmagnetic ﬁeld in theform857
the systems
consisting ofand
three(triple), molecularfour (quadruple),among individual members, are importantinmultipleformastrength (Nakano & Nakamura 1978; Tomisaka
et al.
1988)
two (binary),and onerespectively.
lines,
withanis angularresolutionof (single)∼0 1stars,(∼14
au at 140 pc).
1/2
tion, since protostars grow their masses
mostinμG for our models. The
According to our MHD simulations, the natal cloud of IRAS
Σ = 20estimatedas Baggressively
cr = 2πG
IRAS
04239+2436
is an embedded
close
protobinary
witha a04239+2436
must
have
a
very
weak
magneticﬁeld, withinitialmagnetic
ﬁeld
strengths
of
1,
5,
and
10
μG
are
therefore
the early formation stage, when most of the mass still remains
projected toseparationof 42 au atmass-to-ﬂux
a distanceoflarger140thanpc μ(Reipurthnondimensionalratio= 2–4, to form
the0.05, 0.25, and 0.5Bcr, respectively (corresponding
in the envelope. Therefore, simultaneous
imaging of the innera multiple system.has
Recentobservationssimilarof the submillimeter
nondimensional mass-to-ﬂux ratios of μ = 20,2). Thus,
et 4,al.and2000).
IRAS 04239+2436
a luminosity
to thatenvelope and the central compactonly
multiplesystemis alsodust emission polarization also show that the strong magnetic
magneticallysupercritical
clouds are considered here. For
of
the
Sun
and
very
rich
near-infrared(NIR)
spectrum(Greeneﬁeld in large ﬁlamentary structures is weakened down to the
an MHDscheme,adopted Boris-HLLD (Matsumoto et al.
required. One prominent structure thatshows
theweinteractions& evolution,
Lada 1996).
high–spectral
resolution
observation
supercriticalconditionNIR
at smallcore scales (Doi2019a), whichallowsisusthe
to follow
a long-termeven Amagneticallybetween forming binary and circumbinarymaterial
spiralin cases of strong magnetic ﬁelds. The ohmic
dissipation is revealed
et al. 2021).
The present observation providesgaseous
more strong
kinematicallya sub-astronomical-unit-scalearm structure (Takakuwa et al. 2017; Alves et al. 2019).
disk
with
a
KeplerianrotationaroundacentralmassofFigure 14. (a) Infall and rotation velocity distribution along the arms as a function of radius for the simulation. For comparison, black lines are shown for the freefall
10<10 auB
2

Initial conditions of star formation at .2000 au: physical structure
and NH3 depletion of three early-stage cores
Y. Lin1 , S. Spezzano1 , J. E. Pineda1 , J. Harju2 , A. Schmiedeke3 , S. Jiao4 , H. B. Liu5 , and P. Caselli1
1

Max-Planck-Institut für Extraterrestrische Physik, Giessenbachstr. 1, D-85748 Garching bei München, Germany
e-mail: ylin@mpe.mpg.de
2
Department of Physics, P.O. Box 64, FI-00014, University of Helsinki, Finland
3
Green Bank Observatory, PO Box 2, Green Bank, WV 24944, USA
4
National Astronomical Observatories, Chinese Academy of Sciences, Beijing 100101, China
5
Department of Physics, National Sun Yat-Sen University, No. 70, Lien-Hai Road, Kaohsiung City 80424, Taiwan, R.O.C.

ISM

rXiv:2308.12835v1 [astro-ph.GA] 24 Aug 2023

D
JVLA, K-band, NH3 1,1 & 2,2, dv=0.05 km/s. GBT

NH3Received ; accepted
Y. Lin et al.: Initial conditions of star formation at .2000 au: physical structure and NH3 depletion of three early-stage cores

10 57 40

8 13 20

L694-216

ABSTRACT 4.5

30 38 30

Herschel 250, 350, 500 µm
T .
JCMT 850 µm
ττ .
Spitzer 8 µm
N(H ).
𝑛𝑛 𝑟 =(𝑟⁄𝑟 )' +1, σ , T , T , N(NH ).

4.0 star formation. Characterisations of the physical condiL429
L1517B
Context. Pre-stellar cores represent
a criticalevolutionary phase in low-mass
14

Dec. (J2000)

NH3 1,1

NH

1,1

2,2

NH2 ( 1022 cm 2 )tions of pre-stellar cores provide important
constraints
159
pc to star and planet3.5formation theory, and are pre-requisite for establishing the20436
pc evolution and the12 related chemical
40
dynamical
processes.
3.0
dust
203 pc00
6
10 detailed thermal structure and density distribution of three early-stage cores, starless core L1517B, andAims. We aim to unveil the
2.5
00
prestellar core L694-2 and L429,
with the high angular resolution observations of the NH3 (1,1) and (2,2) inversion transitions obtained
8
14 00
850µm
with VLA and GBT. In addition, we explore where/if NH3 depletes in the2.0central regions of the cores.
4
37 30
6
1.5
Methods. We calculate the physical parameter maps of gas kinetic temperature,
NH3 column density, line-width and centroid velocity
56 40
8µm method to the Spitzer 8 µm map to20of the three cores, utilising4 the NH3 (1,1) and (2,2) lines. We apply the mid-infrared
extinction
1.0
2
obtainahighangularresolutionhydrogencolumndensitymap.
Weexaminethecorrelation
between the derived parameters, and
20
2
0.5
au
4000 au
00 derive the gas density 4000
23
au
the properties of4000
individual
cores. We
profile from the column density maps and assess the variation of NH
40
h m s
s
s
s
s
4
55
22
20
18
16
14
h m abundancessss a function
s
s
as
of
gas
volumedensity.
19h 41m 06s
04s
02s
18 17 08 07
06
05
04
03
Results.
The measured temperature profiles of thecores
L429 and L1517B show a minor decrease towards the core center, dropping (R.A.
(J2000)R.A. (J2000)
R.A. (J2000)
from ⇠9 K to below 8 K, and ⇠11 K to 10 K, while L694-2 has a rather uniform temperature distribution around ⇠9 K. Among thethree cores,
L429 has thecentraldensity, close
to sonic
line-width, and largest localised velocity gradient, allFig. 1. Hydrogen column density maps (N(H2 )) derived from Spitzer
8 µm extinction.
Thehighest
contours
show gas
the integrated
intensity
mapvelocityof the main
7&89
an advanced
stage. We resolvethat theabundanceofNH3 becomes 2 times lower in the central region ofcomponent of NH3 (1,1) line. The void regions in L694-2 and indicative
L1517B ofofvery
low N(Hevolutionary
brightpoint
8µmsources.
The2 ) correspond to nearby
4 1
3
1
3
LSR
V
kin
ex
3
L429, occurringa gas startdensityof14.4⇥10
ComparedOphiuchus/H-MM1 (Pineda et al. 2022) which shows an evenbeam size is indicated as a white filled circle in each figure. For L694-2
and L429aroundthe contours
from
K km
s cm
and. 1.5
K km s to, respectively,
5 1 3stronger
NHa3spacing
abundanceat 2⇥10with a spacing of 0.5 K km s 1 . For L1517B the contours start from
0.5 drop
K kmofs 1thewithof 0.3K km scm. , the abundance variations of the three cores plus Ophiuchus/H-MM1 suggest a
progressive NH3 depletion with increasing central density in pre-stellar cores.
8V

Key words. ISM: pre-stellar core – ISM: L429, L694-2, L1517B – ISM: structure – stars: formationand 500 µm emission maps4 from the Herschel/SPIRE tele- the interpolated background image. We therefore also adopt thescope5 (Griffin et al. 2010) to facilitate the derivation of the high small-scale median filter method (Butler & Tan 2009) to estiangular resolution, hydrogen column density maps (N(H2 )) from mate the emission background, which first e↵ectively smooths
1. Introduction
NH3 with for
deuterated ions in the gas phase (Rodgers & Charnley
the Spitzer/IRAC 8 µm data.
out the local inhomogeneity and then conducts interpolation
2001, Roue↵the masked region. The detailed procedure is elaborated
in Ap- et al. 2005, Sipilä et al. 2015), the process of which
With a modified black-body emission model, Pre-stellar
we first emcores represent a critical stage in the process of low- also a↵ects the abundance of NH3 , composing another source ofpendix D. At any rate, the process of background estimation isployed Herschel/SPIRE maps to derive a dust temperature
map
mass
star
formation:
the molecular gas has reached adequate
3 depletion. However, among various molecules, nitrogensubject to the difficulty of distinguishing betweenNH
small-scale(T dust ) at ⇠2100 resolution (with the high angularresolution
density for self-gravity to balance or even surpass the outward bearing species (e.g., N2 H+ , NH3 ) show relative longevity in the
backgroundemissionvariations
andthegenuineabsorbingcomSPIRE maps in the archive) from Spectral Energyforces
Distribution(thermal and turbulent pressure, rotation and magnetic gas (Caselli et al. 1999, Caselli et al. 2002b, Aikawa et al. 2005,
ponents (Butler & Tan 2009).
(SED) fitting. The temperature map is later used toderivethefield, see e.g. Myers & Benson 1983, Bergin & Tafalla 2007, Flower et al. 2006, Bergin & Tafalla 2007, Sipilä et al. 2015)
dust opacity map at 850 µm (⌧850 µm ). We adopted dust
emissivof its
theprior
foregroundlevel Ifg is guided by
Pinedaet al. 2023). Estimation
Compared tostage ofemission
starless cores,
andhave been important molecular tracers of the inner core reities for unprocessed (not coagulated) dust grains with
thin
ice
the
dust
opacitymapobtainedfromSED
fitting.
foregroundthese cores are denser and at the verge of forming protostars. TheThegion(e.g., inBergin et al. 2006; Friesen et al. 2009; Pineda et al.
mantles from Ossenkopf & Henning (1994) (OH4 model);inthis
to conditions
be spatiallyofconstantand is estimatedsignificance of emission
studying istheassumed
physicalpre-stellar
2010,2011; Chitsazzadeh et al. 2014; Pineda et al. 2015). Whilemodel, the dust emissivities in the submm regime can be approx- an iterative way such that the resultant peak ⌧aftersmoothing

L429 N(H2)~8e22 cm 2dN(NH3)/dN(H2)Initial conditions of star formation at .2000 au: physical structure
and NH3 depletion of three early-stage cores
Y. Lin1 , S. Spezzano1 , J. E. Pineda1 , J. Harju2 , A. Schmiedeke3 , S. Jiao4 , H. B. Liu5 , and P. Caselli1

rXiv:2308.12835v1 [astro-ph.GA] 24 Aug 2023

1

Max-Planck-Institut für Extraterrestrische Physik, Giessenbachstr. 1, D-85748 Garching bei München, Germany
e-mail: ylin@mpe.mpg.de
2
Department of Physics, P.O. Box 64, FI-00014, University of Helsinki, Finland
3
Green Bank Observatory, PO Box 2, Green Bank, WV 24944, USA
4
National Astronomical Observatories, Chinese Academy of Sciences, Beijing 100101, China
5
Department of Physics, National Sun Yat-Sen University, No. 70, Lien-Hai Road, Kaohsiung City 80424, Taiwan, R.O.C.

L694-2
ß
prototypicalReceived ; accepted

RADMC3D

L1544

Tkin

ABSTRACT

Context. Pre-stellar cores represent a critical evolutionary phase in low-mass star formation. Characterisations of the physical condi-

tions of pre-stellar cores provide important constraints to star and planet formation theory, and are pre-requisite for establishing the
dynamical evolution and the related chemical processes.
Aims. We aim to unveil the detailed thermal structure0and density distribution of three early-stage cores, starless core L1517B, andprestellar core L694-2 and L429, with the high angular resolution observations of the NH3 (1,1) and (2,2) inversion transitions obtainedwith VLA and GBT. In addition, we explore where/if NH3 depletes in the central regions of the cores.
Methods. We calculate the physical parameter maps of gas kinetic temperature, NH3 column density, line-width and centroid velocity
32
3
of the three cores, utilising the NH3 (1,1) and (2,2) lines. We apply the mid-infrared extinction method to the Spitzer 8 µm map to
3 map. We examine the correlation between the derived parameters, andobtain a high angular resolution hydrogen column density
2 the properties of individual cores. We derive the gas density profile from the column density maps and assess the variation of NH3
abundance
3as a function of gas volume density.
Results. The measured temperature profiles of the cores L429 and L1517B show a minor decrease towards the core center, droppingfrom ⇠9 K to below 8 K, and ⇠11 K to 10 K, while L694-2 has a rather uniform temperature distribution around ⇠9 K. Among thethree cores, L429 has the highest central gas density, close to sonic velocity line-width, and largest localised velocity gradient, allindicative of an advanced evolutionary stage. We resolve that the abundance of NH3 becomes 2 times lower in the central region of
L429, occurring around a gas density of 4.4⇥104 cm 3 . Compared to Ophiuchus/H-MM1 (Pineda et al. 2022) which shows an evenstronger drop of the NH3 abundance at 2⇥105 cm 3 , the abundance variations of the three cores plus Ophiuchus/H-MM1 suggest a
progressive NH3 depletion with increasing central density in pre-stellar cores.

G =1-2 (

ISRF)

L429
dN(NH )/dN(H )~1
NH
L429 dN(NH3)/dN(H )4e4 cmH-MM1
depletion
2e5 cm
(
)
Tkin

à

Tkin
<8 K à 9 K

Key words. ISM: pre-stellar core – ISM: L429, L694-2, L1517B – ISM: structure – stars: formation

1. Introduction

L429, L694-2

(?)

L429,
L694-2 L1517B
NH with deuterated ions in the gas phase (Rodgers & CharnleyPre-stellar cores represent a critical stage in the process of lowmass star formation: the molecular gas has reached adequate
density for self-gravity to balance or even surpass the outward
forces (thermal and turbulent pressure, rotation and magnetic
field, see e.g. Myers & Benson 1983, Bergin & Tafalla 2007,
Pineda et al. 2023). Compared to its prior stage of starless cores,
these cores are denser and at the verge of forming protostars. The
significance of studying the physical conditions of pre-stellar

3

2001, Roue↵ et al. 2005, Sipilä et al. 2015), the process of which
also a↵ects the abundance of NH3 , composing another source ofNH3 depletion. However, among various molecules, nitrogenbearing species (e.g., N2 H+ , NH3 ) show relative longevity in thegas (Caselli et al. 1999, Caselli et al. 2002b, Aikawa et al. 2005,
Flower et al. 2006, Bergin & Tafalla 2007, Sipilä et al. 2015)and have been important molecular tracers of the inner core region (e.g., Bergin et al. 2006; Friesen et al. 2009; Pineda et al.
2010, 2011; Chitsazzadeh et al. 2014; Pineda et al. 2015).
While
Fig. 11. Radial temperature profiles (solid lines) calculated by

F
v

7
Lupi

INAF-Osservatorio Astronomico di Capodimonte, via Moiariello 16, 80131 Napoli, Italy
INAF-Osservatorio Astronomico di Roma, Via di Frascati 33, 00078 Monte Porzio Catone, Italy
10
Brightness
and mass
accretionrate
evolution
duringthe 2022
SUPA, Schoolof Science and Engineering,
Universityof Dundee,
Nethergate,
Dundee DD1
4HN, UK
11?Marcos, TX 78666, USA
Department of Physics, Texas State University, 749
N Comanche
San
burst
of EXStreet,
Lupi
12
Max Planck Institute for Radioastronomy, Auf dem Hügel 69, 53121, Bonn, Germany
13
1,2,3 , Á. Kóspál
1,2,4,5 ,of
1,2,4 North
6 , C.
6 , J. Wendeborn
SUPA, Schoolof Physics deandMiera
Astronomy,
University
Andrews,
Andrews,
9SS, UK7 ,F. Cruz-Sáenz
P.StÁbrahám
, R. A. Haugh,
B. ClaesSt
F. ManaraKY16
9

E. Fiorellino8 , T. Giannini9 , B. Nisini9 , A. Sicilia-Aguilar10 , J. Campbell-White6 , J. M. Alcalá8 , A. Banzatti11 ,Received MayZs.
31,M.2023;
August4, ,2023
1,2,12,13 , F.
Szabóaccepted
Lykou1,2
S. Antoniucci9 , J. Varga1,2 , M. Siwak1,2 , S. Park1,2 , Zs. Nagy1,2 , and M. Kun1,2
1

Konkoly Observatory, Research Centre for AstronomyABSTRACT
and Earth Sciences, Eötvös Loránd Research Network (ELKH), KonkolyThege Miklós út 15–17, 1121 Budapest, Hungary
e-mail: fernando.cruz-saenz @ irap.omp.eu
Context. EX Lupi2 is
the MTA
prototypewhich EXor-typewere1121,
defined.
It hasexperienced multiple accretion-related bursts
CSFK,
Centre ofby
Excellence,
Konkoly Thegeoutbursts
Miklós út 15–17,
Budapest,
Hungary3
Institut deen Astrophysiqueet Planétologie,
Universitéde Toulouse,
OMP, CNRS,
av. e↵ects
du Colonel
and outbursts throughout
theRecherche
last decades,
whose study
have greatly
extended
ourUT3-PS,
knowledge
about 9theofRoche,these types of events.
31028 Toulouse Cedex 4, FranceThis star experienced
a
new
burst
in
2022.
4
ELTE Eötvös Loránd University, Institute of Physics, Pázmány Péter sétány 1/A, 1117 Budapest, Hungary
Aims. We aim to 5investigate
whetherrecentKönigstuhl
brightening
wasHeidelberg,
caused bytemporarily increased accretion or by a brief decrease in
Max Planck Institutefor this
Astronomy,
17, 69117
Germany
EuropeanObservatory,Karl-Schwarzschild-Strasse
2, 85748, Garching
bei München, Germany
the extinction, and67 to
studySouthern
the evolutionof theEX Lupi system throughout
this event.
Department of Astronomy & Institute for Astrophysical Research, Boston University, 725 Commonwealth Avenue, Boston, MA
Methods. We used02215,
multi-band
photometry to create color-color and color-magnitude diagrams to exclude the possibility that theUSA
8brightening could be
explained byAstronomico
a decreasein extinction.
We obtained
spectra
INAF-Osservatorio
di Capodimonte,
via Moiariello16, 80131
Napoli, using
Italy the X-shooter instrument of the Very Large
Via di Frascati 33, 00078 Monte Porzio Catone, Italy
Telescope (VLT)109toINAF-Osservatorio
determine theAstronomico
Lacc and diṀRoma,acc during the peak of the burst and after its return to quiescence using two di↵erent
SUPA, School of Science and Engineering, University of Dundee, Nethergate, Dundee DD1 4HN, UK
methods: empirical
relationships
betweenlineUniversity,
luminosity
Lacc , and
slab
model
theUSA
whole spectrum. We examined the 130
11
Departmentof Physics,
Texas State749 NandComanche
Street,aSan
Marcos,TX of78666,
12
year light curve of EXprovidestatistics on Auf
the demnumberof outbursts
during this period of time.
Max Lupi
Planck to
Institutefor Radioastronomy,
Hügel 69,
53121, Bonn,experienced
Germany
13
SUPA,
of Physics
Astronomy,
University
St Andrews, North
St Andrews,
KY16 9SS,
UK responsible for the brightenResults. Our analysisof School
the data
takenandduringthe 2022
burstofconfirmedthat Haugh,a changein extinction
is notrXiv:2308.02849v1 [astro-ph.SR] 5 Aug 2023

arXiv:2308.02849v1 [astro-ph.SR] 5

8ing. Our two approaches
in calculating
the AugustṀacc werein agreement, and resulted in values that are two orders of magnitude above what
Received May
31, 2023; accepted
4, 2023
had previously been estimated for EX Lupi using only a couple of individual emission lines, thus suggesting that EX Lupi is a strong
accretor even when in quiescence. We determined that in 2022ABSTRACT
March the Ṁacc increased by a factor of 7 with respect to the quiescentlevel. We also found hints that even though the Ṁacc had returned to almost its pre-outburst levels, certain physical properties of theContext. EX Lupi is the prototype by which EXor-type outbursts were defined. It has experienced multiple accretion-related bursts
gas (i.e. temperatureand density) had not returned to the quiescent values.
and outbursts throughout the last decades, whose study have greatly extended our knowledge about the e↵ects of these types of events.
Conclusions. WeThis
foundthat the amass
accretedstar experienced
new burst
in 2022.during this three month event was 0.8 lunar masses, which is approximately half of
We aaim
to investigate
whether thisrecentbrighteningthat
was ifcausedtemporarily
increasedaccretiona brief
decreasein past 130 years,
what is accreted Aims.
duringyearof quiescence.
We
calculated
EXbyLupi
remains
as active
as orit byhas
been
for thethe extinction, and to study the evolution of the EX Lupi system throughout this event.
during which it has
experienced
at
least
3
outbursts
and10
bursts,
then
it
will
deplete
the
mass
of
its
circumstellar
material in less
Methods. We used multi-band photometry to create color-color and color-magnitude diagrams to exclude the possibility that thethan 160 000 yr. brightening could be explained by a decrease in extinction. We obtained spectra using the X-shooter instrument of the Very Large
Telescope (VLT) to determine the Lacc and Ṁacc during the peak of the burst and after its return to quiescence using two di↵erent

methods:
empirical
relationships
line luminosityand Lacc ,Herbig
and a slabAe/Be
model ofwhole spectrum.
We examined
130Key
pre-mainsequence– betweenstars: variables:
T Tauri,– theaccretion,accretion
disks –thetechniques:
spectroEXwords.
Lupistars:
2022
year light curve
EX Lupi to provide statistics on the number of outbursts experienced during this period of time.
scopic – stars: individual:
EXofLupi

3

Results. Our analysis of the data taken during the 2022 burst confirmed that a change in extinction is not responsible for the brighten-

ing. Our two approaches in calculating the Ṁacc were in agreement, and resulted in values that are two orders of magnitude above what
had previously been estimated for EX Lupi using only a couple of individual emission lines, thus suggesting that EX Lupi is a strong
accretor even when in quiescence. We determined that in 2022 March the Ṁacc increased by a factor of 7 with respect to the quiescentalmost its pre-outburst
levels,certainphysicalpropertiesof theand references). This
1. Introduction level. We also found hints that even though the Ṁacc had returned tosmaller
bursts (see
Fig.
1 for
the light
curvegas (i.e. temperature and density) had not returned to the quiescent values.
stareventhaswasbeen
taken
to whichbe theprototypehalfofoferuptive young stars
Conclusions. We found that the mass accreted during this three month
0.8 lunar
masses,is approximately
Over the past 130 years,
Lupi
hasa year
experienced
threewhat isEX
accreted
duringof quiescence.at
Weleastcalculated
that if EX
remains as active
as it has been
the magnitudes
past 130 years, lasting from a few
thatLupiexperience
outbursts
of for
2–5
duringwhichandit has2008)
experienced
least 3 outbursts
10 bursts, then it will deplete the mass of its circumstellar material in lesslarge outbursts (1944,
1955
andat morethan aanddozen

⾒160,000

0.8

EX Lupimonths up to a year (EXors, Herbig 1989). Together with theirEX4

Hungary
CSFK, MTA Centre of Excellence, Konkoly-Thege Miklós út 15–17, 1121 Budapest, Hungary

Forming rocky exoplanets around K-dwarf stars

Received 6 March 2023 / Accepted 13 June 2023

P. Hatalova1,2 , R. Brasser3,4 , E. Mamonova1,2 , and S. C. Werner1,2
1

ABSTRACTCentre for Earth Evolution and Dynamics (CEED), University of Oslo, PO Box 1028 Blindern, 0315 Oslo, Norway
e-mail: petra.hatalova@geo.uio.no2
Centrefor Planetary
Habitability
Universityof Oslo,
PO Box 1028
0315andOslo,study
NorwayContext. New spacetelescopes,
such
as the (PHAB),
upcoming
PLATO
mission,aimBlindern,
to detectthousands of exoplanets, especially3
Origins
Research
Institute,
Research
Centrefor
Astronomy
andEarth
Sciences,
Konkoly-Thege
Miklósút 15–17,1121multiple
Budapest, close-in superterrestrial planets around main-sequence stars. This motivates us to study how these planets formed.
How
Hungary
Earths form around4 CSFK,
stars MTA
withCentre
masseslower thanthat of theMiklós
Sun útis 15–17,
still anopen
issue.
Several recent modeling studies have focused onof Excellence,
Konkoly-Thege1121
Budapest,
Hungary

planet formation around M-dwarf stars, but so far no studies have focused specifically on K dwarfs, which are of particular interest in
Received 6 March 2023 / Accepted 13 June 2023
the search for extraterrestrial
life.
Aims. We aim to reproduce the currently known population of close-in super-Earths observed around K-dwarf stars and their system
ABSTRACTcharacteristics. Additionally, we investigate whether the planetarysystems that we form allow us to decide which initial conditions are
the most favorable.
Context. New space telescopes, such as the upcoming PLATO mission, aim to detect and study thousands of exoplanets, especiallyterrestrial planets
around main-sequencestars. This
motivates usof
to study
howformationthese planetsvia
formed.
How multipleaccretion
close-in superMethods. We performed
48 high-resolution
N-bodysimulations
planetplanetesimalusing the existing
Earths form around stars with masses lower than that of the Sun is still an open issue. Several recent modeling studies have focused on
GENGA softwareplanet
running
on
GPUs.
In
the
simulations
wevariedtheinitialprotoplanetary
disk
mass
and
the
solid
formation around M-dwarf stars, but so far no studies have focused specifically on K dwarfs, which are of particular interest in and gas surface
density profiles. Eachsimulation
beganlife.
with 12 000 bodies with radii of between 200 and 2000 km around two different stars, withthe searchfor extraterrestrial
Aims.
We
aim
to
reproducethecurrently
known
of close-in
observedaround K-dwarf
andMyr.
their system
masses of 0.6 and 0.8 M . Most simulations
ran for20 population
Myr, withseveralsuper-Earths
simulations
extended
to 40 orstars
100
characteristics. Additionally, we investigate whether the planetary systems that we form allow us to decide which initial conditions are
Results. The mass
the distributions
most favorable. for the planets with masses between 2 and 12 M show a strong preference for planets with masses
Mp < 6 M and aMethods.
lesser preference
planets withN-bodylargersimulations
masses, ofwhereas
the mass
distribution
for theobserved
sample increases
We performedfor
48 high-resolution
planet formation
via planetesimalaccretionusing
the existing
GENGA
softwarerunning
on
GPUs.
In
the
simulations
wevariedtheinitialprotoplanetary
disk
mass
and
the
solid
and
gas
surfacealmost linearly. However, we managed to reproduce the main characteristics and architectures of the known planetary systems anddensity profiles. Each simulation began with 12 000 bodies with radii of between 200 and 2000 km around two different stars, withproduce mostly long-term
systems,
but wean initial disk mass of 15 M or
masses of 0.6angular-momentum-deficit-stable,
and 0.8 M . Most simulations ran for 20 Myr,nonresonantwith several simulations
extended
to 40requireor 100 Myr.
higher and a gas surface
densityvalue at 1forAU
of 1500
cm-2 orhigher.
producemany
low-massplanets with
Results. The
mass distributions
the planets
withgmasses
between2 andOur
12 Msimulations
show a strongalso
preference
for planetswithmasses
M
<
6
M
and
a
lesserpreferenceforplanets
withlargermasses,whereasthe
mass
distributionfortheobserved
sampleincreasesp
M < 2 M , which are not yet found in the observed population, probably due to the observational biases. Earth-mass planets form
almost linearly. However, we managed to reproduce the main characteristics and architectures of the known planetary systems andquickly (usually within
few long-term
million angular-momentum-deficit-stable,
years), mostly before the gas
disk dispersal.
final ansystems
containsmall number of
produce amostly
nonresonant
systems, butThe
we require
initial disk
mass ofonly
15 M a or
highergasM
surfacedensitycould
value at
1 AU of 1500
g cm-2substantial
or higher. Ouramounts
simulationsofalso
producemany low-massplanets with masses
Mpand> a10, whichpossiblyaccretegas,
and these
formedplanets
after withthe gas had mostly
M < 2 M , which are not yet found in the observed population, probably due to the observational biases. Earth-mass planets form
dissipated.
quickly (usually within a few million years), mostly before the gas disk dispersal. The final systems contain only a small number ofConclusions. Weplanets
mostlyreproduce observed K-dwarf exoplanetary systems from our GPU simulations.
withmanage
masses Mto
p > 10 M , which could possibly accrete substantial amounts of gas, and these formed after the gas had mostly
Key K
words. planetsand satellites:– planets
and satellites:
dynamical
evolutionand simulations.
stability – planet-disk interactions –Conclusions.
We mostlyformation
manage to reproduce
observed
K-dwarf exoplanetarysystemsfrom our GPU
dissipated.

2
1

protoplanetary disks
Key words. planets and satellites: formation – planets and satellites: dynamical evolution and stability – planet-disk interactions –48protoplanetary disks

1. Introduction1. Introduction

1500g cm-2

of February
5200 exoplanets
been
As of February As
2023,
more 2023,
than more
5200thanexoplanets
havehavebeen
1
confirmed
from stellar observations, and several thousand plan1
2
M_E
confirmed from stellar
observations, andseveralsuch
thousand
etary candidates
discoveredby missions
as Kepler,planK2, andetary candidates discoveredby missions
suchas
Kepler,K2, planets
andTESS are awaiting
confirmation.
Themajorityof these10
M_E
are confirmation.
so-called super-Earths
mini- or sub-Neptunes,
which fall
TESS are awaiting
Theandmajorityof these planets
between Earth and Neptune in radius and/or mass; such planets

6M_E

0.6

0.8

15M_E

space
telescopes,
such
as the upcoming
PLATO mission,
space telescopes,
such asthe upcoming
ESA PLATOESA
mission,aim totodetectand andstudystudy
thousandsof exoplanets.
The ESA-statedThe ESA-statedaim
detectthousandsof exoplanets.
primary
goal
of
the
PLATO
missionisthe“detection
and“detection
characprimary
goal of the PLATO mission is theand characterization of terrestrial exoplanets around bright solar-type stars,
terizationofterrestrial
exoplanets
aroundbrightsolar-typestars,
with emphasis on planets orbiting in the habitable zone”. F,
G, andemphasis
K dwarfs areon
considered
be “solar-like
This
with
planetstoorbiting
in stars”.
the habitablezone”. F,plannedamongareothers,will provide exciting
opporG,
andmission,
K dwarfs
considered
to be new“solar-likestars”. This⾒

9 Center for Astrophysics and Space Astronomy, Astrophysical and Planetary Sciences Department, University of Colorado, UCB 389 Boulder, Colorado

The VLT MUSE NFM viewof outflows and externally photoevaporating
80309, USA⋆
10 Instituto de Fı́sica Fundamental (CSIC). Calle Serrano 121-123,
28006 Madrid, Spain
discs near the orion bar

1 Astronomy Unit, School of Physics and Astronomy, Queen Mary University of London, London E1 4NS, UK

ABSTRACT

2 Department of Physics and Astronomy, Rice University, 6100 Main St - MS 108, Houston, TX 77005, USA
3 Department of Physics and Astronomy, Vanderbilt University, Nashville, TN 37235-1807, USA
4

Downloaded from https://academic.oup.com/mnras/article/525/3/4129/7251483 by Dankook University / Library user on 10 O

Thomas J. Haworth ,1 ‹ † Megan Reiter ,2 C. Robert O’Dell,3 Peter Zeidler,4 Olivier Berne,5Accepted 2023 August
23. Received
2023,6August
in original,7form
2023 July
17 8 John Bally,9 Javier R. Goicoechea,10Carlo
F. Manara
Giulia23;
Ballabio
Jinyoung
S. Kim,6
6
6
Mari-Liis Aru, Aashish Gupta and Anna Miotello

AURA for the European Space Agency (ESA), ESA Office, Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA
We present Very 5Large
Telescope/Multi-Unit Spectroscopic Explorer Narrow Field Mode observations of a pair of disc-bearing
Institut de Recherche en Astrophysique et Planétologie, Université de Toulouse, CNRS, CNES, UPS, 9 Av. du colonel Roche, 31028 Toulouse Cedex 04, France
6 European Southern Observatory, Karl-Schwarzschild-Strasse 2, 85748 Garching bei München, Germany
young stellar objects
towards the Orion Bar: 203–504 and 203–506. Both of these discs are subject to external photoevaporation,
7 Astrophysics Group, Imperial College London, Blackett Laboratory, Prince Consort Road, London SW7 2AZ, UK
8where winds are Steward
launchedfromtheirof Arizona,outer 933
regions
dueTucson,
to environmental
Observatory,
University
N. Cherry Ave,
AZ 85721-0065, USA irradiation. Intriguingly, despite having projected
9 Center for Astrophysics and Space Astronomy, Astrophysical and Planetary Sciences Department, University of Colorado, UCB 389 Boulder, Colorado
separation from one
another of only 1.65 arcsec(660 au at 400 pc), 203–504 has a classic teardrop shaped ‘proplyd’ morphology
80309, USA
2
10 Instituto
Fundamental (CSIC).
Calle Serrano
28006 Madrid,
Spainstar, rather than θ 1 Ori C) but 203–506 has no ionization
pointing towards θ OrideAFı́sica(indicating
irradiationby121-123,
the EUVof thatfront, indicating it is not irradiated by stellar EUV at all. However, 203–506 does show [C I] 8727 Å and [O I] 6300 Å in emission,indicating irradiation
by2023stellar
FUV.
This
demonstrates
the importance of FUV irradiation in driving mass loss fromAccepted
August 23.
Received
2023 explicitly
August 23; in originalform 2023 July 17
discs. We conclude that shielding of 203–506 from EUV is most likely due to its position on the observers side of an ionized layer
lying in the foregroundA B S T R Aof
C Tthe Huygens Region. We demonstrate that the outflow HH 519, previously thought to be emanatingfrom 203–504 is We
actually
cloud edge andidentify Explorer
a new compact
outflowfrom thatofobject
along ourpresent an
VeryirradiatedLarge Telescope/Multi-Unit
Spectroscopic
Narrow Field
Mode observations
a pair ofapproximately
disc-bearing−1
youngstellarobjects
towards
the
Orion
Bar:
203–504and
203–506.
Both
of
these
discs
are
subjecttoexternal
photoev
aporation,line of sight with a velocity ∼130 km s .

where winds are launched from their outer regions due to environmental irradiation. Intriguingly, despite having projected
separation
one anotherof only 1.65–arcsec(660
au at 400 pc),
203–504jets
has –a classic
shaped–‘proplyd’morphology
Key words: planets
andfromsatellites:
formationstars: formation– stars:
stars:teardrop
protostars
ISM: HII
regions – ISM: kinepointing towards θ 2 Ori A (indicating irradiation by the EUV of that star, rather than θ 1 Ori C) but 203–506 has no ionization
matics and dynamics.
front, indicating it is not irradiated by stellar EUV at all. However, 203–506 does show [C I] 8727 Å and [O I] 6300 Å in emission,indicating irradiation by stellar FUV. This explicitly demonstrates the importance of FUV irradiation in driving mass loss fromdiscs. We conclude that shielding of 203–506 from EUV is most likely due to its position on the observers side of an ionized layer
insight into the gas and dust composition and conditions, revealing
1 I N T RO D U C lying
T I OinNthe foreground of the Huygens Region. We demonstrate that the outflow HH 519, previously thought to be emanatingfrom 203–504 is actually an irradiated cloud edge and identify a newcompact outflow
frommaythat object
approximately
along ourwith the planet
sub-structuresthatpossiblybe associated
−1
line
of
sight
with
a
velocity
∼130
km
s
.
With ever increasing numbers of exoplanets being discovered (Gaudi,
formation process, or planets themselves (e.g. ALMA Partnership2

203-504[C I] 8727 Å[O I] 6300 Å

θ2 Ori A
203-506

203-504
660 auEUV
EUV

203-506//academic.oup.com/mnras/article/525/3/4129/7251483 by Dankook University / Library user on 10 October 2023

8 Steward Observatory, University of Arizona, 933 N. Cherry Ave, Tucson, AZ 85721-0065, USA

FUV

Meyer & Christiansen), therea needto understandhow formation
theDipierroal. 2015; Andrews
et al.
2018kine; Zhang et al. 2018;
Key 2021
words:
planetsisandsatellites:
formation – stars:–2015
stars:; jets
– stars: etprotostars– ISM:HII regions– ISM:
FUV
203-506
EUV
matics anddynamics. to exoplanet diversity. Thereplanet formation processcontributes
Pinte et al. 2020; Öberg et al. 2021). While those observations of
is now overwhelming evidence that planets form from circumstellar
nearby bright discs are extremely important, recently there is a shift‘planet-forming’ discs, to the extent that we can even directly observe
towards
tryingstudy
more typicalinsight into
the gastoand
dust compositionand planet-forming
conditions, revealingdiscs. For nearby
1 I N T RO D U C T I O N
HH 519130km/s
203-504
sub-structuresthatmaypossiblybeassociatedwith ofthemoreplanettypical radii (R ≤young planets within them in some cases (Keppler et al. 2018; Pinte
regions, this entails studying in detail thoseWith ever increasing numbers of exoplanets being discovered (Gaudi,
formation process, or planets themselves (e.g. ALMA Partnershipet al. 2018; TeagueMeyer
et al. &2018; Haffert et al. 2019; Pinte et al. 2022).
au, Zhang et al. 2023; Miotello et al. 2021, Miley et al. submitted,
50
Christiansen 2021), there is a need to understand how the2015; Dipierro et al. 2015; Andrews et al. 2018; Zhang et al. 2018;
Our most detailed studiesof theseprocess
have contributes
typicallytobeen
targeting
large
and
DECO
ALMA
program:2022.1.00875.L,
PI: Cleeves).
planet formation
exoplanet
diversity.
There
Pintetheet al.
2020; Öberg
et al.large2021).
While thoseobservations ofTwenty-Five Years of Accretion onto TW Hya
Twenty-Five Years of Accretion onto the Classical T Tauri Star TW Hya

Herczeg, Chen, et al.

arXiv:2308.14590v1 [astro-ph.SR] 28 Aug 2023

41 Konkoly Observatory, Research Centre for Astronomy
45
7
and Earth
Sciences,
Eötvös Loránd(ELKH),
HungarianAcademy of Sciences,
Gregory J. Herczeg (à˜L),1, 2, 3 Yuguang
Chen (H1I),
Jean-Francois
Donati,Research
Andrea K. Network
Dupree,6 Frederick
M. Walter,
8
9
10
11
12, 13
Miklós útCarlo
15–17,
1121 Hans
Budapest,
HungaryMin Fang (?O),
Lynne A. Hillenbrand,Konkoly-Thege
Christopher M. Johns-Krull,
F. Manara,Moritz Günther,P. Christian Schneider,14 Jeff A. Valenti,15 Silvia H.P. Alencar,16 Laura Venuti,17 Juan Manuel Alcalá,18 Antonio Frasca,1942 CSFK,MTA Centreof Excellence,
Budapest,
Miklós15-17,
H-1121,
Hungary
15
21
6
22
23
Nicole Arulanantham,
Jeffrey
L. Linsky,20 Jerome
Bouvier,Konkoly
Nancy S. Thege
Brickhouse,
NuriaútCalvet,Catherine
C. Espaillat,
1024
25
26, 27, 282930
31
43 Department
Justyn Campbell-White,
John M.
Carpenter,Seok-JunSchool
Chang, of
Kelle
L. Cruz,OsakaS.E.
Dahm, Jochen
Eislöffel, Osaka
Suzan Edwards,of Earth andSpace
Science,
Graduate
Science,
University,
Toyonaka,
560-0043,
Japan
15
32, 33, 34, 35
36
1, 2
37
38
36
William J. Fischer, Zhen Guo (Ì⌥),
Thomas Henning, Tao Ji (™õ), Jessy Jose, Joel H. Kastner, Ralf Launhardt,
David A. Principe,11 Connor E. Robinson,39 Javier Serna,40 Michał Siwak,41, 42 Michael F. Sterzik,10 and Shinsuke Takasao43

Abstract

1 Kavli Institute for Astronomy and Astrophysics, Peking University, Beijing 100871, China
2 Department of Astronomy, Peking University, Beijing 100871, China

3 Visiting
Accretion plays a centralrole
in thephysics
that governs
theofevolution
andCAdispersal
astronomer,
Departmentof Astronomy;
California InstituteTechnology; Pasadena,
91125, USA of protoplanetary disks.
4 Department of Physics & Astronomy, University of California, Davis, 1 Sheilds Avenue, Davis, 95616, CA, USA
The primary goal of this paper5 IRAP-UMR
is to analyze
the stability over time of the mass accretion rate onto TW Hya, the5277, CNRS & Univ. de Toulouse, 14 Av. E. Belin, F-31400, Toulouse, France
6 Centerfor Astrophysics| Harvard &Smithsonian
nearest accreting solar-mass young star. We
measureveilingacrossthe optical spectrum in 1169 archival high7 Department of Physics & Astronomy, Stony Brook University, Stony Brook NY 11794-3800, USA
resolution spectra of TW Hya,8 Department
obtainedfrom 1998–2022. The veiling is then converted to accretion rate usingof Astronomy; California Institute of Technology; Pasadena, CA 91125, USA
9 Department of Physics and Astronomy, Rice University, 6100 Main Street, Houston, TX 77005, USA
26 flux-calibrated spectra 10that cover the Balmer jump. The accretion rate measured from the excess continuum
European Southern Observatory,Karl-Schwarzschild-Strasse 2, 85748 Garching bei München, Germany
1
has an average of 2.51 ⇥10Kavli9Institute
M foryrAstrophysics
andanda Space
Gaussiandistribution
witha FWHM
11 MIT
Research, 77
Massachusetts Avenue,
Cambridge,
MA 02139,of
USA0.22 dex. This accretion
12 Purple Mountain Observatory, Chinese Academy of Sciences, 10 Yuanhua Road, Nanjing 210023, China
rate may be underestimated by a factor
of up to 1.5 because of uncertainty in the bolometric correction and
13 University of Science and Technology of China, Hefei 230026, China
14 Hamburger Sternwarte,
another factor of 1.7 because of excluding
the fractionof accretion
energythat escapes in lines, especially Ly↵.
Gojenbergsweg
112, 21029, Hamburg,
Germany
15 Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA
The accretion luminosities
are well correlated with He line luminosities but poorly correlated with H↵ and H
16 Departamento de Fisica - ICEx - UFMG, Av. Antonio Carlos 6627, 30270-901 Belo Horizonte, MG, Brazil
SETI Institute,
339 Bernardo Ave.,
Suitehours
200, Mountain
94043, USA
luminosity. The accretion rate is 17always
flickeringoverbutView,
onCAlonger
timescales has been stable over 25
18 INAF-Osservatorio Astronomico di Capodimonte, via Moiariello 16, 80131 Napoli, Italy
years. This level of variability is19consistent
with previous measurements for most, but not all, accreting young
INAF - Osservatorio Astrofisico di Catania, via S. Sofia, 78, 95123 Catania, Italy
20 JILA, University of Colorado, Boulder CO 80309-0440
stars.
21

Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France

22 Department of Astronomy, University of Michigan, 1085 South University Avenue, Ann Arbor, MI 48109, USA

TW Hya23 Institute for Astrophysical Research, Department of Astronomy, Boston University, 725 Commonwealth Avenue, Boston, MA 02215, USA
24 Joint ALMA Observatory, Avenida Alonso de Córdova 3107, Vitacura, Santiago, Chile

25 Max-Planck-Institut
Unified Astronomy Thesaurus
concepts:
Classical
T Tauri(252);
Protoplanetary(1300);
fürAstrophysik, Karl-Schwarzschild-Strasse
1, 85748 Garchingb. München
2.51 ×10^(-9)M starsyr^(-1)
FWHM disks
0.22
dexStellar ac26 Dept of Physics and Astronomy, Hunter College, City University of New York, NY, NY, USAcretion disks (1579); Variable stars27 (1761); High resolution spectroscopy (2096)

1.5

Physics, Graduate Center of the City University of New York, NY, NY, USALyα
1.7
Avenue,AZ, 85719, USA
HeTucson,measurements
1. Introduction Gemini Observatory/NSF’s NOIRLab, 950 N. Cherry
continuum(e.g. Natta et al. 2004; Fang et Hαal.
Thüringer Landessternwarte, Sternwarte 5, D-07778 Tautenburg, Germany
HβThe evolution of protoplanetary disksFiveand
2009;
Alcalá
et25
al.USA2017).
College Astronomy Department, Smith College,
Northampton,
MA 01063,
the final out28 Department of Astrophysics, American Museum of Natural History, Ny, NY, USA

29

30

31

32

Instituto de Física y Astronomía, Universidad de Valparaíso, ave. Gran Bretaña, 1111, Casilla 5030, Valparaíso, Chile
Each1111,
ofCasilla
the5030,
observablediagnostics of accretion rate has
33 Núcleo Milenio de Formación Planetaria (NPF), Chile, ave. Gran Bretaña,
Valparaíso, Chile
come of any planet formation depends
in part on how gas
34 Centre for Astrophysics Research, University of Hertfordshire, Hatfield AL10 9AB, UK
its own advantages and disadvantages. Photometry covers
flows through the disk (see review
by Manara
et al. 2022).
35 Departamentode Física, Universidad Tecnicá Federico Santa María, Avenida España 1680, Valparaíso, Chilelarge samples on short and long timescales, though changes
36E-mail: 223s415s@gsuite.kobe-u.ac.jp

Tenri JINNO,1 ,∗ Takayuki R. SAITOH,1 Yota ISHIGAKI,2 and Junichiro MAKINO1

Abstract1

Department of Planetology, Graduate School of Science, Kobe University, 1-1 Rokkodai-cho, Nada-ku,
Kobe, Hyogo 657-8501, Japan
2
National Astronomical
Observatoryof Japan,
2-21-1 Osawa,
Tokyo 181-8588,that
Japan protoplanetary disks
In the conventional
theory ofplanetformation,
itMitaka,
is assumed∗
E-mail: and223s415s@gsuite.kobe-u.ac.jp
are axisymmetrichave a smooth radial profile. However, recent radio observations
Received 2023 May 31; Accepted 2023 August 2
of protoplanetary disks have revealed that many of them have complex radial structures.
In this study, we
perform a series of N-body simulations to investigate how planets are
Abstract
In the conventional
theorywithof planet
formation,
it is assumed
protoplanetary
formed in protoplanetarydisks
radialstructures.
Forthatthis
purpose,disks
we consider theare axisymmetric and have a smooth radial profile. However, recent radio observations
effect of continuous
pebble accretion onto the discontinuity boundary within the terresof protoplanetary disks have revealed that many of them have complex radial structures.
In this study,
we perform
a series
N-bodysimulations
to investigate howgrow
planetsefficiently
aretrial planet-forming
region(∼0.6au).ofWe
findthat protoplanetsat theformed in protoplanetary disks with radial structures. For this purpose,
4 we consider the
discontinuity boundary,reaching the Earth mass within ∼10 yr. We confirm that giant
effect of continuous pebble accretion onto the discontinuity boundary within the terrescollisions of protoplanets
occur
uni
verau).
sally
Moreover,weatfind
trial planet-forming
region(∼0.6
We infindourthat model.
protoplanetsgrow efficiently
the that multiple
4
discontinuity
boundary,
reaching
the Earth mass
within∼10 yr.of
Wethe
confirmthat giant
planet-sized bodies
form
at regularintervals
in thevicinity
discontinuity
boundary.
collisions of protoplanets occur universally in our model. Moreover, we find that multiple
These results planet-sizedindicatebodies
the form
possibilityof thein theformation
solar system-likeplanetaryat regular intervals
vicinity of theof
discontinuity
boundary.
These results
indicate protoplanetary
the possibility of the disks.
formation of solar system-like planetarysystems in radiallystructured

systems in radially structured protoplanetary disks.

Key words: methods:
numerical
planetsand satellites:
formation— planets
andterrestrial
satellites:
Key words:
methods:—
numerical— planets
and satellites: formation— planets
and satellites:
planetsterrestrial planets
1 Introduction

N(~0.6 au)

Downloaded from https://academic.oup.com/pasj/article/75/5/951/7262729 by Dankook University / Library user on 10 October 2023

through pebble accretion in a radially structured
protoplanetary disk
Received 2023 May 31; Accepted 2023 August 2∗

when the mass of a protoplanet reaches a certain criticalvalue (!10 M⊕ ), disk gas flows into the planet, forming
One of the pioneering studies of planetary formation theory
when the mass of a protoplanet reaches a certain criticala gas giant (Mizuno 1980). Since this model successfully
was conducted by Chushiro Hayashi’s group at Kyoto Univalue(!10 of
Mterrestrial
gas gasflows⊕ ), disk
described
the sequenceplanets,
giants,into
and the planet, forming
One of the pioneeringstudies1981of planetary
formation
theoryversity (Hayashi; Hayashi et al.
1985). In theirsceicy planets,
it
came
to
be
regardedasthestandard
theorya gas giant (Mizuno 1980). Since this model successfully
nario, the planetary formation process starts from a disk
was conducted by
Chushiro Hayashi’s group at Kyoto Uni- of planet formation, even though many difficulties werethat is stable against self-gravity. This disk is assumed todescribed the sequence of terrestrial planets, gas giants, andknown.
versity (Hayashi have
1981
; Hayashi
et al. 1985). current
In theira dustcomponent consistent
with themass sceofSomethe most criticalproblems
haveregarded
been resolved.
icy ofplanets,
it cameto beas the standard theory
planets,
and the dust-gas
ratio comesfrom thesolar aabunnario, the planetary
formation
process
startsfromdisk For example, the problem of the formation time of Jupiter

1 Introduction 10,000

https://academic.oup.com/pasj/article/75/5/951/7262729 by Dankook University / Library user on 10 Octob

Department of Planetology, Graduate School of Science, Kobe University, 1-1 Rokkodai-cho, Nada-ku,
Kobe, Hyogo 657-8501, Japan
2
N-bodysimulationof2-21-1
planetary
formation
National Astronomical
Observatoryof Japan,
Osawa, Mitaka,
Tokyo 181-8588, Japan

1

Institut für Astronomie und Astrophysik, Universität Tübingen, Auf der Morgenstelle 10, D-72076, Germany
Astronomy Unit, Department of Physics and Astronomy, Queen Mary University of London, London E1 4NS,UK
numerical
protoplanetary3
Leibniz-InstitutMeasuring
für Astrophysikthe
Potsdam(AIP), Anviscosity
der Sternwartein16,simulations
14482 Potsdam, ofGermany.
disks in Cartesian grids g.a.turpin, r.p.nelson}@qmul.ac.uk
e-mail: jjoseph@aip.de,lucas.jordan@uni-tuebingen.de,{a.ziampras,
Accepted date: 24 July 2023

The viscously spreading ring revisited

Jibin Joseph1, 3 , Alexandros Ziampras1, 2 , Lucas Jordan1 , George A. Turpin2 , Richard P. Nelson2
1

ABSTRACT

Institut für Astronomie und Astrophysik, Universität Tübingen, Auf der Morgenstelle 10, D-72076, Germany
2
Astronomy Unit, Department of Physics and Astronomy, Queen Mary University of London, London E1 4NS,UK
3
Context. Hydrodynamical
simulations
solve Potsdam
the governing
equationson14482a discrete
grid of space and time. This discretization causes
Leibniz-Institut
für Astrophysik(AIP), An derSternwarte 16,
Potsdam, Germany.
e-mail:jjoseph@aip.de,lucas.jordan@uni-tuebingen.de,{a.ziampras,
g.a.turpin,
r.p.nelson}@qmul.ac.uk
numerical di↵usion similar to a physical viscous di↵usion, whose magnitude is often unknown
or poorly constrained. With the current

rXiv:2308.03881v1 [astro-ph.IM] 7 Aug 2023

308.03881v1 [astro-ph.IM] 7 Aug 2023

2

trend of simulating accretion
disks
with
Accepted date:24 July2023no or very low prescribed physical viscosity, it becomes essential to understand and quantify
this inherent numerical di↵usion, in the form of a numerical viscosity.
ABSTRACT
Aims. We study the behavior of the viscous spreading ring andthe spiral instability that develops in it. We then use this setup toquantify the numericalviscosity
in Cartesian
grids
study its
properties.
Context.
Hydrodynamical
simulations
solveandthe governing
equationson a discrete grid of space and time. This discretization causes
numerical
di↵usion
similartoaphysicalviscousdi↵usion,
whose
is oftenor poorly constrained.
With thecurrent
Methods. We simulate the viscous spreading ring and the related magnitude
instability
onunknown
a two-dimensional
polar
grid
using PLUTO as well as
trend of simulating accretion disks with no or very low prescribed physical viscosity, it becomes essential to understand and quantify
FARGO, and ensure thethisconvergenceof
ourresultswitharesolution
study.
WethenrepeatourmodelsonaCartesiangrid and measure
inherent numerical di↵usion, in the form of a numerical viscosity.
Aims.
study the behaviorof thetoviscous
spreading analyticalring and the spiralinstability
that develops
it. WeAthena++.
then use this setup to
the numerical viscosity
byWecomparing
results
the knownsolution,using
PLUTOin andquantify the numerical viscosity in Cartesian grids and study its properties.
2
Results. We find thatMethods.
the numericalin a Cartesianscaleswith
aspolarapproximately
We simulateviscosity
the viscous spreading
ring and thegrid
relatedinstability
on aresolution
two-dimensional
grid using PLUTO⌫asnum
well/as x and is equivalent
4
FARGO,forandaensure
the convergenceof our results
withWea resolution
study. We
thentherepeatour models
on a Cartesian
grid and measure
to an e↵ective ↵ ⇠ 10
commonnumerical
setup.
also showthatspiralinstability
manifests
as a single leading spiral
the numerical viscosity by comparing results to the known analytical solution, using PLUTO and Athena++.
throughout the wholeResults.
domainon polar grids. This is contrary to previous results and indicates that sufficient resolution is necessary
We find that the numerical viscosity in a Cartesian grid scales with resolution as approximately ⌫num / x2 and is equivalent
in order to correctly resolvethe↵instability.
to an e↵ective⇠ 10 4 for a common numerical setup. We also show that the spiral instability manifests as a single leading spiralthroughout
the
whole
domaingrids. This
contrary towhereprevious results
and indicatesthat sufficient
resolutionin
is necessary
Conclusions. Our results
are relevant
in on
thepolarcontext
ofismodels
the origin
shouldbe included
the computational domain,in order to correctly resolve the instability.
or when polar grids cannot
be Ourused.
Examplesofthesuchcases
includemodelsof diskaccretion
onto a central
Conclusions.
resultsare relevant incontext
of models
where theorigin should
be includedin the computational
domain,binary and inherently
or
when
polar
grids
cannotbeused.
Examples
of
such
cases
includemodelsofdiskaccretionontoacentralbinaryandinherently
Cartesian codes.
Cartesian codes.
Key words. Hydrodynamics
accretion, accretion– protoplanetary disksdisks– Methods:
numerical numerical
Key words. Hydrodynamics– accretion,– accretion
disksdisks– protoplanetary– Methods:

νnum
1. Introductionx^2

1. Introduction

Hydrodynamical simulations are a useful tool in studying a variety of astrophysical processes. An example of such a process isthe accretion flow of a protoplanetary disk around a star. AccreHydrodynamical simulations
are a useful tool in studying a varition is thought to be achieved via turbulence operating in parts
ety of astrophysical processes.
exampleofthatsuch
a processisof the disk (LyraAn& Umurhan2019),
transports
angular momentum
radially
outwards
similartoaviscosity(Lynden-Bellthe accretion flow of&a Pringle
protoplanetary
disk around a star. Accre1974; Balbus & Papaloizou 1999), or via a magnetic
tion is thought to betorqueachieved
via turbulencein parts
that is appliedon the material operatingclose to the surfaceof thedisk
asastellarwind
expels
gas
from
the
disk
surfacelayers(Baiof the disk (Lyra & Umurhan 2019), that transports angular
mo-

α

x
10^(-4)winds are becoming a favored means of interpreting accretion in
protoplanetary disks, and numerical models of such disks tend
windsthearebecoming
favoredmeans of interpreting accretion in
towards
inviscid
limit (e.g.,aLega
et al. 2022).
protoplanetarydisks,
numerical
Here, however, we run
into a anddi↵erent
problem: themodels
numeri- of such disks tendcalschemes
used
to
model
protoplanetarydisks
introducecer-2022).
towards the inviscid limit (e.g., Lega et aal.
tain amount of numerical di↵usion. The magnitude of this di↵usion, isHere,often unknown
or poorly
calculatingproblem: the numerihowever,
we constrained.
run into Thus,a di↵erent
an upper limit to this non-physical di↵usion, or equivalently acal schemes
usedto modelprotoplanetary"numericalviscosity",
is essential
in ensuringthat results are notdisks introduce a cera↵ectedby the e↵ects
can induce. Thisis particularly
tain amountof itnumericaldi↵usion.
Theimpormagnitude of this di↵utant for models using Cartesian grids, which are primarily used2∗

College of Science, Ibaraki University, 2-1-1 Bunkyo,
Mito, Ibaraki
310-8512, Japan
Kinematic
signaturesof a low-mass planet with
School of Natural Sciences, Institute for Advanced Study, Princeton, NJ 08544, USA
E-mail: kanagawa.k@gmail.com

Received 2022 November 10; Accepted 2023 August 16

Abstracta moderately inclined orbit in a protoplanetary
disk
Kazuhiro D. KANAGAWA,1 ,∗ Tomohiro ONO,2 and Munetake MOMOSE1
1
2

College of Science, Ibaraki University, 2-1-1 Bunkyo, Mito, Ibaraki 310-8512, Japan
School of Natural Sciences, Institute for Advanced Study, Princeton, NJ 08544, USA

A planet embedded in a protoplanetary∗ disk produces a gap by disk–planet interaction.
E-mail: kanagawa.k@gmail.com
It also generates velocity perturbation ofgas, which can also be observed as deviations
Received 2022 November 10; Accepted 2023 August 16
from the Keplerian rotation in the channel map of molecular line emission, called kineAbstract
matic planetary features. These observed
signatures provide clues to determine the massA
planetembeddedby
in aaprotoplanetarya gap
by disk–planet interaction.
of the planet. We investigated the features induced
planet withdisk
anproduces
inclinedorbit
It also generates velocity perturbation of gas, which can also be observed as deviations
through three-dimensional hydrodynamic
simulations. We found that a smaller planet,from the Keplerian rotation in the channel map of molecular line emission, called kine◦◦with an inclination of ∼10 –20 , can produce
kinematic
featuresas prominent
as those
matic planetaryfeatures.
These observed
signatures provide
clues to determine the massoftheplanet.
We
investigatedthefeaturesinducedbyaplanetinduced by a massive coplanar planet. Despite the kinematic features being similar, the with an inclined orbit
through three-dimensional hydrodynamic simulations. We found that a smaller planet,
gap is shallower and narrower comparedwith the case in which the kinematic featureswith an inclination of ∼10◦ –20◦ , can produce kinematic features as prominent as thoseare formed by a coplanar planet. We also
found
the kinematicfeaturesby featuresan
inducedbythata massive
coplanar planet.
Despiteinduced
the kinematic
being⾒similar, theis shallower and because
narrower comparedwith thecase in which the kinematic features
inclined planet were fainter for rarer COgap
isotopologues
the velocity
perturbation
are formed by a coplanar planet. We also found that the kinematic features induced by anis weaker at the position closer to the inclined
midplane,
which was different in the case with aplanet were fainter for rarer CO isotopologues because the velocity
CO perturbation
coplanar massive planet. This dependence
on attheis midplane,
distinguished
if different
theis weaker
the isotopologues
position closer to thewhich wasin the case with acoplanarplanet.
This dependence
on theplanet has an inclined orbit. We discussed
twomassive
observedkinematic
featuresin isotopologues
the disk of is distinguished if theplanet has an inclined orbit. We discussed two observed kinematic features in the disk of
HD 163296. We concluded that the kink observed at 220 au can be induced by an inclined
HD 163296. We concluded that the kink observed at 220 au can be induced by an inclined
planet, while the kink at 67 au is consistent
towhilethattheinducedcoplanarplanet,kink at 67by
au isaconsistent
to planet.
that induced by a coplanar planet.
Key words: protoplanetary
disk — planet— disk interaction
Key words: protoplanetary disk — planet — disk interaction— submillimeter:planetarysystems— submillimeter: planetary systems

1 Introduction

revealed a number of protoplanetary disks with gap structures (e.g.,disks
ALMA
Partnership
2015; Tsukagoshi et al.
a numberofgravitaprotoplanetarywithgap strucA planet is formed in revealed
a protoplanetary
disk and1 Introduction
2016
;
Long
et
al.
2018
;
Andrews
et al. 2018; van dertionally interacts with thesurroundinggas
via spiral
waves
tures
(e.g.,
ALMA
Partnership
2015;
Tsukagoshi
etal.
Marel et al. 2019, Cieza et al. 2019). The planetary gap is
A planet is formed in a protoplanetary disk (Linand&gravitaPapaloizou 1979; Goldreich & Tremaine 1980).
2016
;
Long
et
al.
2018
;
Andrews
al. 2018; vander
one of theetconvincing
originswhich
explain these observed
When a planet is sufficiently large, a gas-depleted region

tionally interacts with the surrounding gas via spiral waves

Downloaded from https://academic.oup.com/pasj/advance-article/doi/10.1093/pasj/psad059/7273622 by Dankook University / Library user on 10 Octo1

academic.oup.com/pasj/advance-article/doi/10.1093/pasj/psad059/7273622 by Dankook University / Library user on

Kazuhiro D. KANAGAWA,1 ,∗ Tomohiro ONO,2 and Munetake MOMOSE1

10°

20°

CO

Jeremy Karam

ABSTRACT‹

and Alison Sills

Department of Physics & Astronomy, McMaster University, 1280 Main Street West, Hamilton, ON L8S 4M1, Canada

Key words: stars: formation – stars: kinematics and dynamics – galaxies: star clusters: general.

Key words: stars: formation – stars: kinematics and dynamics – galaxies: star clusters: general.

GMC

Downloaded from https://academic.oup.com/mnras/article/521/4/5557/7083154 by Dankook University / Library user on 10 O

The formation of star clusters involves
the growth of smaller, gas-rich subclusters through accretion of gas from the giant
Accepted 2023 March 17. Received 2023 February 9; in original form 2022 December 1
molecular cloud (GMC) within which the subclusters are embedded. The two main accretion mechanisms responsible for thisare accretion of gas from dense filaments
A B S T Rand
A C Tfrom the ambient background of the cloud. We perform simulations of both of these
The formation
star clusters
involvessmoothed
the growth ofparticle
smaller, gas-rich
subclusters throughaccretion theof gasgas
fromandthe giantaccretion processes onto gas-rich starclustersofusing
coupledhydrodynamics
to modelN-bodymolecular cloud (GMC) within which the subclusters are embedded. The two main accretion mechanisms responsible for thisdynamics to model the stars. We findforofbothaccretion
processes,
accretingstar
clusterloses
someof itsof original
arethat,accretion
gas from
dense filaments
and from thethe
ambient
backgroundof thecloud. We
perform
simulations
both of thesemass
processes
onto gas-richor
startheclusters
using filament.
coupled smoothedhydrodynamics
modelfromthe gasbothand N-bodywhile gaining mass from either theaccretion
ambient
background
dense
Theparticle
amountof massto lostof these
dynamics to model the stars. We find that, for both accretion processes, the accreting star cluster loses some of its original massprocesses is small compared to the while
totalgaining
massmassof thecluster. However, in the case of accretion from a background medium, thefrom either the ambient background or the dense filament. The amount of mass lost from both of these−1
net effect can be a decrease in the total
massis ofthecomparedclusterto if
is travelling
fast Hoenoughambient
medium medium,(>4 kmprocesses
small
theittotal
mass of the cluster.
wever, inthroughthe case of theaccretionfrom a backgroundthes ).
−1
net effect can be a decrease in the total mass of the cluster if it is travelling fast enough through the ambient medium (>4 km s ).
We find that the amount of mass lostfrom the cluster through filamentary accretion is independent of the density, width, orWe find that the amount of mass lost from the cluster through filamentary accretion is independent of the density, width, or
number of filaments funnelling gas number
into theclusterfunnellingand isgasalwaysthatis althe
massof thethemassclusteris constantlyof filaments
intothe suchcluster andways
such thatof the cluster
is constantly increasing
increasing with withtime. We compare
our results toidealizedprescriptions
to modelformation
star cluster formationin larger
scaleGMC
GMC simulations
time. We compare our results to idealizedprescriptions
usedto modelstarusedclusterin larger
scale
simulations
and find that such prescriptions act as an upper limit when describing the mass of the star cluster they represent.
and find that such prescriptions act as an upper limit when describing the mass of the star cluster they represent.

The sizes of these star clusters can vary, but numerical simulations
1 I N T RO D U C T I O N
of GMC evolution have shown that smaller star clusters (subclusters)
Star-forming regions throughout galaxies are found to comprise of
can evolve into more massive ones. Two important mechanisms for
dense molecular gas that can collapse to form stars. These stars form
this growth are subcluster mergers and accretion of the backgroundThe sizesof thesevary, butnumerical
clusteredgroupsthatcancontainmany
stars
gravitationallyboundgas instar
whichclusters
the clustercan
is embedded(Howard,Pudritz simulations& Harris
1 I N T RO D U C T I O N
to one another (called star clusters). Such clusters
areembedded2018;
hereafter
H18
).
We
focus
on
the
latterinthiswork(subclusters)
and pointof GMC evolution have shown that smaller star clusters
inside these giant molecular clouds (GMCs) in their early stages of
the reader to Karam & Sills (2022) (hereafter Paper I) for a study ofStar-forming regions throughout galaxies
are found to comprise of
can evolve into more
massive ones. Two important mechanisms forevolution (Lada & Lada 2003). This embedded phase of the star
subcluster mergers.
life tends
to last stars
only a few
Myr, after thiswhichgrowth
the clusterare subcluster
Simulating starformationandinsideaccretion
GMCs is computationallychaldense molecular gas that can collapse tocluster’s
form stars.
These
form
mergersof the backgroundis able to remove surrounding gas through stellar feedback effects
lenging because of the small timesteps required to properly modelclustered groups that can contain many(Pelupessy
stars gravitationallyboundgas in which thethecluster
is embedded (Howard, Pudritz & Harris& Portegies Zwart 2012; Li et al. 2019). Interactions
high-density gas. Therefore, in order to model star cluster
to one another (called star clusters). Such
arestarembedded
2018;molecularhereafter H18). We
onalongside
the latterin ofthisbetweenclusters
the embeddedcluster and the surrounding
formationand focus
evolution
the rest
the work
GMC, aandsinkpoint
gas
at
this
time
can
have
an
effect
on
the
o
v
erall
evolution
of
the
star
particleprescription
is
often
used.
A
commonlyusedprescription
inside these giant molecular clouds (GMCs) in their early stages of
the reader to Karam & Sills (2022) (hereafter Paper I) for a study ofcluster.
is that described in Federrath et al. (2010), which was used in H18.
evolution (Lada & Lada 2003). This embedded
phase ofthe starof the gassubcluster
mergers.
Because of turbulence,the distribution
surrounding starThis approach models a subcluster of gas and stars as a sink withclusters
in
GMCs
can
be
very
complicated.
For
example,observations
parametersthat describe
overall behaviour(i.e. the total mass, chalcluster’s life tends to last only a few Myr, after which the cluster
Simulating starformation
insideits GMCs
is computationally

4 km/s⾒2

Downloaded from https://academic.oup.com/mnras/article/521/4/5557/7083154 by Dankook Universi

Modelling
starformclusterformation:
Gas accretion
Accepted 2023 March 17. Received 2023 February
9; in original
2022 December
1Shinichi. W. Kinoshita1,2and Fumitaka Nakamura1,2,3

1
Department of Astronomy, the University of Tokyo 7-3-1 Hongo, Bunkyo, Tokyo 113-0033, Japan
2 Unveiling the Dynamics of Dense Cores in Cluster-forming Clumps: A 3D

3

National Astronomical Observatory of Japan NINS, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan

Magnetohydrodynamics Simulation Study of Angular Momentum and Magnetic Field
The Graduate University for Advanced Studies (SOKENDAI) 2-21-1 Osawa, Mitaka, Tokyo 181-0015, Japan; kinoshita.shinichi@grad.nao.ac.jp
Properties
Received 2023 May 25; revised 2023 August 17;
accepted 2023 August 17; published 2023 September 26
Shinichi. W. Kinoshita1,2and Fumitaka Nakamura1,2,3

1
Department of Astronomy, the University of Tokyo 7-3-1 Hongo, Bunkyo, Tokyo 113-0033, Japan
2

3

Abstract

National Astronomical Observatory of Japan NINS, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, JapanThe Graduate University for Advanced Studies (SOKENDAI) 2-21-1 Osawa, Mitaka, Tokyo 181-0015, Japan; kinoshita.shinichi@grad.nao.ac.jp
Received 2023 May 25; revised 2023 August 17; accepted 2023 August 17; published 2023 September 26

We conducted isothermal magnetohydrodynamics simulations with self-gravity to investigate the properties ofdense cores in cluster-forming clumps. Two different Abstract
setups were explored: a single rotating clump and colliding
clumps. We focusedon determining
the extent to which
the rotationand magneticof theofparental clump are
We conducted
isothermal magnetohydrodynamics
simulations
with self-gravity
to investigate ﬁeld
the properties
dense
cores
in
cluster-forming
clumps.
Two
different
setupswere
explored:asinglerotatingclump
andcollidinginherited by the formed dense cores. Our statistical analysis revealed that the alignment
between the angularclumps. We focused on determining the extent to which the rotation and magnetic ﬁeld of the parental clump are
Lcore
, and
rotationalaxis
of theclump
is inﬂuenced
by thetheangularstrength of turbulence
momentum of dense
cores,
inheritedby theformeddensethecores.
Our statisticalanalysis
revealed
that thealignment between
Lcore , and clumps,
the rotational
axisfoundof the clumpby the strength
turbulenceof dense
cores,align ofwith
the clump’s rotational
and the simulationmomentum
setup. In
singlerotating
wethat isLinﬂuenced
core tendstoand the simulation setup. In single rotating clumps, we found that Lcore tendsto align with the clump’s rotational
axis if the initial turbulenceis weak. In colliding clumps, however, this alignmentdoes not occur, regardless of theaxis if the initial turbulence is weak. In colliding clumps, however, this alignmentdoes not occur, regardless of the
turbulence
strength.
This misalignment
colliding clumps
is dueisto due
the induced
turbulencefrom turbulencetheinitial turbulence initialstrength.
This
misalignment
in incolliding
clumpsto theinducedfrom thecollision and the isotropic gas inﬂow into dense cores. Our analysis of colliding clumps also revealed that thecollision and the magnetic
isotropic
gas inﬂow into dense cores. Our analysis of colliding clumps also revealed that theﬁeld globally bends along the shock-compressed layer, and the mean magnetic ﬁeld of dense cores, Bcore ,magnetic ﬁeld globally
bends
along
the
shock-compressed
layer,weandtheand Lofaligns with
it. Bothin single
rotatingclumps and colliding clumps,
foundthatmeanthe anglemagnetic
between Bcoreﬁeldcore dense cores, Bcore ,isgenerally
random,regardless
of
the
clump
properties.
We
also
analyzedthedynamical
statesoftheformedcoresaligns with it. Bothandinfoundsinglerotating clumps and colliding clumps, we found that the angle between Bcore and Lcorea higher proportion of unbound cores in colliding clumps. In addition, the contribution of rotationalis generally random,
regardlessof the clump
properties.
also
analyzed
the dynamical
statesof the formed coresenergywas only approximately
5% of thegravitational Weenergy,regardlessof the model
parameters for both
singleandcollidingcases.
and found a higher proportion of unbound cores in colliding clumps. In addition, the contribution of rotational
Uniﬁed Astronomy Thesaurus
concepts:
Molecular gasenergy,(1073); Starregardlessformation (1569);
Magnetohydrodynamicalenergy was only approximately
5% of thegravitationalof themodel parameters for both single
simulations (1966)
and colliding cases.
Uniﬁed Astronomy Thesaurus
concepts: Molecular gas (1073);
Star formation(1569);
Vázquez-Semadeni
et al. 2011; Padoan
et al. Magnetohydrodynamical2014), allowing us to
1. Introduction
explorethecomplexinterplay
betweenenvironmental conditions,
simulations (1966)
Gravitationally bound dense cores undergo collapse to form

2

protostellar systems, which eventually evolve into stars. The
angular momentum of these cores is a crucial factor in thecreation of protostellar systems, as it plays a key role in theformationof protoplanetary disks (e.g., Belloche et al. 2002;
1. Introduction
Sai et al. 2023). The interplay between a rotating accretion disk
and a dense
magneticcores
ﬁeld isresponsiblecollapse
for the launching
Gravitationally boundundergoto formof
protostellar outﬂows (e.g., Tomisaka 2002; Matsumoto &protostellar systems,Tomisaka
which 2004;
eventually
The
Banerjee evolve& Pudritzinto2006).stars.
The angularthe protoplanetary
disk is factor
inherited in
fromthetheangular momentum momentumof theseofcoresis a crucial
dense core, and previous studies have shown that outﬂows tend
creation of protostellar
systems, as it plays a key role in theto be launched parallel to the disk’s angular momentum (e.g.,
formation of protoplanetarydisksMatsumoto(e.g., Belloche
al. Launhardt
2002;
Tomisaka 2002;& Tomisaka et
2004;
et al. 2009). In addition, the magnetic ﬁeld within collapsingangular momentum, and magnetic ﬁelds in cluster-forming
regions. By simulating the collapse of cluster-forming clumps
and the subsequent formation of dense cores, these simulations
can shed light on the origin and properties of the angularVázquez-Semadeni et al. 2011; Padoan et al. 2014), allowing us tomomentum and magnetic ﬁelds of cores, as well as their
explore fortheprotostellar
complexinterplay
betweenconditions,
implications
outﬂows
and disk
formation environmental(e.g.,
Chen&
Ostriker
2018;
Kuznetsova
et
al.
2019,
2020;
Arroyoangular momentum, and magnetic ﬁelds in cluster-forming
Chávez & Vázquez-Semadeni 2022; Misugi et al. 2023). For
regions.
simulating
the collapseof cluster-forming
clumpsexample,ChenBy& Ostriker(2018) investigatedthe properties
of
dense
cores
in
magnetohydrodynamics(MHD)
simulationsofand the subsequent formation of dense cores, these simulations
large-scale converging ﬂows. They suggest that the internal andcan magnetic
shed lightthe andoriginandmomentum
properties of the angularexternalﬁelds areon
correlated
the angularofmomentum
cores is acquired fromambient
turbulence.
Misugietal.
and magnetic ﬁelds of (2023)
cores, as well as their

5

NRC Herzberg Astronomy and Astrophysics, 5071 West Saanich Road, Victoria, BC, V9E 2E7, Canada
9Department of Physics and Astronomy, University of Victoria, Victoria, BC, V8P 5C2, Canada
10Complex
Organic
Molecules
in a520
Very
Young
HotCharlottesville,
Corino, HOPS
373SW
National
Radio
Astronomy
Observatory,
Edgemont
Road,
VA 22903,
USA
Received 2022 September
18; revised
2023 June
June 29; published
2023 October
5
1,2
1,3
4 27; accepted 2023
1
1
5
Jeong-Eun Lee

1, Giseon Baek
, Seokho Lee , Jae-Hong Jeong , Chul-Hwan Kim
Gregory J. Herczeg6,7 , Doug Johnstone8,9 , and John J. Tobin10, Yuri Aikawa,

Abstract
We present the spectra of complex organic molecules (COMs) detected in HOPS 373SW with the Atacama Large
Millimeter/submillimeter Array (ALMA). HOPS 373SW, which is a component of a protostellar binary with aseparation of 1500au, has been discovered as a variable protostar by the JCMT transient monitoring survey with amodest (∼30%) brightness increase at submillimeter wavelengths. Our ALMA target-of-opportunity observation at∼345 GHz for HOPS 373SW revealed extremely young chemical characteristics with strong deuteration of
methanol. The dust continuum opacity is very high toward
obscuring line emission from within
Abstractthe source center,17
component,
HOPSorganic
373NE,
was(COMs)
detected
only
in C373SW
O inwithour
0 03. The other binary
We present
the spectra of complex
molecules
detectedin HOPS
the observation,
Atacama Large implying a cold
Millimeter/submillimeter
Array(ALMA).
HOPS
373SW,
whichisacomponentofaprotostellar
binary witha373SW
with those ofand quiescent environment. We compare the COM abundances relative to CH3OH in HOPS
separation of 1500au, has been discovered as a variable protostar by the JCMT transient monitoring survey with aV883 Ori, which modest
is an (∼30%)
eruptive
diskincrease
object,as well asother hotto demonstrate
the chemical
evolutionbrightness
at submillimeter
wavelengths.
Our corinos,
ALMA target-of-opportunity
observation
at∼345 GHz
for abundances
HOPS 373SW revealed
extremely
young chemical
characteristics
with strong
deuterationfrom envelope to disk.
Highof singly,
doubly,and triply
deuterated
methanol(CHof2DOH, CHD2OH,
methanol. The dust continuum opacity is very high toward the source center, obscuring line emission from within
low
CH3CN
abundance in HOPS 373SW compared to other hot corinos suggest a very earlyand CD3OH) and 0a 03.
The otherbinary component, HOPS 373NE, was detected only in C17O in our observation, implying a cold
evolutionary stageandofquiescent
HOPSenvironment.
373SW Wein compare
the hotthecorino
phase.relative
Sinceto CH
the3OH
COMs
HOPS
373SW wouldin HOPSdetected
373SW withinthoseof
COM abundances
V883
Ori,
which
is
an
eruptivedisk
object,
as
well
as
other
hot
corinos,todemonstratethechemicalevolution
have been sublimatedvery recently from grain surfaces, HOPS 373SW is a promising place to study the surface
from envelope to disk. High abundances of singly, doubly, and triply deuterated methanol (CH2DOH, CHD2OH,
chemistry of COMs
prestellar
phase inbefore
sublimation.
OH) coldand a low
CH3CN abundance
HOPS 373SW
compared to other hot corinos suggest a very early
and inCD3theDepartment of Physics and Astronomy, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea; lee.jeongeun@snu.ac.kr
2
SNU Astronomy Research Center, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea
3
Research Institute of Basic Sciences, Seoul National University, Seoul 08826, Republic of Korea
4
Korea Astronomy and Space Science Institute, 776 Daedeok-daero, Yuseong, Daejeon 34055, Republic of Korea
5
Department of Astronomy, University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan
6
Kavli Institute for Astronomy and Astrophysics, Peking University, Yiheyuan 5, Haidian Qu, 100871 Beijing, Peopleʼs Republic of China
7
Department of Astronomy, Peking University, Yiheyuan 5, Haidian Qu, 100871 Beijing, Peopleʼs Republic of China
8
NRC Herzberg Astronomy and Astrophysics, 5071 West Saanich Road, Victoria, BC, V9E 2E7, Canada
9Department of Physics and Astronomy, University of Victoria, Victoria, BC, V8P 5C2, Canada
10
National Radio Astronomy Observatory, 520 Edgemont Road, Charlottesville, VA 22903, USA
Received 2022 September 18; revised 2023 June 27; accepted 2023 June 29; published 2023 October 5

evolutionary stage of HOPS 373SW in the hot corino phase. Since the COMs detected in HOPS 373SW would

Uniﬁed Astronomy
Thesaurus
concepts:
Star
Astrochemistry
(75)
have
been sublimatedvery recently
fromformation
grain surfaces,(1569);
HOPS 373SW
is a promising place
to study the surface
HOPS 373SW
1500auCOMschemistry of COMs in the cold prestellar phase before sublimation.

ALMA

Uniﬁed Astronomy Thesaurus concepts: Star formation (1569); Astrochemistry (75)

1. Introduction

HOPS
373SW
Cometsand asteroidsComets
revealthereveal
solar
was wasrichrichinin
and that
asteroids
thatnebula
the solar nebulawater andCOMs
organic
recent Rosetta
missionalso
water and organic molecules.
Themolecules.
recent The
Rosettamissionalso
showed that many complex organic molecules (COMs; carbon1. Introduction

showed that many complex
organic molecules (COMs; carbonbearing molecules made of six or more atoms), as well as sobearing molecules madesix molecules,
or morewhichatoms),
asprecursors
well asof amino
socalled ofprebiotic
could beacids
andsugars,
exist
in
the
comet
67P/C-G(Altwegg
et al.
called prebiotic molecules,
which could be precursors of amino2017). These organics, together with water, could have been
acids and sugars, exist
in the67P/C-G(Altwegg
et al.
brought
to thecomet
young Earth’s
surface by
comets and asteroids,
whichtogether
are believed
to preservesomecompositions2017). These organics,
withwater,
could
have frombeentheearly stage of star formation. Therefore, the investigation ofbrought to the youngCOMs
Earth’s
by environments
comets andasteroids,
in thesurface
circumstellarof Sun-likeyoungstellarobjectsisimportanttobetterinfer
thechemicalwhich are believed to preserve some compositions from the345 GHz

temperature (100 K) of the water-ice mantles of dust grains, areCH CNtemperature(100inK) COMsof the water-ice
mantles of dust
arevery rich(Ceccarelli
et grains,al. 2007),
although not allvery rich in COMs (Ceccarelli et al. 2007), although not all
3
protostars
are observed
tocorinos.
develop
hot corinos.
Hot corinos areprotostars
are observed
to develop hot
Hot corinos
arediscovered
predominantly
toward the youngest
mostyoungest and most
discovered
predominantly
towardandthe
embedded Class 0 protostars with bolometric luminosities
embeddedClass 0et al.protostars
luminosities
greaterthan 5 Le (Belloche
2020; Yang et withal. 2020;bolometric
Hsu
etgreateral. 2022). thanOn the5otherCOMs as etrichal.
as thosein hotLehand,(Belloche
2020;
Yang et al. 2020; Hsu
corinos have been rarely detected in disks, except for special
et al.
On
theofotherCOMs
as rich as those in hot
cases,
such2022).
as the disk
surface
HH 212hand,(Lee et al.
2019a) andtheeruptivedisk
midplane
of
V883
Ori
(van
’t
Hoff
et
al.
2018;
corinos have been rarely detected in disks, except for special
Lee et al. 2019b).
cases,
such
asforthediskrequire
surfaceof andHHice-phase
212 (Lee et al. 2019a) andChemicalmodelsCOMs
both gasprocesses
in protostellar
ice chemistry
the eruptive
diskenvironments,
midplane althoughof V883
Ori (van ’t Hoff et al. 2018;is essential to the formation of COMs. During the warm-up

Hanpu Liu

ABSTRACT

1,2 ‹and Xue-Ning Bai

3,4 ‹1 Kavli Institute for Astronomy and Astrophysics, Peking University, 100871 Beijing, People’s Republic of China
3

Downloaded from https://academic.oup.com/mnras/article/526/1/80/7259167 by Dankook University / Library user on 10 Octo

2 Department of Astronomy, Peking University, 100871 Beijing, People’s Republic of China

Institute for Advanced Study, Tsinghua University,100084 Beijing , China
Recent numerical simulations
revealed that dust clumping and planetesimal formation likely proceed in ring-like disc
4 Departmenthaveof Astronomy, Tsinghua University, 100084 Beijing, China
substructures, where dust gets trapped in weakly turbulent pressure maxima. The streaming instability has difficulty operating insuch rings with external turbulenceandno pressure
To2023
exploreAccepted 2023 August
24. Received
2023 Augustgradient.
24; in original form
July 14 potential paths to planetesimal formation in this context,
we analyse the stability of turbulent dust-trapping ring under the shearing sheet framework. We self-consistently establish theABST
R A C ringTpressure maximum and thedust
in equilibrium, the former via a balance of external forcing versus viscosity and the latter
Recent numerical simulations have revealed that dust clumping and planetesimal formation likely proceed in ring-like discvia dust drift versus turbulent
diffusion. We find two types of ! H-scale instabilities (H being the pressure scale height), which
substructures, where dust gets trapped in weakly turbulent pressure maxima. The streaming instability has difficulty operating in
we term the dusty Rossbysuch
wave
instability(DRWI).
Type
I isgradient.
generalized
thestandard
Rossbywave
instability(RWI, whichrings withexternal turbulence
and nopressure
To explorefrom
potentialpaths
to planetesimalformationin thiscontext,
we analyse the stability of turbulent dust-trapping ring under the shearing sheet framework. We self-consistently establish theis stationary at the pressure
maximum and dominates in relatively sharp pressure bumps. Type II is a newly identified travellingpressure maximum and the dust ring in equilibrium, the former via a balance of external forcing versus viscosity and the lattermode that requires the presencedust.
It canoperatein tworelatively
many scale
thatheight),
are stable
via dust driftof
versusturbulent
diffusion.
We findtypes of ! mild
H-scalebumps,
instabilitiesincluding(H being the pressure
which to the standard
we term the dusty Rossby wave instability (DRWI). Type I is generalized from the standard Rossby wave instability (RWI, which
RWI, and its growth rate is largely determined by the equilibrium gas and dust density gradients. We further conduct two-fluid
is stationary at the pressure maximum and dominates in relatively sharp pressure bumps. Type II is a newly identified travelling
simulations that verify thetypestheofpresence
the DRWI.
While
I leads
strongincluding
to dustconcentrationa large gas vortex
modetwothat requiresof dust. It canoperateTypein relatively
mild bumps,
manythat are stable to theinto
standard
RWI, andgrowthrate isdeterminedthe equilibrium
gas and dust density
gradients. We
further conduct
two-fluid within the ring.
similar to the standard RWI,
theitsdust
ringislargelypreserved
in by
Type
II, and meanwhile
exhibiting
additional
clumping
simulations that verify the two types of the DRWI. While Type I leads strong to dust concentration into a large gas vortex
The DRWI suggests a promising
path towards
formationof planetesimals/planetary
embryos
and azimuthally
similar to the standard
RWI, the dust
ring is preservedin Type II, and meanwhile exhibiting
additional clumping
within the ring.asymmetric dust
DRWI suggests arings.
promising path towards formation of planetesimals/planetary embryos and azimuthally asymmetric dust
structure from turbulent The
dust-trapping

structure from turbulent dust-trapping rings.

Key words:
hydrodynamics
instabilities – methods:
analytical– methods: numerical– planets
and satellites:–
Key words: hydrodynamics– instabilities––methods:
analytical– methods:
numerical– planets
and formation
satellites:
formation –protoplanetary discs.
protoplanetary discs.

Iabundance (vertically integrated dust-to-gas mass ratio) exceeds acertain threshold (depending on dust size, typically ! 0.02; Bai &Stone 2010; but see Li & Youdin 2021), the SI is found in simulations
It has recently been established that ring-like substructures areabundance(vertically integrated dust-to-gas mass ratio) exceeds ato lead to efficient dust clumping, with clumps dense enough to
1 I N T RO D U C T I O N ubiquitous among extended protoplanetary discs (PPDs), as revealed
threshold directly(dependingon dustsize,
typically ! 0.02; Bai &by ALMA [ALMA Partnership et al. 2015; Andrews et al. 2018; for certain
form planetesimals
by gravitational
collapse(e.g. Johansen,
a review, see Andrews (2020)]. While the formation mechanisms of
Youdin & Mac Low 2009; Carrera, Johansen & Davies 2015; Yang,
Stone
2010; but see Li & Youdin 2021), the SI is found in simulations
It has recently been establishedthat ring-like substructures are
such ring-like substructures are debated [see e.g. Bae et al. (2022)
Johansen & Carrera 2017). However, if turbulent dust-trapping ringubiquitous among extendedforprotoplanetary
discs
astrapping
revealed
to efficientdust clumping,with clumps
a review], they are believed(PPDs),
to reflect dustin turbulent to lead
like substructures
are common
as found in observations,
the SI dense enough topressure bumps (e.g. Dullemond et al. 2018; Rosotti et al. 2020).
paradigm for planetesimal formation in such dust rings faces two
by ALMA [ALMA Partnershipet al. 2015; Andrews et al. 2018; forform
planetesimals directly by gravitational collapse (e.g. Johansen,
Such dust-trapping sites not only retain the dust by preventing or
challenges. First, most existing simulations did not include external
1 I N T RO D U C T I O NIIloaded from https://academic.oup.com/mnras/article/526/1/80/7259167 by Dankook University / Library u

The dusty Rossby wave instability (DRWI): linear analysis andAccepted 2023 August 24. Received
2023 August
in original form
2023 July 14 rings in protoplanetary discs
simulationsof24;
turbulentdust-trapping