Spiral arm instability -- III. Fragmentation of primordial protostellar discs
MMon. Not. R. Astron. Soc. , 1– ?? (2014) Printed 22 October 2019 (MN L A TEX style file v2.2)
Spiral arm instability - III. Fragmentation of primordialprotostellar discs
Shigeki Inoue , , , (cid:63) & Naoki Yoshida , , Center for Computational Sciences, University of Tsukuba, Ten-nodai, 1-1-1 Tsukuba, Ibaraki 305-8577, Japan Chile Observatory, National Astronomical Observatory of Japan, Mitaka, Tokyo 181-8588, Japan Kavli Institute for the Physics and Mathematics of the Universe (WPI), UTIAS, The University of Tokyo, Chiba 277-8583, Japan Department of Physics, School of Science, The University of Tokyo, Bunkyo, Tokyo 113-0033, Japan Research Center for the Early Universe, School of Science, The University of Tokyo, Bunkyo, Tokyo 113-0033, Japan
22 October 2019
ABSTRACT
We study the gravitational instability and fragmentation of primordial protostellardiscs by using high-resolution cosmological hydrodynamics simulations. We follow theformation and evolution of spiral arms in protostellar discs, examine the dynamicalstability, and identify a physical mechanism of secondary protostar formation. Weuse linear perturbation theory based on the spiral-arm instability (SAI) analysis inour previous studies. We improve the analysis by incorporating the effects of finitethickness and shearing motion of arms, and derive the physical conditions for SAI inprotostellar discs. Our analysis predicts accurately the stability and the onset of armfragmentation that is determined by the balance between self-gravity and gas pressureplus the Coriolis force. Formation of secondary and multiple protostars in the discsis explained by the SAI, which is driven by self-gravity and thus can operate withoutrapid gas cooling. We can also predict the typical mass of the fragments, which isfound to be in good agreement with the actual masses of secondary protostars formedin the simulation.
Key words: instabilities – methods: numerical – methods: analytical – stars: Popu-lation III – stars: formation.
The first generations of stars affects the thermal and chem-ical state of the inter-galactic medium in the early Universeby photo-ionisation and heating, and by dispersing heavyelements synthesized within them. The first stars also affectthe subsequent formation and evolution of galaxies and mas-sive black holes through a variety of feedback effects (Bromm& Yoshida 2011). The initial stellar mass function is perhapsthe most important quantity that characterises the first starsand their overall impact to cosmic structure formation. Sofar, primordial star formation in a cosmological context hasbeen studied extensively, and detailed protostellar evolutioncalculations suggest that the first stars have a variety ofmasses with 10–1000 M (cid:12) (e.g. Hirano et al. 2014).Formation of binaries or multiple stars has been alsofound in recent simulations that follow the post-collapse evo-lution of primordial proto-stellar systems (Greif et al. 2012;Stacy et al. 2016; Susa 2019). A primordial proto-stellardisc is embedded within an accreting envelope (Hosokawa (cid:63)
E-mail: [email protected] et al. 2016), and often becomes unstable to fragment andyield multiple protostars, but the exact mechanism(s) thatdrives the disc fragmentation has not been identified. Thefrequency of disc fragmentation and the masses of the frag-ments are critically important because low-mass primordialstars, if they are formed, survive until the present day andmight be found in the local Universe.Greif et al. (2012) run cosmological simulations of pri-mordial star formation to show that multiple protostars areformed in the circumstellar disc around the main protostar.In their simulation, spiral arms are excited in the disc, whichoccasionally fragment to generate secondary protostars thathave small masses. If the secondary stars are gravitationallybound within the disc, and do not get accreted by the centralstar, they will form a binary or multiple stellar system. Itis important to understand the physical mechanism of discfragmentation, and it appears that fragmentation of spiralarms is a critical step for the formation of multiple systemsand star clusters.Fragmentation of a spiral arm is often discussed in termsof gravitational instability, such like non-linear evolution ofToomre-unstable regions (Safronov 1960; Toomre 1964). The c (cid:13) a r X i v : . [ a s t r o - ph . GA ] O c t S. Inoue& N. Yoshida
Toomre analysis itself describes the linear stability of a two-dimensional region within a thin disc, and it remains unclearwhat is the end product of Toomre instability – spiral armsor spherical clumps. Furthermore, the Toomre analysis isno longer applicable if there is significant density contrastbetween perturbations (arms) and the disc. Another impor-tant issue may be the effect of radiative cooling. It is oftenargued that rapid gas cooling is responsible for fragmen-tation of arms (e.g. Gammie 2001; Rice et al. 2003). Sincethe gas in primordial protostellar discs evolve approximatelyisothermally, it is not clear if fragmentation can be drivenby radiative cooling. Tsukamoto et al. (2015) and Takahashiet al. (2016) argue that rapid cooling is not essential fordisc fragmentation. While there are a variety of physicalprocesses suggested, and related criteria are proposed, thephysical mechanism of primordial gas disc fragmentation hasnot been fully explored.In our previous studies (Inoue & Yoshida 2018, 2019,hereafter Papers I and II), we present linear perturbationanalysis for spiral arms in galactic discs and derive the phys-ical conditions for spiral arm instability (SAI). We use theoutput of high-resolution simulations of isolated disc galax-ies to show that SAI induce fragmentation of galactic spiralarms. We also derive the characteristic mass of the fragmentsformed via SAI. In this
Letter , we apply the SAI analysis toprotostellar discs formed in early primordial gas clouds. Weinvestigate the physical mechanism of spiral-arm fragmen-tation through which secondary protostars are formed.
We use the cosmological hydrodynamics simulations of Greifet al. (2012), which are performed with the moving-mesh hy-drodynamics code
Arepo (Springel 2010). The simulationsfollow the collapse of primordial gas clouds in small-mass”minihaloes”. Primordial gas chemistry and associated cool-ing processes are followed for the species of H, H + , H − , H +2 ,H , He, He + , He ++ , D, D + , HD and free electrons. Details ofthe simulations are found in Greif et al. (2012). We use theirMH1 and MH4 runs. In the moving-mesh hydrodynamcissimulations, adoptive refinement is performed for gas cellssuch that a Jeans length is resolved with 32 cells. In Papers I and II, we derive the physical conditions for grav-itational instability of a spiral arm in a thin disc, based onlinear analysis for azimuthal perturbations along the armthat is approximated to be tightly wound (see also Taka-hashi et al. 2016). The instability parameter for perturba-tions with wavenumber k is given as S ≡ σ k + κ π G f ( kW )Υ k , (1)where G is the gravitational constant, W and Υ are the halfwidth and the line mass of the spiral arm defined as themass per unit length, and f ( kW ) ≡ [ K ( kW ) L − ( kW ) + K ( kW ) L ( kW )] with K i and L i denoting the modifiedBessel and Struve functions of order i . We compute theepicyclic frequency κ directly from the rotation velocity v φ as κ = 2 v φ R (cid:16) d v φ d R + v φ R (cid:17) , (2)as in Inoue et al. (2016). The velocity dispersion is calculatedas σ = c + σ φ , where c snd and σ φ are sound speed andthe azimuthal turbulent velocity dispersion. We note that σ , Υ, W and κ can be evaluated locally from the simulationoutput (see Section 3.2).The instability parameter S is expressed as a functionof k . If S ( k ) > k , the spiral arm is stable againstall perturbations. The instability condition for a spiral armand for perturbations with k is given by min[ S ( k )] <
1. Ourprevious analyses in Papers I, II and Takahashi et al. (2016)assume rigid rotation within the arm, i.e., κ = 2Ω, but herewe improve the analysis by considering the velocity shear. The above analysis using equation (1) still ignores the ver-tical thickness of a spiral arm. Since we consider spiral armsexcited in primordial protostellar discs, they may have alarge thickness to affect significantly the stability analysis. Athick arm has low mass-concentration, and the gravitationalforce on the disc plane is effectively reduced. Accordingly,using equation (1) overestimates the strength of the gravita-tional instability. To take this effect into account, we adopta correction proposed by Toomre (1964), in which the armhas a vertical height of h and density fluctuations within thearm are assumed to be independent of height from the discplane. In this case, fluctuations of gravitational potential onthe disc plane are weakened by a factor of [1 − exp( − kh )] /kh .With this thickness correction, equation (1) is modified toyield S T ≡ ( σ k + κ ) h π Gf ( kW )Υ k [1 − exp( − kh )] . (3)The instability condition for a spiral arm is given bymin[ S T ( k )] <
1. We define the vertical height of the armas h ≡ ( σ z + c ) / ( π G Σ), where σ z is the vertical veloc-ity dispersion (Elmegreen 2011). Although there are otherstudies that propose slightly different correction factors (e.g.Goldreich & Lynden-Bell 1965; Vandervoort 1970; Romeo1992), they give essentially the same instability parameter . To compute the instability parameter S T , we use verticallyaveraged values of the physical quantities such as Σ, σ φ , σ z , c snd and v φ . In practice, we first apply two-dimensionalGaussian smoothing with r smooth = 0 . R, φ ) for each snapshot.In order to examine the SAI, we first need to detectspiral arms and measure their half widths W . To this end, These analyses consider the instability of local disc regionsagainst radial perturbations, rather than that of a spiral armagainst azimuthal ones. However, the finite thickness only affectsthe estimate of the strength of gravity on the disc plane, and theinstability analysis does not differ between radial and azimuthalperturbations as long as the direction of the perturbations is in-dependent of time (Toomre 1964).c (cid:13) , 1– ?? AI in protostellar discs Figure 1.
The polar-map representation of our analysis for the run that shows arm fragmentation. The left and central panels show theface-on gas distributions in the Cartesian and polar coordinates. Times elapsed since the first snapshot (top panel) are indicated in theleft panels. The right panels show the thickness-corrected instability parameters, min( S T ), where unstable regions with min( S T ) < we perform one-dimensional Gaussian fitting along the ra-dial direction at a given φ using the polar map of thesurface density Σ( R, φ ). The fitting function is defined as˜Σ(
R, ξ, φ ) = Σ(
R, φ ) exp[ − ξ / w ], where ξ represents ra-dial offset from R . The fitting is performed in the range of − . w < ξ < . w with varying w . We search for w thatminimises the goodness-of-fit χ . Then the arm half widthis given by W = 1 . w . With this definition of W , the es-timated Gaussian surface density at the edge of the spiralarm corresponds to 30 percent of the peak value (Takahashiet al. 2016), i.e., ˜Σ( R, ± W, φ ) = 0 . R, φ ). The line-massof an arm is obtained as Υ(
R, φ ) = 1 . W Σ (see Paper I).We define spiral arms to be regions where log χ < − . S T is independent of this choice. Note thatthe above procedures assume a spiral arm to have a Gaus-sian density distribution and a pitch angle of θ = 0. OurPapers I and II confirm that this assumption neverthelessallows us to detect spiral arms robustly and to derive theinstability parameter accurately.We compute the epicyclic frequency κ in the followingmanner. The rotation velocity v φ ( R, φ ) is calculated fromthe smoothed velocity fields. The local velocity gradientd v φ / d R at ( R, φ ) is computed by linear least-square fittingfor v φ in the radial range of − W < ξ < W , and κ ( R, φ )is computed from equation (2). All the quantities necessaryfor computing S are obtained as functions of ( R, φ ). Spiral arms are excited in protostellar discs in all our simu-lations. Although some arms occasionally fragment to formsecondary protostars, not all the spiral arms do so. As illus-trative cases, we focus on two characteristic runs; one showsarm fragmentation and the other with stable arms. We ap-ply our analysis described in Section 3 and show that thefragmentation in the former case is explained by SAI.
Fig. 1 shows the results of our polar-map analysis for therun where two spiral arms fragment to yield secondary pro-tostars. The central star has a mass of 2 .
13 M (cid:12) within1 . We find that the spiral arms and the inter-armregions have similar temperatures T (cid:39) T > ∼ τ cool ∼ τ cool Ω ∼ At t = 2 .
68 yr, the spiral arms are located from R (cid:39) . (cid:13) , 1– ?? S. Inoue& N. Yoshida t = 0, when a pair of prominent spiral arms have developed.In the next few years, both the arms fragment to yield twosecondary stars. Arm A has initially a smooth density dis-tribution, and the density increases until t = 1 .
66 yr (thethird snapshot). Finally at t = 2 .
20 yr, Arm A fragmentsto form a secondary star. Similarly, Arm B bears anothersecondary star at t = 1 .
66 yr.The right panels of Fig. 1 indicate the instability param-eters min( S T ) in the spiral-arm regions with log χ < − . S T ) < t = 2 .
68 yr. Our SAI analysis accurately predicts the frag-mentation in the simulation, and can also quantify the de-gree of gravitational instability by the value of min( S T ). This suggests that the physical mechanism of the spiral-armfragmentation is linear dynamical instability driven by self-gravity of the arm.In Paper I, we successfully predict the characteristicmass of clumps formed by SAI. Here, we aim at predictingthe mass of fragments (protostars) in the disc by followingessentially the same procedure. The wavelength of the mostunstable perturbation can be obtained as λ MU = 2 π /k MU where k MU is the wavenumber that gives the minimum S T ( k ). In the snapshot at t = 0 in Fig. 1, for both thetwo arms, we find that λ MU (cid:39) W (cid:39) . λ MU > ∼ W ), the most un-stable perturbation with λ MU is expected to collapse in thedirection along the spiral arm. Then, the mass of a collaps-ing clump is estimated to be M cl ∼ Υ λ MU . The masses ofthe fragments (secondary stars) can be already predicted at t = 0 to be M cl (cid:39) . . (cid:12) for Arms A and B, respec-tively. In the snapshot at t = 2 .
68 yr (Fig. 1), we find thatthe actual masses of the simulated secondary stars are 0 . .
15 M (cid:12) for A and B , which are in good agreementwith the estimated values.We note that the regions with min( S T ) < t = 0, both of whichare already shorter than λ MU (cid:39) S T ) ∼ S T ) ∼ t = 2 .
68 yr, but afterwards they strongly inter-act with each other and exchange their masses. The region at (
R, φ ) (cid:39) (1 AU , S T ) < t = 0. However, this region corresponds to the tip of theelongated envelope of the main star. Our analysis assuming ananulus is not applicable to such a bar structure. In the mass estimations of the secondary stars, we define theirclump centres to be the positions of the highest Σ. The masses ofthe secondary stars are computed as the enclosed masses within0 . (cid:39) Figure 2.
Our polar-map analysis for the non-fragmenting case.The top panels show gas distributions in the face-on (left) andedge-on (right) views. The second panel from the top shows theface-on density map but in the polar coordinates. The bottomtwo panels indicate the instability parameters with and withoutthe thickness correction, i.e., min S T and min S . Fig. 2 shows the result for another run in which spiral armsdo not fragment. The central star has a mass of 1 .
10 M (cid:12) within 0 . T (cid:39) τ cool Ω ∼ τ cool Ω ∼
1. It appears that the coolingcriterion of Gammie (2001) is not readily applicable to thefragmentation of the spiral arms. Our SAI analysis givesmin( S T ) > should be stable. The disc has a large thickness as canbe seen in the edge-on view, where we also see significantflaring around the position of the spiral arm. Comparingthe bottom two panels, we find that the thickness correc-tion is critically important. Without the correction, thereare portions with min( S ) < Q < c (cid:13) , 1– ?? AI in protostellar discs no fragmentation is seen in the subsequent evolution. Thisagain supports our argument that Toomre instability is notthe direct cause of secondary protostar formation. We have successfully applied the SAI analysis developedin our Paper I/II to primordial protostellar discs. Our im-proved analysis with finite thickness correction and shearingmotions accurately describe the evolution and stability ofspiral arms. The physical process of secondary protostar for-mation can be summarized as follows. First, fragmentationis triggered by gravitational instability of the arms, ratherthan by non-linear processes following the two-dimensionalcollapse owing to Toomre instability. Second, fragmentationis not driven by gas cooling but by self-gravity of the arm.When the self-gravity exceeds the stabilizing effects due togas pressure and the Coriolis force, the spiral arm is subjectto fragment and bears a secondary protostar.We confirm a two-step process for secondary protostarformation as suggested by Takahashi et al. (2016). Namely,the circumstellar disc around a primordial protostar can ex-cite spiral arms through Toomre instability and/or swingamplification, and then SAI operates in some of the spiralarms when they satisfy min( S T ) < S T ) < N -body simulations where there is no dissipation (Paper I).With this result, we therefore conclude that rapid cooling isneither sufficient nor necessary conditions for fragmentationin protostellar discs.From the instability parameter given by equation (3),it is expected that SAI tends to operate in massive (largeΥ), thin (small W and h ) and cold (small c snd and σ φ ) spiralarms around a less massive central star (small κ ). Therefore,we expect SAI-driven fragmentation occurs preferentially inthe early stages of protostellar evolution when the gas massaccretion rate is large (Yoshida et al. 2008; Hosokawa et al.2016). Further studies are clearly warranted to investigatethe formation of protostars via arm/disc fragmentation dur-ing later evolutionary stages of protostellar systems.Formation of binary or multiple stars has been found ina number of numerical simulations (e.g. Clark et al. 2011;Stacy & Bromm 2014) but the exact reason for the discfragmentation has not been identified. Often, complexities in technical implementation of sink particle generation andstar-gas interaction hampered proper understanding of thephysical mechanism. We have shown, for the first time, thatthe SAI drives the formation of secondary protostars. It isimportant to explore other mechanisms such as break-up orfission of protostars (Tohline 2002) and filament fragmenta-tion (Chiaki et al. 2016). ACKNOWLEDGEMENTS
We are grateful to the reviewer for his/her useful com-ments. We thank Volker Springel for kindly providing thesimulation code
Arepo and Thomas Greif for providingthe outputs of his simulations. This study was supportedby World Premier International Research Center Initia-tive (WPI), MEXT, Japan and by SPPEXA through JSTCREST JPMHCR1414. SI receives the funding from KAK-ENHI Grant-in-Aid for Young Scientists (B), No. 17K17677.Final note: This paper is dedicated to our friend and col-league, Dr. Thomas H. Greif (1981-2019).
REFERENCES
Bromm V., Yoshida N., 2011, ARA&A, 49, 373Chiaki G., Yoshida N., Hirano S., 2016, MNRAS, 463, 2781Clark P. C., Glover S. C. O., Smith R. J., Greif T. H.,Klessen R. S., Bromm V., 2011, Science, 331, 1040Elmegreen B. G., 2011, ApJ, 737, 10Gammie C. F., 2001, ApJ, 553, 174Goldreich P., Lynden-Bell D., 1965, MNRAS, 130, 97Greif T. H., Bromm V., Clark P. C., Glover S. C. O., SmithR. J., Klessen R. S., Yoshida N., Springel V., 2012, MN-RAS, 424, 399Hirano S., Hosokawa T., Yoshida N., Umeda H., OmukaiK., Chiaki G., Yorke H. W., 2014, ApJ, 781Hosokawa T., Hirano S., Kuiper R., Yorke H. W., OmukaiK., Yoshida N., 2016, ApJ, 824, 119Inoue S., Dekel A., Mandelker N., Ceverino D., BournaudF., Primack J., 2016, MNRAS, 456, 2052Inoue S., Yoshida N., 2018, MNRAS, 474, 3466Inoue S., Yoshida N., 2019, MNRAS, 485, 3024Omukai K., 2001, ApJ, 546, 635Rice W. K. M., Armitage P. J., Bate M. R., Bonnell I. A.,2003, MNRAS, 339, 1025Romeo A. B., 1992, MNRAS, 256, 307Safronov V. S., 1960, Annales d’Astrophysique, 23, 979Springel V., 2010, MNRAS, 401, 791Stacy A., Bromm V., 2014, ApJ, 785, 73Stacy A., Bromm V., Lee A. T., 2016, MNRAS, 462, 1307Susa H., 2019, ApJ, 877, 99Takahashi S. Z., Tsukamoto Y., Inutsuka S., 2016, MN-RAS, 458, 3597Tohline J. E., 2002, ARA&A, 40, 349Toomre A., 1964, ApJ, 139, 1217Tsukamoto Y., Takahashi S. Z., Machida M. N., InutsukaS., 2015, MNRAS, 446, 1175Vandervoort P. O., 1970, ApJ, 161, 87Yoshida N., Omukai K., Hernquist L., 2008, Science, 321,669 c (cid:13) , 1–, 1–