The Limited Role of the Streaming Instability during Moon and Exomoon Formation

The Limited Role of the Streaming Instability during Moon and Exomoon Formation
Miki Nakajima1,2
, Jeremy Atkins2
, Jacob B. Simon3
, and Alice C. Quillen2
1
Department of Earth and Environmental Sciences, University of Rochester, P.O. Box 270221, Rochester, NY 14627, USA
2
Department of Physics and Astronomy, University of Rochester, P.O. Box 270171, Rochester, NY 14627, USA
3
Department of Physics and Astronomy, Iowa State University, Ames, IA 50010, USA
Received 2023 September 20; revised 2024 March 27; accepted 2024 April 17; published 2024 June 17
Abstract
It is generally accepted that the Moon accreted from the disk formed by an impact between the proto-Earth and
impactor, but its details are highly debated. Some models suggest that a Mars-sized impactor formed a silicate
melt-rich (vapor-poor) disk around Earth, whereas other models suggest that a highly energetic impact produced a
silicate vapor-rich disk. Such a vapor-rich disk, however, may not be suitable for the Moon formation, because
moonlets, building blocks of the Moon, of 100 m–100 km in radius may experience strong gas drag and fall onto
Earth on a short timescale, failing to grow further. This problem may be avoided if large moonlets (?100 km)
form very quickly by streaming instability, which is a process to concentrate particles enough to cause gravitational
collapse and rapid formation of planetesimals or moonlets. Here, we investigate the effect of the streaming
instability in the Moon-forming disk for the first time and find that this instability can quickly form ∼100 km-sized
moonlets. However, these moonlets are not large enough to avoid strong drag, and they still fall onto Earth quickly.
This suggests that the vapor-rich disks may not form the large Moon, and therefore the models that produce vapor-
poor disks are supported. This result is applicable to general impact-induced moon-forming disks, supporting the
previous suggestion that small planets (<1.6 R⊕) are good candidates to host large moons because their impact-
induced disks would likely be vapor-poor. We find a limited role of streaming instability in satellite formation in an
impact-induced disk, whereas it plays a key role during planet formation.
Unified Astronomy Thesaurus concepts: Earth-moon system (436)
Supporting material: tar.gz file
1. Introduction
1.1. Review of the Moon-formation Hypotheses
The origin of the Earth’s Moon is a long-standing open
problem in planetary science. While it is accepted that the
Moon formed from a partially vaporized disk generated by a
collision between the proto-Earth and an impactor approxi-
mately 4.5 billion yr ago (the “giant impact hypothesis”;
Hartmann & Davis 1975; Cameron & Ward 1976), the details,
such as the impactor radius and velocity, are actively debated.
According to the canonical model, the proto-Earth was hit by a
Mars-sized impactor (Canup & Asphaug 2001). This model can
successfully explain several observations, such as the mass of
the Moon, the angular momentum of the Earth–Moon system,
and the small lunar core. Moreover, this model could
potentially explain the observation that the Moon is depleted
in volatiles, if they escaped during the impact or subsequent
processes (see reviews by Halliday 2004; Canup et al. 2023).
However, this model fails to explain the Earth and Moon’s
nearly identical isotopic ratios (e.g., Si, Ti, W, and O;
Armytage et al. 2012; Zhang & Dauphas 2012; Touboul
et al. 2015; Thiemens et al. 2019, 2021; Kruijer et al. 2021).
This is because the disk generated by an oblique Mars-sized
impactor is primarily made of the impactor materials
(Canup 2004; Kegerreis et al. 2022; Hull et al. 2024), which
likely had different isotopic ratios from Earth, unless the inner
solar system was well mixed and homogenized (Dauphas 2017)
or equilibrated in the disk phase (Pahlevan & Stevenson 2007;
Lock et al. 2018). In contrast, more energetic impact models
may solve this problem by mixing the proto-Earth and
impactor. These energetic models include an impact between
two half-Earth-sized objects (Canup 2012), an impact between
a rapidly rotating proto-Earth and a small impactor (Cúk &
Stewart 2012), and general impact events that involve high
energy and high angular momentum, which create a doughnut-
shaped vapor disk connected to the Earth (a so-called synestia;
Lock et al. 2018). In these scenarios, the disk and the proto-
Earth compositions could have been mixed well, potentially
solving the isotopic problem. Other proposed scenarios include
(a) the Moon formation by multiple small impacts (Rufu et al.
2017), (b) hit-and-run collisions that mix the proto-Earth and
the Moon more than the canonical model does (Asphaug et al.
2021), and (c) an impact between the proto-Earth covered by a
surface magma ocean and an impactor, which leads to more
contribution from Earth to the disk than the canonical scenario
(Hosono et al. 2019).
All of these models have potential explanations for the
isotopic similarity, but each model faces challenges to explain
other constraints. For example, the half-Earth and synestia
models predict a much higher angular momentum of the Earth–
Moon system than that of today (Canup 2012; Cúk &
Stewart 2012). The excess of the angular momentum may or
may not be removed easily (Cúk & Stewart 2012; Cúk et al.
2016, 2021; Rufu et al. 2017; Rufu & Canup 2020; Ward et al.
2020). Another constraint may come from the Earth’s interior;
the deep portion of the Earth still retains solar neon (Yokochi &
Marty 2004; Williams & Mukhopadhyay 2019) and helium
isotopic ratios (Williams et al. 2019), which are distinct from
the rest of the mantle. These noble gases may have been
The Planetary Science Journal, 5:145 (12pp), 2024 June https://doi.org/10.3847/PSJ/ad4863
© 2024. The Author(s). Published by the American Astronomical Society.
Original content from this work may be used under the terms
of the Creative Commons Attribution 4.0 licence. Any further
distribution of this work must maintain attribution to the author(s) and the title
of the work, journal citation and DOI.
1

delivered to Earth’s mantle during the planet formation phase.
This suggests that Earth’s mantle developed chemical hetero-
geneity during Earth’s accretion, when the protoplanetary disk
was still present, and that the Earth’s mantle has never been
completely mixed, even by the Moon-forming impact. Previous
work suggests that the preservation of the mantle heterogeneity
can be explained by the canonical model but not by the
energetic models, because these energetic impacts tend to mix
the mantle (Nakajima & Stevenson 2015). It should be noted,
however, the possibility that the neon and helium were
delivered to the Earth's core and have been slowly leaking
into the mantle over time has been discussed (Bouhifd et al.
2020), which would eliminate the need to preserve the mantle
heterogeneity. However, it is not clear if delivery of these noble
gases to the core was efficient.
Moreover, due to recent analytical capability, small isotopic
differences between Earth and the Moon have been observed
(e.g., K, O, W, V, Cr, and others; Wiechert et al. 2001; Kruijer
et al. 2015; Touboul et al. 2015; Wang & Jacobsen 2016;
Young et al. 2016; Thiemens et al. 2019; Nielsen et al. 2021;
Sossi et al. 2018). Some isotopic difference, such as the
Moon’s enrichment in heavy K isotopes compared to those of
Earth, can be explained by isotopic fractionation due to liquid–
vapor phase separation during the Moon accretion process
(Wang & Jacobsen 2016; Nie & Dauphas 2019; Charnoz et al.
2021). Observed small oxygen isotopic differences between
Earth and the Moon may suggest that the Moon still records the
impactor component (Cano et al. 2020). This suggestion may
not be compatible with the energetic impact models because
these impacts are so energetic that Earth and the Moon would
be efficiently mixed. Whether the proposed models for the
lunar origin can explain the observation that the Moon is
depleted in volatiles is actively debated (Canup et al. 2015,
2023; Dauphas et al. 2015, 2022; Nakajima & Stevenson 2018;
Sossi et al. 2018; Nie & Dauphas 2019; Charnoz et al. 2021;
Halliday & Canup 2022; Canup et al. 2023).
1.2. Gas Drag Problem in a Vapor-rich Disk
Another key constraint that has not been discussed until
recently is the vapor mass fraction (VMF) of the disk, which
sensitively depends on the impact model. Less energetic
impacts, such as the canonical model and the multiple impact
model, produce relatively small VMFs (VMF ∼ 0.2 for the
canonical model, Nakajima & Stevenson 2014, and ∼0.1–0.5
for the multiple impact model, Rufu et al. 2017), while more
energetic models, such as the half-Earth model and synestia
model, produce nearly pure vapor disks (VMF ∼0.8–1,
Nakajima & Stevenson 2014). The VMF of the Moon-forming
disk can significantly impact the Moon accretion process; if the
initial Moon-forming disk is vapor-poor (liquid-rich), moonlets
can quickly form by gravitational instability (GI) from the
liquid layer in the disk midplane outside the Roche radius
(Thompson & Stevenson 1988; Salmon & Canup 2012). These
moonlets eventually accrete to the Moon in tens to hundreds of
years (Thompson & Stevenson 1988; Salmon & Canup 2012;
Lock et al. 2018). In contrast, an initially vapor-rich disk needs
to cool until liquid droplets emerge before moonlet accretion
begins. Growing moonlets in the vapor-rich disk experience
strong gas drag from the vapor (Nakajima et al. 2022). The gas
drag effect is strongest when the gas and moonlet are coupled
to an adequate degree and is sensitive to the moonlet radius. It
is strongest with a moonlet radius Rp of a few kilometers
(Nakajima et al. 2022). In contrast, much smaller moonlets are
completely coupled with the gas, and much larger moonlets are
completely decoupled from the gas, experiencing a weaker gas
drag effect. As a result, ∼kilometer-sized moonlets lose their
angular momentum and inspiral onto the Earth within 1 day
(Nakajima et al. 2022), a much shorter timescale than that
required for lunar formation.
This same problem was a major challenge to planet
formation in the protoplanetary disk (Adachi et al. 1976;
Weidenschilling 1977). Here, we quickly review the gas drag
problem in the protoplanetary disk. The radial velocity of a
particle in a disk is (see Equation (A6) for the derivation; see
also Armitage 2010; Takeuchi & Lin 2002)
( )
v
v v
2
, 1
r
r
f
1
,g K
f f
1
t h
t t
=
-
+
-
-
where τf = Ωtf is the dimensionless stopping time (Armitage 2010),
Ω is the Keplerian angular velocity, tf is the friction time
(see further descriptions below), vr,g is the gas velocity, and η is the
pressure gradient parameter, described as
⎜ ⎟
⎛
⎝
⎞
⎠
( )
c
v
p
r
1
2
ln
ln
, 2
s
K
2
g
h = -
¶
¶
where cs is the sound speed, vK is the Keplerian velocity, pg is
the gas pressure, and r is the radial distance. vr is largest when
τf = 1. The friction time is the time until the particle and gas
reach the same velocity and is defined as tf = mvrel/FD, where
m is the particle mass, vrel is the relative azimuthal velocity
between the gas and particle, and FD is the drag force. In the
protoplanetary disk, the gas drag for a particle is generally
written as F R v
C
D 2 g,0 p
2
rel
2
D
pr
= - , where CD is the gas drag
coefficient, ρg,0 is the initial gas density at the midplane, and Rp
is the particle radius. The gas drag coefficients are roughly
(Armitage 2010)
⎧
⎨
⎩
( )
( )
( )
( )
C
24 Re , Re 1 Stokesregime
24 Re , 1 Re 800 transitionregime
0.44, Re 800 Newtonregime
, 3
D
1
0.6
=
<
< <
>
-
-
where Re is the Reynolds number. Assuming that vrel ∼ ηvK
(see Equation (A2)), ν ∼ csλ, where ν is the kinematic
viscosity, λ is the mean free path (
k T
d p
2
B
2
g,0
l =
p
, where kB is
the Boltzmann constant, T is the temperature, d is the molecular
diameter, and pg,0 is the gas pressure at the midplane), Re =
4.26
v Rp
rel
~
n
, which is in the transition regime, in the
protoplanetary disk at 1 au, assuming T = 280 K, η = 0.002,
d = 289 pm, Rp = 0.52 m (see discussion below), p
RT
m
g,0
g,0
m
=
r
,
and mm = 0.002 kg mol–1
, where mm is the mean molecular
weight and R is the gas constant. Here we also assume the
following vertical distribution of the gas density ρg,
⎜ ⎟
⎛
⎝
⎞
⎠
( )
z
H
exp
2
, 4
g g,0
2
2
r r
= -
where z is the vertical coordinate and H is the gas scale height.
This leads to
H
g,0 2
g
r =
p
S
. We assume Σg = 17,000 kg m−2
and
use the relationships of c RT m
s m
= , H = cs/Ω,
p
r
ln
ln
g
=
¶
¶
2
The Planetary Science Journal, 5:145 (12pp), 2024 June Nakajima et al.

( )
3 2
b
- + , and
1
2
b = . The orbital angular frequency is
GM r3
W = , where G is the gravitational constant, M is the
stellar mass and r is the distance from the star (for the Moon-
forming disk, M is the Earth mass and r is the distance from
Earth). These parameters are taken from previous work
(Carrera et al. 2015; other values of β have been discussed,
such as
3
7
b = ; Chiang & Youdin 2010).
The dimensionless stopping time becomes
( )
t
C
R
v C
R
r
8
3
8
3
. 5
f f
D
p p
g,0 rel D
p p
g,0
t
r
r h
r
r
= W = W =
τf = 1 when Rp = 0.52 m. The particle density ρp = 3000 kg
m−3
is assumed. The residence time of the particle at 1 au is
1 au/vr = 75 yr, where vr is the radial fall velocity of the
particle (see Equation (A7); for simplicity, the radial gas
velocity vr,g = 0 is assumed). Thus, an approximately meter-
sized particle falls toward the Sun at 1 au within 80 yr, a
timescale much shorter than the planet formation timescale
(several tens of Myr). This is the so-called “meter barrier”
problem, which was a major issue to explain planetary growth
in a classical planet formation model (Adachi et al. 1976;
Weidenschilling 1977).
In contrast, for the vapor-rich Moon-forming disk, τf = 1
(Re 1010
> ) is achieved when Rp = 1.8 km, assuming
ρg,0 = 40 kg m–3
, ρp = 3000 kg m–3
, r = 3 R⊕, η = 0.04,
T = 5200 K (see Section 2.3 for justifications for these parameters),
and d = 300 pm. The residence time, r/vr, is approximately 1 day
(1.21 days), which is much shorter than the Moon-formation
timescale of several tens of yr to 100 yr. This formation timescale
is ultimately determined by the radiative cooling timescale, but it is
model-dependent. Here we provide a very simple estimate; the
timescale for radiative cooling is
( )
10
M L C T
r T
4
p
disk
2
ph
4 ~
p s
+ D
yr (this is
also consistent with numerical work; Lock et al. 2018), where
Mdisk is the disk mass, L is the latent heat, Cp is the specific heat,
ΔT is the temperature change over time, σ is the Stefan–
Boltzmann constant, and Tph is the photosphere temperature. Here,
Mdisk ∼ 0.015 M⊕, L = 1.2 × 107
J kg–1
(Melosh 2007), Cp =
103
J K–1
kg–1
, ΔT = 2000 K, and Tph = 2000 K (Thompson &
Stevenson 1988) are assumed. However, the actual Moon-
formation timescale can be longer than this for several reasons,
including additional heating due to viscous spreading (Thompson
& Stevenson 1988; Charnoz & Michaut 2015) and radial material
transport efficiency (e.g., Salmon & Canup 2012). A short Moon-
formation timescale has been proposed (e.g., Mullen & Gam-
mie 2020; Kegerreis et al. 2022), but the typical Moon-formation
timescale has been estimated to range in 10–102
yr.
The vertical settling velocity of condensing particles is
v C
R z
settle
8
3 D
p p
2
g
=
r
r
W
(Armitage 2010). Assuming z ∼ H, the
settling time for a particle with a radius off 2 is 0.22 days.
Determining the collision history of moonlets requires
conducting orbital dynamics simulations, but here we provide
a rough estimate. A rough estimate of the mass-doubling time
of the largest moonlet in the Moon-forming disk is typically
∼1 day (Salmon & Canup 2012). If a 2 km-sized moonlet mass
doubles in a day, the radius change is very small (2 km ×
21/3
= 2.5 km), and the gas drag effect remains strong on such
a moonlet; this change does not prevent the moonlet from
falling into Earth. However, the actual collision time can vary
depending on the local concentration of particles, which needs
a detailed future study.
This gas drag effect is a problem with forming the Moon
from an initially vapor-rich disk. The Moon can still form after
most of the vapor condenses, but by that time a significant
portion of the disk mass could be lost (Ida et al. 2020;
Nakajima et al. 2022), which would fail to form a large moon
from an initially vapor-rich disk (here we use a ”large” moon
when its mass is ∼1 wt% or larger of the host planet). If this is
the case, an initially vapor-rich disk may not be capable of
forming a large Moon (Nakajima et al. 2022). In contrast, the
gas drag effect is weak for vapor-poor disks, which are
generated by less energetic models, such as the canonical and
multiple impact models.
1.3. Streaming Instability in the General Impact-induced
Moon-forming Disk
A potential solution to this vapor drag problem is forming a
large moonlet very quickly (much larger than 2 km), so that the
moonlet would not experience strong gas drag. This is the
accepted solution for the gas drag problem in the protoplane-
tary disk. The proposed mechanism is the streaming instability
(Youdin & Goodman 2005; Johansen et al. 2007), where
particles spontaneously concentrate in the disk, gravitationally
collapsing and forming a large clump (∼100 km in size;
Johansen et al. 2015). If this mechanism works for the Moon-
forming disk, an initially vapor-rich disk may be able to form
the Moon despite the gas drag issue. If this mechanism turns
out not to work for the Moon-forming disk, it is an interesting
finding as well, given that a Moon-forming disk is often treated
as a miniature analog of the protoplanetary disk. Understanding
what makes these two disks differ would deepen our under-
standing of planet and satellite formation processes.
Moreover, knowing whether the streaming instability can affect
moon formation processes informs our understanding of moon
formation in the solar and extrasolar systems. Moon formation in
an impact-induced disk is common in the solar system (e.g., the
moons of Mars (Craddock 2011), the moons of Uranus (Slattery
et al. 1992), and Pluto-Charon (Canup 2005)). While there are no
confirmed exomoons (moons around exoplanets) to date (e.g.,
Cassese & Kipping 2022), impact-induced exomoons should be
common because impacts are a common process during planet
formation (Nakajima et al. 2022) and because impacts in
extrasolar systems may have already been observed (e.g., Meng
et al. 2014; Bonomo et al. 2019; Thompson et al. 2019;
Kenworthy et al. 2023). If streaming instability operates in these
disks, it can affect what types of planets can host exomoons,
which can be compared with future exomoon observations.
It should also be noted that other instability, such as secular
GI (e.g., Youdin 2011; Takahashi & Inutsuka 2014; Tominaga
et al. 2018) and two-component viscous GI (TVGI; Tominaga
et al. 2019), have been discussed as a mechanism for
planetesimal and dust ring formation. The secular GI occurs
when dust–gas interaction reduces the rotational support of the
rotating disk, which leads to dust concentration. The TVGI is
an instability caused by dust–gas friction and turbulent gas
viscosity. These instabilities can lead to clump formation even
if the disk is self-gravitationally stable. The implications of
these instabilities are beyond the scope of this paper.
3

1.4. Motivation of This Work
The goal of this work is to investigate whether the streaming
instability can form moonlets that are large enough to avoid
strong gas drag from the vapor-rich disk. The result could
constrain the Moon-formation model as well as general impact-
induced models in the solar and extrasolar systems. We
conduct hydrodynamic simulations using the code Athena
(Stone et al. 2008; Bai & Stone 2010a). We first conduct 2D
simulations to identify the section of parameter space that leads
to streaming instability. Subsequently, we conduct 3D simula-
tions with self-gravity to identify the size distribution of the
moonlets. Finally, we identify the lifetime of the moonlets
formed by streaming instability to investigate whether it is
possible to form a large moon from an initially vapor-rich disk.
For the general impact-induced disks, we consider ”rocky” and
”icy” disks, where these disks form by collisions between
rocky planets (with silicate mantles and iron cores) and
between icy planets (with water-ice mantles and silicate cores),
respectively.
2. Method
2.1. Athena
We use the Athena hydrodynamics code, which solves the
equations of hydrodynamics using a second-order accurate
Godunov flux-conservative approach (Stone et al. 2008). We
use the configuration of Athena that couples the dimensionally
unsplit corner transport upwind method (Colella 1990) to the
third-order in-space piecewise parabolic method by Colella &
Woodward (1984) and calculates the numerical fluxes using the
HLLC Riemann solver (Toro 1999). We also integrate the
equations of motion for particles following Bai & Stone
(2010a) and include particle self-gravity for 3D simulations
following the particle-mesh approach described in previous
work (Simon et al. 2016). Orbital advection is taken into
account following previous work (Bai & Stone 2010a, 2010b).
Our setup is the local shearing box approximation in which a
small patch of the disk is corotating with the disk at the
Keplerian velocity (Stone & Gardiner 2010). The local
Cartesian frame is defined as (x, z) for 2D and (x, y, z) for
3D simulations, with x as the radial coordinate from the planet,
z parallel to the planetary spin axis, and y in the direction of the
orbital rotation.
In the Athena code, we solve the following equations for our
simulations (Bai & Stone 2010a; Simon et al. 2016; Li et al.
2019):
· ( ) ( )
t
u 0, 6
g
g
r
r
¶
¶
+  =
· ( )
( )
t
p
t
u
uu I
x z u
v u
3 2 , 7
g
g
g
g
2
g
2
g p
f
r
r
r r r r
W
¶
¶
+  +
= W - W - ´ +
-
ˆ ˆ
ˆ ( )
d
dt
v x
z
t
v
x x
z v
v u
a
2 3
2 . 8
i
i
i i
i
K
2
2
f
g
h
W
= - W + W
-W - ´ -
-
+
The first two equations specify mass and momentum
conservation for the gas, respectively, while the third equation
represents the motion of a particle i coupled with the gas. Here,
u is the velocity of the gas and I is the identity matrix. v is the
mass-weighted averaged particle velocity in the fluid element,
assuming that particles can be treated as fluid (Bai &
Stone 2010a). The terms of the right-hand side of
Equation (7) are radial tidal forces (gravity and centrifugal
force), vertical gravity, and the Coriolis force and the feedback
from the particle to the gas. In Equation (8), the first term on
the right-hand side describes a constant radial force due to gas
drag. vi is the particle velocity, x̂ and ẑ represent the unit
vectors in the x- and z-directions, and xi and zi are the values of
x and z for the particle i. u is the gas velocity interpolated from
the grid cell centers to the location of the particle. The second,
third, and fourth terms are radial and vertical tidal forces and
the Coriolis force. ag is the acceleration due to self-gravity,
which is considered only in 3D simulations. ag = −∇Φp, where
Φp is the gravitational potential and satisfies the Poisson’s
equation, ∇2
Φp = 4πGρp. To reduce computational time, the
particles are organized into “superparticles,” each representing
a cluster of individual particles of the same size. In the code
units, we normalize Ω = cs = H = ρg,0 = 1. The gas and
particle initial distributions are described in Equation (4) and as
⎛
⎝
⎜
⎞
⎠
⎟ ( )
H
z
H
2
exp
2
, 9
p
p
p
2
p
2
r
p
=
S
-
respectively, where Hp is the scale height of the particles and is
set to 0.02H (Bai & Stone 2010a) and Σp is the particle surface
density. The system uses an isothermal equation of state
P c
g s
2
r
= , and the particles are distributed uniformly in the x-
and y-directions and normally in the z-direction (Equation (9)).
The computational domains are 0.2H × 0.2H in 2D simulations
and 0.2H × 0.2H × 0.2H in 3D simulations, with all-periodic
boundary conditions. The resolution for 2D is 512 × 512, and
each grid cell has one particle (512 × 512 = 262,144 particles).
The resolution for 3D is 10 × 1283
(= 21 million particles in
total). Previous work shows that the resolution of 1283
can
produce the large clump size distribution well compared to
those with 2563
and 5123
(Simon et al. 2016); therefore, this
resolution of our 3D simulations is sufficient to resolve large
clumps, which are the focus of this work. The initial particle
size is assumed to be constant. Variable initial particle sizes
could affect the growth speed (Krapp et al. 2019) and the
concentrations of particles depending on their sizes (Yang &
Zhu 2021), but it is not known to affect the largest clump size
formed by streaming instability.
2.2. Clump Detection in 2D and 3D
After we conduct 2D simulations, we identify filaments
forming in the disk in order to constrain the parameter space
favorable for the streaming instability. We use the Kolmo-
gorov–Smirnov (KS) method, which is based on previous work
(Carrera et al. 2015). For each time step in a given time
window (25/Ω in our case), we compute the particle surface
density ( ) ( )
x x z dz
,
p p
0.1
0.1
ò r
S =
-
and average it over the time
window to give p t
áS ñ . We then sort the values in p t
áS ñ from
highest to lowest and compute the cumulative distribution
of this sorted data set. Finally, we use the KS test to output a
4

p-value that measures the likelihood that the underlying
cumulative distribution was linear:
( ) ( )
( ) ( )
Q z e
p Q D n
2 1
, 10
j
j z j
1
1 2 2 2
å
= -
=
=
¥
- -
where D is the maximum distance between the data and linear
cumulative distributions and n is the number of data points. The
higher the p is, the more homogeneous the system is, and
therefore filament formation is unlikely, whereas a small p
indicates filament formation is more likely. Here we assume
that filament formation is very likely when p < 0.10, the same
as in previous work (Carrera et al. 2015).
For self-gravitating clump detection in our 3D simulations,
we use PLanetesimal ANalyzer (PLAN; Li et al. 2019). This is
a tool to identify self-gravitating clumps specifically made for
Athena output. The density of each particle is assessed based
on the nearest particles. Particles with densities higher than a
threshold are associated with neighboring dense particles until
a density peak is achieved. In contrast, particle groups with a
saddle point less than a threshold remain separated. Detailed
descriptions are found in previous work (Li et al. 2019).
2.3. Model Parameters
The main parameters for the Athena simulations are (1) the
dimensionless stopping time τf (the Newton regime; see
Section 1.2 and Nakajima et al. 2022), (2) the normalized
pressure parameter Δ = ηvK/cs, (3) the ratio of the particle surface
density to the gas surface density Z = Σp/Σg (VMF = Z
1
1 +
), and
(4) the normalized gravity G̃ G
4 g,0
2
p r
º W for 3D simulations,
where G is the gravitational constant. A small value of τf indicates
a small particle radius Rp, where the particle is well coupled with
the gas, whereas a large value of τf corresponds to a large value of
Rp, which is more decoupled from the gas. A large value of Δ
corresponds to a large pressure gradient and quicker radial infall.
G̃ represents the strength of the self-gravity. The parameter space
we are exploring for our 2D simulations is τf = 10−3
, 10−2
, 10−1
,
100
, 101
; Δ = 0.1, 0.2, 0.3, 0.4, 0.5; and Z = 0.05 and 0.1. For 3D
simulations, we use Z = 0.1 and 0.3 and G̃ 0.1788
= and 0.5898,
where the former G̃ value corresponds to a slightly cooler
temperature (4700 K), while the latter corresponds to a hotter
temperature (5200 K). Here, we justify the choice of these
parameters. This range of τf corresponds to the particle radii
between 2 m and 20 km (see Equation (5)). The global disk
structures have been calculated based on hydrodynamic simula-
tions in previous work. The overall disk mass of the Moon-
forming disk is typically a few percent of Earth, depending on the
impact model (MD/ML = 1.35–2.80, where MD is the disk mass
and ML is the lunar mass; Nakajima & Stevenson 2014). The
midplane disk temperature ranges from 3000 K to 7000 K, and the
radial range of the disk is ∼1–8 R⊕ (see Figure 5 in Nakajima &
Stevenson 2014). The disk temperature is ∼4000–5500 K at r = 3
R⊕. For the general rocky and icy impact-induced disks, the disk
temperature can vary, typically in the range of thousands of K for
vapor-rich disks (Nakajima et al. 2022). The pressure gradient can
vary, but the typical value of η is ∼0.02–0.06 based on impact
simulations (Nakajima et al. 2022). In the Moon-forming disk, the
value of Z increases as the disk cools (on a timescale of
10–102
yr). In other words, Z can initially be zero in an energetic
Moon-forming impact model (e.g., Lock et al. 2018) and
eventually becomes infinity as the disk materials condense and
the gas disappears. Thus, picking a value of Z means that we are
seeing physics at a specific time. Since we are primarily interested
in high-vapor disks (i.e., an early phase of the disk), we explore
the mass ratio range Z ä [0.05, 0.1] for our 2D simulations and
Z ä [0.1, 0.3] for our 3D simulations. We use the large value of
Z = 0.3 as a sensitivity test. Higher values of Z can be achieved as
the disk cools, but we focus on small values of Z („0.3) for
several reasons. First, at a larger value of Z (>0.3), the
conventional GI in the liquid part of the disk can form moonlets.
At Z 0.76
Z
1
1
=
+
, and the total (Σp + Σg) surface density at
r ∼ 3 R⊕ is ∼108
kg m−2
in energetic models, which means that
Σp = 0.76 × 108
kg m−2
= 7.6 × 107
kg m−2
. Previous work
suggests that when the liquid (melt) layer’s thickness reaches
5–10 km (or the equivalent surface density of a few 107
kg m−2
),
GI can happen in the melt layer (Machida & Abe 2004). Whether
streaming instability occurs at the same time or whether
streaming instability affects the GI in the Moon-forming disk
are unknown and have not been explored in previous studies.
Under these circumstances, the important of streaming instability
becomes unclear. Additionally, these streaming instability
simulations have been conducted at low-Z values („0.1) in
previous work to reproduce conditions in the protoplanetary
disks (e.g., Abod et al. 2019; Li & Youdin 2021), and
simulations with high-Z (>0.3) values have not been fully
tested. For these reasons, we focus on relatively small values of Z
in this study.
For the set of Athena simulations, we focus on reproducing
two sets of the Moon-forming disk thermal profile; assuming
M = M⊕ and r = 3 R⊕, where M⊕ is the Earth mass and R⊕ is
the Earth radius, the gas pressure at the midplane around 3 R⊕
is ≈12 MPa, T = 4700 K, and mm = 30 g mol−1
. This makes
the density ρg,0 = 12.13 kg m−3
for an ideal gas, cs = 1140 m
s−1
, vK = 4562 m s−1
, Ω = 2.38 × 10−4
s−1
, the scale height
H = cs/Ω = 4.78 × 106
m, and G̃ 0.1788
= . For the higher-
temperature scenario, T = 5200 K, ρg,0 = 40.01 kg m−3
,
cs = 1200 m s−1
, H = 5.03 × 106
, and G̃ 0.5898
= . These
temperatures are motivated by previous hydrodynamic calcula-
tions of the Moon-forming disk formation (Nakajima &
Stevenson 2014), and the other parameters are calculated
based on an equation of state of dunite assuming that the vapor
and liquid phases are in equilibrium (MANEOS; Thompson &
Lauson 1972; Melosh 2007). Since we are assuming that the
disk is in liquid–vapor equilibrium, a higher disk temperature
leads to a higher vapor density at the midplane. For example,
the gas density at liquid–vapor equilibrium is 150 kg m−3
at
6000 K and 1.8 kg m−3
at 4000 K for dunite, according to
MANEOS (Thompson & Lauson 1972; Melosh 2007). Assum-
ing η ∼ 0.02–0.06 for the Moon-forming disk, Δ ∼ 0.1–0.2.
However, higher values are possible (Nakajima et al. 2022);
therefore, we explore the range of Δ ∼ 0.1–0.5.
The parameter values of Δ and G̃ are significantly different
from values in the protoplanetary disk, where the typical values
used are Δ ∼ 10−3
–10−2
and G̃ 0.05
~ in the protoplanetary
disk (Carrera et al. 2015; Simon et al. 2016). We use a fixed
value of Z in hydrodynamic calculations because the
condensation timescale (∼years) is longer than the simulation
timescale (we run 2D simulations for ∼100 orbits, which
correspond to 123 hr. A steady state is reached by this time).
5

We also assume that the disk does not evolve on this short
timescale.
3. Results
3.1. 2D Athena Simulations
Figure 1(A) shows four examples of spacetime plots of our
2D simulations, where (1) (τf, Δ) = 10−3
, 0.3, (2) 10−2
, 0.1,
(3) 10−1
, 0.2, and (4) 10−1
, 0.2. Z = 0.1 for the four cases.
These simulations represent four characteristic regimes. The
horizontal axis is x normalized by the scale height H, and the
vertical axis is the number of orbits. The color shows particle
concentration. In case 1, τf is small, which means that gas and
particles are well coupled, and this does not lead to filament
formation. In case 2, after ∼20 orbits, filaments start forming,
and these filaments are stable during the rest of the simulations,
which indicate that streaming instability occurs with this
chosen parameter. A radially shearing periodic boundary
condition is used in the x-direction; therefore, the filament that
appears at x/H = −0.1 at ∼100 orbits reappears at x/H = 0.1.
Two distinct filaments form in this simulation. In case 3,
filaments form, but they are not stable because their radial
movements are high due to the high Δ value. Thus, this is not
an ideal parameter space for streaming instability. Similar
behaviors have been seen in simulations for the protoplanetary
disk with other parameter combinations (e.g., Δ = 0.05, τf = 3,
Z = 0.005; Carrera et al. 2015). In case 4, the filaments are not
as clearly defined as those in case 2, but a coherent filament
structure is found after several orbits. Thus, this is also
considered as a streaming instability regime.
Figure 1(B) shows the summary of our 2D simulations for
Z = 0.1 (top) and Z = 0.05 (bottom). We use the KS test to
identify the extent of particle concentration (see Section 2.2).
Concentration is measured by the p-value, and when the value of
p is small (p < 0.1), streaming instability is likely (Carrera et al.
2015). The color bar indicates the p-value. For Z = 0.1, p < 0.1
is achieved when Δ = 0.1 and 10−3
„ τf „ 100
. For large Δ
values, the radial motions are large, and stable filaments do not
form, as discussed above. For Z = 0.05, this trend remains the
same, but the streaming instability regime is smaller
(10−2
„ τf „ 100
) than that for Z = 0.1. This is because the
larger Z leads to the presence of more particles in the disk and
therefore to more filament formation (Carrera et al. 2015). Thus,
our 2D simulations show that the most filaments form at
relatively small Δ (Δ = 0.1) and with both Z = 0.1 and Z = 0.05,
but the higher Z has more favorable conditions. This general
trend is consistent with previous work in the protoplanetary disk
(e.g., Carrera et al. 2015), which investigates smaller Δ values
for the protoplanetary disk (Δ „ 0.05).
3.2. 3D Athena Simulations
Now that we have identified the parameter space for the
streaming instability (Δ = 0.1, Z = 0.1), we perform 3D
simulations with self-gravity to estimate the size distributions
of streaming-instability-induced moonlets. Figure 2 shows
Figure 1. (A) Spacetime diagram for the four cases. The input parameters are (1) (τf, Δ) = 10−3
, 0.3, (2) 10−2
, 0.1, (3) 10−1
, 0.2, and (4) 1, 0.1, all at Z = 0.1. The
horizontal axis is x normalized by the scale height H. The vertical axis indicates the number of orbits. The color shows p p
S áS ñ, where Σp is the particle surface
density and p
áS ñ is the average along the x-axis. Filament formation occurs in cases 2 and 4, while no filament formation occurs in cases 1 and 3. (B) Summary of
results for Z = 0.1 and Z = 0.05. The colors show the p-value, and the clumping regime (p < 0.1) is shown in the sky blue shade. Parameters for cases 1–4 are
indicated. This shows that clumping occurs only at small Δ values (Δ = 0.1).
6

snapshots of one of our 3D simulations ( G̃
1,
f
t = =
Z
0.5898, 0.1
= , and Δ = 0.1) at four different times
(t = 0.32, 2.87, 3.39, and 3.55), where t is the number of
orbits. The self-gravity is not initially included and is turned on
at t = 3.18. This is a general practice to avoid artificial
clumping before streaming instability takes place (Simon et al.
2016). The color shows the particle density. The circles
indicate the Hill spheres of self-bound clumps, which are
identified using PLAN (Li et al. 2019). At t = 0.32, no clump is
identified, but by t = 2.87, the streaming instability fully
develops and reaches its steady state. After self-gravity is
turned on at t = 3.18, self-gravitating clumps form right away
(t = 3.39). This general behavior is similar to previous work on
the streaming instability in the protoplanetary disk (Abod et al.
2019), but the time it takes to form clumps by streaming
instability is shorter, probably because of the large Z value
(Simon et al. 2022) compared to those of the protoplane-
tary disk.
Figure 3(A) shows the cumulative distribution of moonlet
mass Mp, normalized by the characteristic self-gravitating mass
MG, defined as (Abod et al. 2019)
⎛
⎝
⎞
⎠
˜ ( ) ( )
M
G
Z G H
2
4
2
2
, 11
G
G
p
p
2
5
2 3
4
9 2 3 2
g,0
3
p
l
p
p r
= S =
S
W
=
where λG is an instability wavelength, which originates from
the Toomre dispersion relation, equating the tidal and
gravitational forces (Abod et al. 2019). The parameters are
τf = 0.1, 1, G̃ 0.1788, 0.5898
= , and Z = 0.1, 0.3. For the
G̃ 0.1787
= cases, MG = 5.19 × 1018
kg, and for the
G̃ 0.5898
= cases, MG = 2.17 × 1020
kg. For all 3D simula-
tions, Δ = 0.1 is assumed. The solid lines indicate the moonlet
mass distribution shortly after the onset of the self-gravity
(t = 3.66 for G̃ 0.1788, 0.1
f
t
= = , and t = 3.34 for all the
other cases). The dashed–dotted lines indicate the same
parameter after an additional time (t + Δt, where Δt = 1/
Ω = 0.16 orbits, except the Z = 0.3 case, where t
1
2
D = W
). The
Figure 2. Snapshots of a 3D simulation (τf = 1, Δ = 0.1, Z = 0.1, G̃ 0.5898
= ). The horizontal and vertical axes are x and y, normalized by the scale height H. The
color indicates the particle density, normalized by the average density. t indicates the number of orbits. The circles represent the Hill radius of each moonlet formed by
streaming instability. At first (t = 0.32), no concentration of particles is observed, but streaming instability clearly develops by t = 2.87. After self-gravity is turned on
at t = 3.18, moonlets form by GI (t = 3.39, 3.55). The self-gravitating clumps are identified using PLAN (see the main text).
7

reason for the shorter Δt for the higher-Z case is that the
streaming instability develops faster for higher-Z values (Simon
et al. 2022). The maximum Mp/MG = 0.254 (G̃ 0.1788,
=
1
f
t = at t = 3.50). This is broadly consistent with previous
studies that suggest that the maximum clump masses formed by
the streaming instability are characterized by MG in the range of
10−1
−101
MG (Abod et al. 2019; Li et al. 2019). Our result lies
on the lower end of this mass range. This indicates that the
moonlet mass distribution based on our numerical simulations
is consistent with the analytical mass model generated for the
protoplanetary disk. Thus, this also indicates that the processes
of streaming instability are similar between the Moon-forming
disk and protoplanetary disk despite the different input
parameters. We also ﬁnd that higher Z does not necessarily
lead to higher Mp.
Here we consider the best-case scenario of forming a large
moonlet. Assuming that the density of the clump is 3000 kg m−3
,
MG = 2.17 × 1020
kg is equivalent to 258 km in radius and 0.254
MG = 163 km. Based on Equation (1), the residence times of
258 km and 163 km moonlets are 92 and 58 days, respectively,
assuming vr,g = 0. These timescales are much shorter than the
Moon-formation timescale (10 s–100 yr, Thompson & Steven-
son 1988). Thus, even though the streaming instability can occur
in an initially vapor-rich Moon-forming disks, it does not help
increase the residence time of moonlets.
3.3. Exomoon Formation by Streaming Instability
The streaming instability likely occurs in impact-induced
moon-forming disks in extrasolar systems. The pressure gradient
parameter η (∼0.02–0.06) is similar regardless of the composition
of the disk (Nakajima et al. 2022). Figure 3(B) shows the disk
residence time for a clump with mass MG in a moon-forming disk
as a function of the planetary mass (Mplanet = 1–6 M⊕). “Rocky
planet” corresponds to disks formed by a collision between
terrestrial planets, while “icy planet” indicates disks formed by a
collision between icy planets whose mantles are made of water ice
(70 wt%) and cores are made of iron (30 wt%). This is produced
by calculating r/vr (see Equation (1)), where for rocky planets, we
use ρg,0 = 40.01 kg m−3
, T = 5200 K, ρp = 3000 kg m−3
, and
r = 3 Rplanet, where Rplanet is the planetary radius, and for icy
planets we use ρg,0 = 10.0 kg m−3
, T = 2000 K, ρp =
1000 kg m−3
, and r = 3 Rplanet. These are the temperatures when
the impact-induced disks reach complete vaporization based on
impact simulations (Nakajima et al. 2022). A higher temperature
is possible, but the disk would need to cool down to this
temperature to form particles (dust); therefore, this is the most
relevant temperature to assess the effect of streaming instability.
The thermal state of the disk formed by icy planetary collisions is
estimated based on an equation of state of water (Senft &
Stewart 2008). Additionally, ( )
R R M M
planet planet
1 3.7
= Å Å for
rocky planets and ( )
R R M M
1.2
planet planet
1 3
= Å Å for icy planets
are assumed (Kipping et al. 2013; Mordasini et al. 2012). At
Mplanet = 1–6 M⊕, these parameters make G̃ 0.42 0.67
= - for
rocky planets and G̃ 0.25
= for icy planets, which are similar to
the ranges covered in our hydrodynamic simulations (Section 2.3).
In both icy and rocky planet cases, the disk residence time is
a few tens of days to several months, which is short compared
to the satellite formation timescale (10 s–100 s yr, Nakajima
et al. 2022). This suggests that the streaming instability likely
plays a limited role in impact-induced moon-forming disks.
4. Discussion
4.1. Streaming Instability in the Moon-forming Disk
Figure 4 shows a schematic view of our hypothesis. An
energetic impact would generate a vapor-rich disk (the disk is
made of silicate vapor for the Moon-forming disk and rocky
moon-forming disks and water vapor for icy moon-forming
disks). Over time, the disk cools by radiation and small droplets
(<cm) emerge. These small droplets would grow by accretion
and streaming instability. However, once these moonlets reach
100 m–100 km in size, gas drag from the vapor is so strong that
they fall onto the planet on a short timescale (days–weeks).
This continues to occur until the VMF of the disk decreases so
that the gas drag effect is no longer strong (we assume the
silicate vapor is in equilibrium with the silicate liquid. As the
disk temperature decreases, the gas density decreases). Once
this condition is reached, a liquid layer emerges, and moonlets
can stay in the disk. However, by this time, a signiﬁcant disk
mass could have been lost (>80 wt%; Ida et al. 2020; Nakajima
et al. 2022). For this reason, only a small moon (or moons)
Figure 3. (A) Cumulative mass distribution of moonlets formed by streaming instability for Δ = 0.1, Z = 0.1, 0.3. The purple, dark blue, light blue, dark yellow, and
light yellow lines indicate parameter values of ( ˜ )
G Z
, ,
f
t = (0.5898, 0.1, 0.1), (0.1788, 0.1, 0.1), (0.5898, 1, 0.1), (0.1788, 1, 0.1), and (0.1788, 0.1, 0.3). The solid and
dashed–dotted lines represent the mass distribution at different times (see Section 3.2). The horizontal axis indicates the moonlet mass (Mp) normalized by MG, and the
vertical axis indicates the number of moonlets whose masses are larger than the given moonlet mass. (B) Residence time of a moonlet whose mass is MG formed by
streaming instability at r = 3 Rplanet as a function of the planetary mass Mplanet normalized by the Earth mass. The blue solid and yellow dashed–dotted lines represent
disks formed by collisions between icy planets and rocky planets, respectively.
8

could form from an initially vapor-rich disk. Thus, we suggest
that an initially vapor-rich disk is not suitable for forming our
Moon, and our result supports the hypothesis that the Moon
formed from an initially vapor-poor disk, including the
canonical model where the proto-Earth was hit by a Mars-
sized impactor (Canup & Asphaug 2001).
The isotopic composition of the Moon may constrain whether
streaming instability played a role during the Moon formation.
Some may argue that the observational fact that the Moon is
depleted in volatiles may be explained by accreting moonlets
formed by streaming instability before volatiles accreted onto the
Moon. However, this process would make the Moon isotopically
light if these isotopes experience kinetic fractionation (Dauphas
et al. 2015, 2022), which is inconsistent with the observation of
the enrichment of heavy isotopes in some elements in the Moon
(such as K; Wang & Jacobsen 2016). The lunar isotopes would be
heavier if they experience equilibrium fractionation, but the
equilibrium fractionation would not produce the observed isotopic
fractionations under the high disk temperature condition. Alter-
natively, some volatiles could be lost from the lunar magma ocean
under an equilibrium condition (Charnoz et al. 2021), but the
efficiency of the volatile loss is not fully known (Dauphas et al.
2022).
4.2. Streaming Instability in General Impact-induced Disks
Our model supports the previous work that suggests that
relatively large rocky (>6 M⊕, >1.6 R⊕) and icy (>1 M⊕, >
1.3 R⊕) planets cannot form impact-induced moons that are
large compared to the host planets (Nakajima et al. 2022).
Larger planets than those thresholds generate completely vapor
disks, because the kinetic energy involved in an impact scales
with the planetary mass. Thus, these large planets are not
capable of forming moons that are large compared to their host
planets. Moons can form by mechanisms other than impacts,
such as formation in a circumplanetary disk and gravitational
capture, but these moons tend to be small compared to the sizes
of their host planets (the predicted moon-to-planet mass ratio is
∼10−4
for moons formed by circumplanetary disk, Canup &
Ward 2006; and gravitationally captured moons are small in the
solar system, Agnor & Hamilton 2006). Thus, fractionally large
moons compared to the host planet sizes, which are
observationally favorable, likely form by impact. So far, no
exomoon has been confirmed despite extensive searches, but
future observations, especially with the James Webb Space
Telescope (Cassese & Kipping 2022), may be able to find
exomoons and test this theoretical hypothesis.
4.3. Comparison between the Moon-forming Disk and
Protoplanetary Disk
It is certainly intriguing that the streaming instability is able to
solve the gas drag problem in the protoplanetary disk but not in
the Moon-forming disk. In both scenarios, the clump sizes
formed by the streaming instability happen to be similar
(MG ∼ 100 km, see Johansen et al. 2015, for the protoplanetary
disk and ∼100 km for the protolunar disk, respectively). This
size is ∼105
times larger than the size of the particle (0.52–1 m;
Adachi et al. 1976) that experiences the strongest gas drag
(τf = 1 and vr ∼ −ηvK) in the protoplanetary disk. Therefore,
once particles become ∼100 km-sized planetesimals by
streaming instability, the radial velocities of the planetesimals
decrease drastically (vr
v
2 K
f
= -
h
t
for large τf; see Equation (A7)),
which helps planetesimal growth. In contrast, the largest moonlet
size (100 km) formed by streaming instability in the Moon-
forming disk is only 50 times larger than the particle size that
experiences the strongest gas drag (2 km). This does not result in
a large τf change or vr change; therefore, streaming instability
does not effectively help moon formation. Therefore, streaming
instability is an effective mechanism to skip the “1 m barrier” in
Figure 4. Schematic view of Moon formation from an initially vapor-rich disk.
9

the protoplanetary disk, whereas it is not an effective mechanism
to skip the “kilometer barrier” in the Moon-forming disk.
4.4. Streaming Instability in the Circumplanetary Disk around
a Gas Giant
The streaming instability may occur during moon formation in
various satellite systems around gas giants. The Galilean moons
around Jupiter and Titan around Saturn likely formed from their
circumplanetary disks. In these disks, the particle-to-gas ratio Z is
typically 10−4
–10−2
, which was thought to be too small for
streaming instability to occur (Carrera et al. 2015; Yang et al.
2017; Shibaike et al. 2017). However, recent streaming instability
calculations show that streaming instability can occur as low as
Z = 4 × 10−3
(at Δ = 0.05 and τf = 0.3; Li & Youdin 2021), and
the required Z depends on the disk conditions (e.g., Carrera et al.
2015; Yang et al. 2017; Sekiya & Onishi 2018). This may
indicate that streaming instability can be potentially important in
circumplanetary disks around gas giants if the dust–gas ratio of
the disk is relatively high. The velocity difference between the gas
and particles vrel is (Canup & Ward 2002)
⎛
⎝
⎞
⎠
( )
v v c
c
r
. 12
rel k s
s
h
= ~
W
Assuming 0.1
c
r
c
vk
s s
= ~
W
(Canup & Ward 2002), this yields
η ∼ 0.01 and 0.1
v
c
k
s
D = ~
h
. If Z is sufficiently large (at least
Z > 4 × 10−3
; Li & Youdin 2021), it is possible that the
streaming instability takes place in a circumplanetary disk
around a gas giant. A rough estimate of the size of a moonlet in
a circumplanetary disk is MG = 1.4 × 1016
kg (Equation (11)),
which is ∼10.3 km in radius if it is a rocky moonlet. Here, we
are assuming Z = 10−2
, Σp = 3 × 104
kg m−2
, cs = 1000 m s−1
,
H ∼ cs/Ω, and ( )
H 2
g g
r p
= S (Canup & Ward 2002). These
moonlets are relatively small compared to large moons around
gas giants (e.g., the Galilean moon masses are in the range of
1022
–1023
kg), and it is unclear if they would significantly impact
these moon-formation processes. Further research is needed to
understand their effects on the satellite formation in a
circumplanetary disk around a gas giant.
4.5. Model Limitations
There are several model limitations that need to be addressed in
our future work. First, the effect of the Roche radius (aR ∼ 3 R⊕)
needs to be taken into account to understand the Moon's
formation. An inspiraling moonlet would not directly reach the
Earth but would be tidally disrupted near the Roche radius, where
it then might be incorporated into an accretion disk (e.g., Salmon
& Canup 2014). Second, our model presented here does not take
into account the evolution of the disk in detail, which is important
to identify the final mass and composition of the resulting moon.
As the disk spreads out, it is likely that Δ decreases, which would
slow down the radial drift of moonlets. This gas drag effect would
disappear once the local vapor condenses, which would occur at
the outer part of the disk first, since that part of the disk can cool
efficiently due to its large surface area. Given that such moonlets
that directly form from the disk would not have time to lose
volatiles in the disk phase, this model requires the Moon to have
lost its volatiles before or after the disk phase, for example, during
the lunar magma ocean phase (Charnoz et al. 2021; see further
discussion in Section 4.1). Nevertheless, such a scenario may be
possible in moon-forming disks in the solar and extrasolar
systems. The effects of the Roche radius and disk evolution will
be addressed in our future work.
5. Conclusions
In conclusion, we conduct hydrodynamic simulations using
Athena and show that the streaming instability can form self-
gravitating clumps (∼102
km) in a vapor-rich moon-forming disk
generated by a giant impact, but the sizes of the clumps it
generates are not large enough to avoid inspiraling due to the
strong gas drag. This is a major difference from processes in the
protoplanetary disk, where the streaming instability can efficiently
form large clumps to avoid the strong gas drag effect. As a result,
growing moonlets in an initially vapor-rich moon-forming disk
continue to fall onto the planet once they reach sizes of 100 m–
100 km. These moonlets could grow further once the disk cools
enough and the VMF of the disk becomes small. However, by this
time, a significant amount of the disk mass is lost, and the
remaining disk could make only a small moon. This result is
applicable to impact-induced moon-forming disks in the solar
system and beyond; we find that the streaming instability is not an
efficient mechanism to form a large moon from an impact-induced
vapor-rich disk in general. As a result, we support previous work
that suggests that fractionally large moons compared to their host
planets form from vapor-poor disks. This supports the Moon-
formation models that produce vapor-poor disks, such as the
canonical model. For exomoons, our work supports the previous
work that suggests that the ideal planetary radii that host
fractionally large moons are 1.3–1.6 R⊕ (Nakajima et al. 2022)
given that rocky or icy planets larger than these sizes would likely
produce completely vapor disks, which are not capable of forming
large moons (Nakajima et al. 2022). The streaming instability may
take place in circumplanetary disks, but their effect on the moon-
formation process needs further investigation.
Acknowledgments
We thank the code developers of Athena (Stone et al. 2008;
Bai & Stone 2010a; Simon et al. 2016) and PLAN (Li et al.
2019). We appreciate discussions with Shigeru Ida and Scott D.
Hull. J.A. was partially supported by the Research Experience
for Undergraduates (REU) Program, National Science Founda-
tion (NSF), under grant No. PHY-1757062. M.N. was supported
in part by the National Aeronautics and Space Administration
(NASA) grant Nos. 80NSSC19K0514, 80NSSC21K1184, and
80NSSC22K0107. Partial funding for M.N. was also provided
by NSF EAR-2237730 as well as the Center for Matter at
Atomic Pressures (CMAP), an NSF Physics Frontier Center,
under award PHY-2020249. Any opinions, findings, conclu-
sions, or recommendations expressed in this material are those of
the authors and do not necessarily reflect those of the National
Science Foundation. M.N. was also supported in part by the
Alfred P. Sloan Foundation under grant G202114194.
Appendix
Appendix A
The momentum equation for gas in the radial direction is
written as
( )
v
r
GM
r
dp
dr
1
, A1
2
2
g
g
r
= +
f Å
10

where vf is the azimuthal velocity of the gas. Using
Equation (2) in the main text and the relationship of p c
g g 2
2
r
= ,
( ) ( )
v v 1 2 . A2
K
1
2
h
= -
f
The equations of motion of the particles in the gas in the radial
and azimuthal directions, vr and vf, are (Armitage 2010;
Takeuchi & Lin 2002)
( ) ( )
dv
dt
v
r
r
t
v v
1
, A3
r
r r
2
2
f
,g
= - W - -
f
( ) ( ) ( )
d
dt
rv
r
t
v v , A4
f
,g
= - -
f f f
where vr,g and vf,g are the gas velocities in the radial and
azimuthal directions. Assuming that the
d
dt
terms are negligible,
one finds that
( )
v
v v
2
. A5
r
r
f
1
,g K
f f
1
t h
t t
=
-
+
-
-
This formulation is slightly different from previous formula-
tions (−2ηvK instead of −ηvK; Armitage 2010; Takeuchi &
Lin 2002) due to the different definition of η by a factor of 2.
This leads to
( )
v
v v
2
. A6
r
r
f
1
,g K
f f
1
t h
t t
=
-
+
-
-
This can be further simplified by the values of τf as
⎧
⎨
⎪
⎪
⎩
⎪
⎪


( )
( )
v
v v
v v
v
2 at 1,
1
2
2 at 1,
2
at 1.
A7
r
r
r
K
,g f K f
,g K f
f
f
ht t
h t
h
t
t
=
-
- ~
-
For a steady disk flow, one could assume vr r
,g
3
2
~ -
n
(Armitage 2010). In the main text, we simply assume vr,g = 0.
Appendix B
An example input file for a 2D Athena simulation is listed in
the supplementary material (athinput.par_strat2d). The resolu-
tion in this file is 512 × 512, and τf = 1, Z = 0.1, and Δ = 0.1.
ORCID iDs
Miki Nakajima https:/
/orcid.org/0000-0001-5014-0448
Jeremy Atkins https:/
/orcid.org/0009-0005-0722-8408
Jacob B. Simon https:/
/orcid.org/0000-0002-3771-8054
Alice C. Quillen https:/
/orcid.org/0000-0003-1280-2054
References
Abod, C. P., Simon, J. B., Li, R., et al. 2019, AJ, 883, 192
Adachi, I., Hayashi, C., & Nakazawa, K. 1976, PThPh, 56, 1756
Agnor, C. B., & Hamilton, D. P. 2006, Natur, 441, 192
Armitage, P. J. 2010, Astrophysics of Planet Formation (Cambridge:
Cambridge Univ. Press)
Armytage, R. M. G., Georg, R. B., Williams, H. M., & Halliday, A. N. 2012,
GeCoA, 77, 504
Asphaug, E., Emsenhuber, A., Cambioni, S., Gabriel, T. S. J., &
Schwartz, S. R. 2021, PSJ, 2, 200
Bai, X.-N., & Stone, J. M. 2010a, ApJ, 722, 1437
Bai, X.-N., & Stone, J. M. 2010b, ApJS, 190, 297
Bonomo, A. S., Zeng, L., Damasso, M., et al. 2019, NatAs, 3, 416
Bouhifd, A. M., Jephcoat, A. P., Porcelli, D., Kelley, S. P., & Marty, B. 2020,
GChPL, 15
Cameron, A. G. W., & Ward, W. R. 1976, Lunar Planet. Sci. VII, 120
Cano, E. J., Sharp, Z. D., & Shearer, C. K. 2020, NatGe, 13, 270
Canup, R. M. 2004, Icar, 168, 433
Canup, R. M. 2005, Sci, 307, 546
Canup, R. M. 2012, Sci, 338, 1052
Canup, R. M., & Asphaug, E. 2001, Natur, 412, 708
Canup, R. M., Visscher, C., Salmon, J., & Fegley, B. 2015, NatGe, 8, 918
Canup, R. M., & Ward, W. R. 2002, AJ, 124, 3404
Canup, R. M., & Ward, W. R. 2006, Natur, 441, 834
Canup, R. M., Righter, K., Dauphas, N., et al. 2023, RvMG, 89, 53
Carrera, D., Johansen, A., & Davies, M. B. 2015, A&A, 579, A43
Cassese, B., & Kipping, D. 2022, MNRAS, 516, 3701
Charnoz, S., & Michaut, C. 2015, Icar, 260, 440
Charnoz, S., Sossi, P. A., Lee, Y.-N., et al. 2021, Icar, 364, 114451
Chiang, E., & Youdin, A. N. 2010, AREPS, 38, 493
Colella, P. 1990, JCoPh, 87, 171
Colella, P., & Woodward, P. R. 1984, JCoPh, 54, 174
Craddock, R. A. 2011, Icar, 211, 1150
Cúk, M., Hamilton, D. P., Lock, S. J., & Stewart, S. T. 2016, Natur, 539, 402
Ćuk, M., Lock, S. J., Stewart, S. T., & Hamilton, D. P. 2021, PSJ, 2, 147
Cúk, M., & Stewart, S. T. 2012, Sci, 338, 1047
Dauphas, N. 2017, Natur, 541, 521
Dauphas, N., Poitrasson, F., Burkhardt, C., Kobayashi, H., & Kurosawa, K.
2015, E&PSL, 427, 236
Dauphas, N., Nie, N. X., Blanchard, M., et al. 2022, PSJ, 3, 29
Halliday, A. N., & Canup, R. M. 2022, NRvEE, 4, 19
Halliday, N. H. 2004, Natur, 427, 505
Hartmann, W. K., & Davis, D. R. 1975, Icar, 24, 504
Hosono, N., Karato, S.-i., Makino, J., & Saitoh, T. R. 2019, NatGe, 12, 418
Hull, S. D., Nakajima, M., Hosono, N., Canup, R. M., & Gassmöller, R. 2024,
PSJ, 5, 9
Ida, S., Ueta, S., Sasaki, T., & Ishizawa, Y. 2020, NatAs, 4, 880
Johansen, A., Mac Low, M. M., Lacerda, P., & Bizzarro, M. 2015, SciA, 1,
e1500109
Johansen, A., Oishi, J. S., Mac Low, M. M., et al. 2007, Natur, 448, 1022
Kegerreis, J. A., Ruiz-Bonilla, S., Eke, V. R., et al. 2022, ApJL, 937, L40
Kenworthy, M., Lock, S., Kennedy, G., et al. 2023, Natur, 622, 251
Kipping, D. M., Forgan, D., Hartman, J., et al. 2013, ApJ, 777, 134
Krapp, L., Benítez-Llambay, P., Gressel, O., & Pessah, M. E. 2019, ApJL,
878, L30
Kruijer, T. S., Archer, G. J., & Kleine, T. 2021, NatGe, 14, 714
Kruijer, T. S., Kleine, T., Fischer-Godde, M., & Sprung, P. 2015, Natur,
520, 534
Li, R., & Youdin, A. N. 2021, ApJ, 919, 107
Li, R., Youdin, A. N., & Simon, J. B. 2019, ApJ, 885, 69
Lock, S. J., Stewart, S. T., Petaev, M. I., et al. 2018, JGRE, 123, 910
Machida, R., & Abe, Y. 2004, ApJ, 617, 633
Melosh, H. J. 2007, M&PS, 42, 2079
Meng, H. Y. A., Su, K. Y. L., Rieke, G. H., et al. 2014, Sci, 345, 1032
Mordasini, C., Alibert, Y., Georgy, C., et al. 2012, A&A, 547, 112
Mullen, P. D., & Gammie, C. F. 2020, ApJL, 903, L15
Nakajima, M., Genda, H., Asphaug, E., & Ida, S. 2022, NatCo, 13, 568
Nakajima, M., & Stevenson, D. J. 2014, Icar, 233, 259
Nakajima, M., & Stevenson, D. J. 2015, E&PSL, 427, 286
Nakajima, M., & Stevenson, D. J. 2018, E&PSL, 487, 117
Nie, N. X., & Dauphas, N. 2019, ApJL, 884, L48
Nielsen, S. G., Bekaert, D. V., & Auro, M. 2021, NatCo, 12, 1817
Pahlevan, K., & Stevenson, D. J. 2007, E&PSL, 262, 438
Rufu, R., Aharonson, O., & Perets, H. B. 2017, NatGe, 10, 89
Rufu, R., & Canup, R. M. 2020, JGRE, 125, e2019JE006312
Salmon, J., & Canup, R. M. 2012, ApJ, 760, 83
Salmon, J., & Canup, R. M. 2014, RSPTA, 372, 20130256
Sekiya, M., & Onishi, I. K. 2018, ApJ, 860, 140
Senft, L. E., & Stewart, S. T. 2008, M&PS, 43, 1993
Shibaike, Y., Okuzumi, S., Sasaki, T., & Ida, S. 2017, ApJ, 846, 81
Simon, J. B., Armitage, P. J., Li, R., & Youdin, A. N. 2016, ApJ, 822, 55
Simon, J. B., Blum, J., Birnstiel, T., & Nesvorný, D. 2022, Comets III, in press
(arXiv:2212.04509)
Slattery, W. L., Benz, W., & Cameron, A. G. W. 1992, Icar, 99, 167
Sossi, P. A., Moynier, F., & van Zuilen, K. 2018, PNAS, 115, 10920
Stone, J. M., & Gardiner, T. A. 2010, ApJS, 189, 142
Stone, J. M., Gardiner, T. A., Teuben, P., Hawley, J. F., & Simon, J. B. 2008,
ApJS, 178, 137
11

Takahashi, S. Z., & Inutsuka, S.-i. 2014, ApJ, 794, 55
Takeuchi, T., & Lin, D. N. C. 2002, ApJ, 581, 1344
Thiemens, M. M., Sprung, P., Fonseca, R. O. C., Leitzke, F. P., & Münker, C.
2019, NatGe, 12, 696
Thiemens, M. M., Tusch, J., O. C. Fonseca, R., et al. 2021, NatGe, 14,
716
Thompson, C., & Stevenson, D. J. 1988, ApJ, 333, 452
Thompson, M. A., Weinberger, A. J., Keller, L. D., Arnold, J. A., &
Stark, C. C. 2019, ApJ, 875, 45
Thompson, S. L., & Lauson, H. S. 1972, Improvements in the Chart D
Radiation-Hydrodynamic Code. III: Revised Analytic Equations of State
Tech. Rep. SC-RR-71-0714, Sandia National Laboratory
Tominaga, R. T., Inutsuka, S.-i., & Takahashi, S. Z. 2018, PASJ, 70, 3
Tominaga, R. T., Takahashi, S. Z., & Inutsuka, S.-I. 2019, ApJ, 881, 53
Toro, E. F. 1999, Riemann Solvers and Numerical Methods for Fluid
Dynamics: A Practical Introduction (Dordrecht: Springer)
Touboul, M., Puchtel, I. S., & Walker, R. J. 2015, Natur, 520, 530
Wang, K., & Jacobsen, S. B. 2016, Natur, 538, 487
Ward, W. R., Canup, R. M., & Rufu, R. 2020, JGR, 125, e2019JE00626
Weidenschilling, S. J. 1977, MNRAS, 180, 57
Wiechert, U., Halliday, A. N., Lee, D. C., et al. 2001, Sci, 294, 345
Williams, C. D., & Mukhopadhyay, S. 2019, Natur, 565, 78
Williams, C. D., Mukhopadhyay, S., Rudolph, M. L., & Romanowicz, B. 2019,
GGG, 20, 4130
Yang, C. C., Johansen, A., & Carrera, D. 2017, A&A, 606, 80
Yang, C.-C., & Zhu, Z. 2021, MNRAS, 508, 5538
Yokochi, R., & Marty, B. 2004, E&PSL, 225, 77
Youdin, A. N. 2011, ApJ, 731, 99
Youdin, A. N., & Goodman, J. 2005, ApJ, 620, 459
Young, E. D., Kohl, I. E., Warren, P. H., et al. 2016, Sci, 351, 493
Zhang, J., Dauphas, N. M. D. A., Leya, I., & Fedkin, A. 2012, NatGe,
1429, 251
12

The Limited Role of the Streaming Instability during Moon and Exomoon Formation

Related slideshows

Recommended for you

Recommended for you

Recommended for you

More Related Content

Similar to The Limited Role of the Streaming Instability during Moon and Exomoon Formation

Similar to The Limited Role of the Streaming Instability during Moon and Exomoon Formation (20)

More from Sérgio Sacani

More from Sérgio Sacani (20)

Recently uploaded

Recently uploaded (20)

The Limited Role of the Streaming Instability during Moon and Exomoon Formation