Effects of turbulence in the Circumnuclear Disk
DDraft version February 22, 2021
Typeset using L A TEX twocolumn style in AASTeX63
Effects of turbulence in the Circumnuclear Disk
Cuc K. Dinh, Jesus M. Salas, Mark R. Morris, and Smadar Naoz
1, 2 Dept. of Physics & Astronomy, University of California, Los Angeles, CA, 90095, USA Mani L. Bhaumik Institute for Theoretical Physics, Department of Physics and Astronomy, University of California, Los Angeles, CA90095, USA
ABSTRACTA Circumnuclear Disk (CND) of molecular gas occupies the central few parsecs of the Galactic Center.It is likely subject to turbulent disruptions from violent events in its surrounding environment, butthe effect of such perturbations has not yet been investigated in detail. Here we perform 3D, N-body/smoothed particle hydrodynamic (SPH) simulations with an adapted general turbulence drivingmethod to investigate the CND’s structural evolution, in particular its reaction to varied scales ofinjected turbulence. We find that, because of shear flow in the disk, transient arcs of gas (streams)naturally arise when turbulence is driven on large scales (up to ∼ INTRODUCTIONThe Galactic Center (GC) is a vibrantly active regionhosting a 4 × M (cid:12) supermassive black hole mani-fested as the radio source Sgr A* (e.g., Ghez et al. 2008).The black hole is orbited by a molecular gas structureknown as the Circumnuclear Disk (CND, Becklin et al.1982; Morris & Serabyn 1996). This moderately dense(10 − cm − ) and massive (few 10 M (cid:12) ) disk has aninner radius of ∼ . Corresponding author: Cuc K. [email protected] arms (e.g., G¨usten et al. 1987; Sanders 1998; Christo-pher et al. 2005; Montero-Casta˜no et al. 2009; Hsiehet al. 2017). Density clumps are also frequently men-tioned in the literature, though their definition variesdepending on the study (e.g., G¨usten et al. 1987; Sut-ton et al. 1990; Jackson et al. 1993; Marr et al. 1993;Shukla et al. 2004; Christopher et al. 2005; Montero-Casta˜no et al. 2009; Mart´ın et al. 2012; Lau et al. 2013;Smith & Wardle 2014).At present, the origin of the CND remains an openquestion; some studies argue that this disk is a transientstructure, while others suggest that it is quasi-stable andlong-lived. The transient argument is largely based onthe disk’s non-uniformity: the presence of clumps andstreams (e.g., Genzel et al. 1985; Christopher et al. 2005;Montero-Casta˜no et al. 2009; Mart´ın et al. 2012; Lauet al. 2013). Because the rotational dynamics of the diskare dominated almost entirely by the central black hole,the disk should experience strong differential rotationthat smooths out non-uniformities such as clumps on anorbital time scale. The observed non-axisymmetric den-sity distribution has therefore been invoked as evidencethat either the CND was recently formed, and that theobserved structure is therefore short-lived (e.g., G¨ustenet al. 1987; Requena-Torres et al. 2012), or that thestructures within the disk are dense enough to be tidallystable (Jackson et al. 1993; Vollmer & Duschl 2001a,b; a r X i v : . [ a s t r o - ph . GA ] F e b Dinh et al.
Shukla et al. 2004; Christopher et al. 2005; Montero-Casta˜no et al. 2009). However, multiple studies havefound that the clumps have densities ranging up to afew times ∼ cm − , which falls well below the criti-cal Roche density of ∼ × cm − (e.g., G¨usten et al.1987; Requena-Torres et al. 2012; Lau et al. 2013; Millset al. 2013). Thus, the observed density concentrationsare tidally unstable, which has led some to concludethat the CND is a transient structure, with a lifetimeof a few dynamical timescales ( ∼ -10 yrs, G¨ustenet al. 1987). A popular model to account for the tran-sient structures within the CND is that a passing cloudwas gravitationally captured by the black hole. This“infalling cloud” was then tidally disrupted into what isnow the CND (Sanders 1998; Bradford et al. 2005; Bon-nell & Rice 2008; Wardle & Yusef-Zadeh 2008; Hobbs& Nayakshin 2009; Alig et al. 2011; Oka et al. 2011;Mapelli et al. 2012; Mapelli & Trani 2016; Trani et al.2018; Goicoechea et al. 2018; Ballone et al. 2019).The long-term evolution of the CND is therefore thequestion before us. Are we seeing it soon after it hasbeen created, and will it persist as a quasi-stable struc-ture well into the future, or will it be dissipated by var-ious forms of energetic activity at the Galactic center?An important factor to consider when addressing thesequestions is that the CND is subject to various sourcesof turbulence, which come from the active environmentof the Galactic Center. These can be in the form of su-pernovae (e.g., Mezger et al. 1989; Mart´ın et al. 2012),outflows from the Nuclear Stellar Cluster (e.g., Morris& Serabyn 1996; Genzel et al. 2010) and Sgr A*, andturbulence due to the CND’s toroidal magnetic field,which has a strength on the order of a few milligauss(e.g., Marshall et al. 1995; Bradford et al. 2005; Hsiehet al. 2018). These perturbers may play a vital role inthe CND’s morphology. However, the general effect ofturbulence driving on the disk’s structural evolution isstill largely uncertain.In this paper, we perform 3D smoothed particle hy-drodynamic (SPH) simulations that include self-gravityand a turbulence driving method to investigate how dif-ferent scales of turbulence affect the CND’s dynamicalevolution. Our aim is to determine whether turbulencecan help us distinguish between the transient and non-transient ideas. This paper is organized as follows: InSection 2 we describe the numerical code and methods.In Section 3 we describe our simulations results. Finally,we summarize our paper and present some insights be-hind our results in Section 4, as well as compare our re-sults with observations to determine whether there areany similar morphological features. NUMERICAL METHODSWe used the N-body/SPH code Gadget2 (Springel2005), which is based on the tree-Particle Mesh methodfor computing gravitational forces and on the SPHmethod for solving the Euler equations of hydrodynam-ics. The smoothing length of each particle in the gas isfully adaptive down to a set minimum of 0.001 pc. Gad-get2 employs an entropy formulation of SPH, as out-lined in Springel & Hernquist (2002), with the smooth-ing lengths defined to ensure a fixed mass (i.e., fixednumber of particles) within the smoothing kernel volume(set at N neigh = 64). The code adopts the Monaghan-Balsara form of artificial viscosity (Monaghan & Gingold1983; Balsara 1995), which is regulated by the parame-ter α MB , set to 0 . Gravitational potential
We include into Gadget2 the gravitational potentialof the inner 10 pc of the Galaxy, which includes thesupermassive black hole (Sgr A*) and the Nuclear StellarCluster (NSC, e.g., Do et al. 2009; Sch¨odel et al. 2009).In other words, the total potential, Φ T is:Φ T = Φ bh + Φ NSC . (1)The gravitational potential of the black hole is includedas a point-mass potential:Φ bh = − GM bh r , (2)where M bh = 4 × M (cid:12) , and r is the radial distancefrom SgrA*. Furthermore, to include the gravitationaleffects of the NSC, we adopt the potential described byStolte et al. (2008):Φ NSC = 12 v c ln ( R c + r ) , (3)where v c = 98 . R c = 2 pc is the core radius.2.2. External pressure
The interstellar medium surrounding the CND is mod-eled via an external pressure term to approximate a con-stant pressure boundary. Following Clark et al. (2011),we modify Gadget2’s momentum equation (Springel & urbulence in the CND dv i dt = − (cid:88) j m j (cid:34) f i P i ρ i ∇ i W ij ( h i ) + f j P j ρ j ∇ i W ij ( h j ) (cid:35) , (4)where v i is the velocity of particle i , m j is the mass ofparticle j , P i is the pressure, ρ i is the density, W ij ( h i )is the kernel function which depends on the smoothinglength h i , and f is a unitless coefficient that dependson ρ i and h i . We replace P i and P j with P i − P ext and P j − P ext , respectively, where P ext is the external pres-sure. The pair-wise nature of the force summation overthe SPH neighbors ensures that P ext cancels for particlesthat are surrounded by other particles. At the bound-ary, where the P ext term does not disappear, it mimicsthe pressure contribution from a surrounding medium(Clark et al. 2011). We set P ext equal to 10 − ergscm − , an approximate value for the GC (Spergel & Blitz1992; Morris & Serabyn 1996).2.3. Initial conditions
For simplicity, we start with a gas disk having initialparameters drawn from observations: an inner radius of1 . . . SPH particles and atotal mass of 4 . × M (cid:12) . The particles are initially incircular orbits, with their velocities calculated using thepotential described in Section 2.1. All simulations wererun using an isothermal equation of state with T = 200K. 2.4. Turbulence driving
We use the turbulence driving method described bySalas et al. (2019), which consists of a Fourier forcingmodule modelled with a spatially static pattern. Thisturbulence module follows the methods by Stone et al.(1998) and Mac Low (1999). For completeness, we re-capitulate the key factors of the algorithm here.First, we create a library of 10 files of turbulence,which our modified version of the Gadget2 code readsin at the start of the simulation. Each file contains aunique realization of a turbulent velocity field (in theform of a 3D matrix) with power spectrum P ( k ) ∝ k − (where k is the wavenumber). This power law is steeperthan Kolmogorov turbulence, but it is suitable for com-pressible gas (Clark et al. 2011). Each 3D matrix iscreated using fast Fourier transforms inside a 128 box,following the techniques described in Rogallo (1981) andDubinski et al. (1995). These 3D turbulence matricescan be visualized as lattice cubes (or grids) containing128 × ×
128 equally spaced lattice points. We set thephysical size of these cubes to L cube per side, dependingon the turbulence model (see Section 2.5). In some studies, the driving module only containspower on the larger scales (e.g., Federrath et al. 2010).This type of driving models the kinetic energy in-put from large-scale turbulent fluctuations, which thenbreak up into smaller structures as the kinetic energycascades down to scales smaller than the turbulence in-jection scale. However, in SPH, the artificial viscos-ity can damp this energy cascade and prevent it fromreaching the smaller scales. Consequently, to createthe different realizations of turbulent velocity fields, weuse a discrete range of k values from k min = 2 to k max = 128, thus effectively injecting energy on scalesbetween λ max = L cube / λ min = L cube / ∼ N of these cubiclattices per side (where N = 4 or 100, depending on themodel, see Section 2.5), each drawn randomly from ourlibrary. Using this method, the spatial resolution of theturbulence in our simulation domain is the same as theresolution of an individual turbulence lattice cube.To drive the turbulence, we follow a method similar tothat described by Mac Low (1999): every N t timesteps(the timestep in all our simulations is fixed to ∆ t s = 100yrs) we add a velocity increment to every SPH particle, i , given by:∆ (cid:126)v i ( x, y, z ) = F ( ρ i ) A(cid:126)I i ( x, y, z ) , (5)where (cid:126)I i ( x, y, z ) is the turbulent velocity interpolatedfrom the turbulence field of the lattice cube that con-tains the particle in question. The amplitude A is cho-sen to maintain a constant kinetic energy input rate˙ E in = ∆ E in / ( N t ∆ t s ), and the term F ( ρ i ) is a factorthat depends on the particle’s density, ρ i , and on theturbulence model (see Section 2.5). Any particle out-side the simulation domain does not receive any turbu-lent energy.Our turbulence implementation contains two free pa-rameters: ∆ E in , the total energy input per instance ofturbulent energy injection, and ∆ t = N t ∆ t s , the timebetween velocity “kicks”. We show below the chosenparameters in our models.2.5. Turbulence models
Here we describe our turbulence models and assump-tions. We consider 1) a large-scale model and 2) a small-scale model, each corresponding to different physical sce-narios. Table 1 summarizes the runs and parameters tobe used based on our models.
MODEL - SMALL TURBULENCE SCALES (STS) We first investigate a small-turbulence-scale model bysetting L cube = 0 . Dinh et al.
Run name Turbulence scale ∆ E in (ergs) ∆ t (yrs)STS-t100-47 Small 10 − pc)STS-t500-47 Small 10 − pc)LTS-t5-49 Large 10 (2-0.3 pc)LTS-t5-50 Large 10 (2-0.3 pc) Table 1.
Parameters of all turbulence runs. ∆ E in corre-sponds to the total energy per injection, and ∆ t = N t ∆ t s isthe time between injections. λ = 0 .
05 pc to 10 − pc (for k = 2 and k = 128, respec-tively). We fill the volume of our simulation domain(10 × ×
10 pc) with N = 100 of these boxes perside. This model can represent, for example, turbulencedriven by stellar winds, occasional novae, and outflowsfrom Wolf-Rayet stars throughout the CND. We expectsuch perturbations to affect the surrounding gas of theCND on few-hundred-year timescales. Therefore, we runtwo injection timescales, ∆ t = 100 yrs and ∆ t = 500yrs. The injection energy was set to ∆ E in = 10 ergs.In comparison, the wind of a standard AGB star emitsabout 10 ergs of energy every 100 yrs, assuming a massloss rate of 10 − M (cid:12) /yr and a wind speed of 15 km/s.The energy injection in this STS model would then beequivalent to the wind of ∼ AGB stars distributedthroughout the disk, assuming all of that energy is trans-ferred to the gas .Furthermore, in Salas et al. (2019) the turbulent veloc-ity kicks have a density dependence factor F ( ρ ) = √ Gρ ,which helps counteract the effects of self-gravity. Thisdensity factor is not needed in this STS model sincea temperature of 200 K creates a sufficient thermalpressure gradient in the gas disk to counteract its self-gravity. Thus, in this model we eliminate this densitydependence by setting F ( ρ ) = 1. MODEL 2 - LARGE TURBULENCE SCALES (LTS)
We next consider a large-turbulence-scale model bysetting L cube = 4 pc. Thus, energy is injected on scales λ = 2 − .
03 pc (for k = 2 and k = 128, respectively).This model represents, for example, turbulence driven The energy injection in the STS model represents an upper limit.However, as we discuss in Section 3, even this large energy injec-tion does not affect the CND significantly. Thus, any lower (andmore realistic) energy injection would have even less discernibleeffects on the CND. by supernovae (SNe). We fill the volume of our sim-ulation domain with N = 4 of these boxes, each boxoccupying a quadrant of the simulation domain (whichin this case has a volume of 8 × × yrs,roughly the average lifetime of a supernova remnant.Furthermore, given the number of massive stars withinthe central few parsecs of the Galaxy (e.g., Do et al.2013; Nishiyama & Sch¨odel 2013) we might expect theirrecurrence timescale to be also on the order of 10 years.SN explosions release ∼ ergs of energy, which drivesa blast wave through the ambient ISM. However, theproperties of the environment within which a SN ex-plodes strongly affect the fraction of energy that is de-posited into the surrounding gas (Dwarkadas & Gruszko2012). For example, Cox (1972) estimates that 10% ofthe initial SN energy is retained as thermal and kineticenergy in the gas. Chevalier (1974) obtained a frac-tion between 4-8%, depending on the ambient density,assumed homogeneous. Other studies have found frac-tions from up to 50% (McKee & Ostriker 1977; Cowieet al. 1981) to just a few percent (Slavin & Cox 1992;Walch & Naab 2015). Therefore, we opted for testingenergy injection values of ∆ E in = 10 and 10 ergs,which correspond to 1% and 10% of a typical SN en-ergy, respectively. Additionally, here we use a densityfactor of F ( ρ ) = ρ − / (see Equation 5) to capture thedensity dependence of the expansion velocity of a SNshockwave during the snowplow phase (e.g., Cioffi et al.1988).Finally, as was shown in Salas et al. (2019), our turbu-lence method enhances the inward angular momentumtransport of gas. In this model, rather than a slow in-ward migration, the large-scale perturbations could po-tentially send gas onto ballistic trajectories towards thecenter. Thus, to counteract the rapid filling of the CNDcavity due to turbulence, we mimic the outflow fromthe young nuclear cluster of massive, windy stars. Weachieve this by adding a radially outwards velocity, v i ( r ),to every SPH particle, i , within 1 . v i ( r ) = f v esc,i ( r ) (cid:114) n n i , (6)where v esc,i ( r ) is the escape velocity of an SPH particleat radius r , n = 10 cm − is the assumed outflow num- urbulence in the CND Figure 1.
Column number density of the Small Turbulence Scale (STS) model featuring its evolution through the simulation.The simulation was run for 3 Myrs with injection energy of 10 ergs at intervals of ∆ t N = 100 years. Top panel: Face-on view;bottom panel: edge-on view. ber density where it strikes the inner edge of the CND(the simulations by Blank et al. 2016 assumed an outflowdensity of n = 100 cm − at r = 0 . /r , n ≈
10 cm − at r = 1 . n i is the number density of an SPHparticle at the inner edge of the CND. Finally, we addthe free parameter f in order to adjust the magnitudeof this radial velocity so that the inner edge maintainsa stable radius. This also avoids launching particles offto infinity, which creates numerical errors in the simula-tions. In our calculations, f = 0 . RESULTSWe tested each turbulence run for 3 Myrs, which isfar greater than the dynamical timescale of the disk ( ∼ . − . Model 1 (STS)
Snapshots of the time evolution of the STS-t100-47run are shown in Figure 1, where we show column num-ber density maps at 1, 2 and 3 Myrs (we present only oneof the runs, since both STS simulations are qualitatively identical). Given that the spatial scales over which theenergy injection is distributed are relatively small, thenet effect of turbulence is manifested on only the small-est of scales, as expected. The small-scale convergingflows create a granular structure throughout the disk,and this texture continues during the entire simulation.Furthermore, our driving method creates a “turbu-lent” viscosity which promotes the transfer of angularmomentum, as discussed in Salas et al. (2019). Theseviscous forces broaden the disk slightly: the disk ex-pands radially outwards, which increases the outer ra-dius, while the inner radius decreases slightly. The cav-ity closes at a faster rate than the disk’s expansion,but the inward migration rate was still relatively slow.We did not impose an artificial radial outflow as in theLTS runs, since we presumed, correctly, that it was notneeded. Overall, the STS models reached steady stateearly on in the simulations.3.2.
Model 2 (LTS)
The evolutionary progressions of the LTS-t5-49 andLTS-t5-50 runs are shown in Figures 2 and 3, respec-tively. The effects of the turbulence are much moreprominent here than in the STS models. The large-scale turbulence perturbations are quickly stretched outinto long arcs by the strong tidal shear, which promotesthe formation of continuous, orbiting spiral stream seg-ments. We note that while these features are appar-ent in both runs, streams in the LTS-t5-49 run are lessdistinguishable in comparison to those in the LTS-t5-
Dinh et al.
Figure 2.
Column number density maps featuring the evolution of the Large Scale Turbulence model (LTS-t5-49). The runwas run for 3 Myrs with injection energy of 10 ergs at intervals of ∆ t = 10 years. Top: face-on view, rotates counter-clockwise. The disk is fully and densely permeated by identifiable spiral streams by ∼ .
50 run. For example, the streams fully encapsulate thedisk by ∼ . − ∼ . .
71 Myrs (left panelin Figure 4). This imprint occurs at every injection,and it was was an expected result, since the size of theturbulence grid is comparable to the scale of the CND.While the injection procedure in this LTS model is ge-ometrically oversimplified, the grid’s imprint, as well as the turbulent perturbations, disappear and are quicklysheared into orbiting streams, leaving no memory of theshape of the initial disturbance domain.The bottom panels of Figures 2 and 3 show thatmaterial was perturbed off the disk as the simulationsevolved. This is because the turbulent velocity kicks inthe LTS simulations are density-dependent: due to theimposed density factor of F ( ρ ) = ρ − / (see Section 2.4),lower density promotes higher turbulent kicks. Thus,the turbulence is very effective in perturbing the lowerdensity gas off the disk plane, creating an “atmosphere”of low-density gas surrounding the higher density disk.As such, the LTS-t5-50 run, due to its higher energy in-jection rate, had more material perturbed off the diskplane (see bottom panel of Figure 3) compared to the urbulence in the CND Figure 3.
Column number density maps depicting the evolution of the Large Turbulence Scale model (LTS-t5-50). Thesimulation was run for 3 Myrs with injection energy of 10 ergs at intervals of ∆ t = 10 years. Top: face-on view, rotatescounter-clockwise. Much more distinct streams with higher densities (especially comparing to LTS-t5-49 at 3 Myrs). Bottom:edge-on view. More prominent perturbations off the disk plane than in the LTS-t5-49 model, as expected. LTS-t5-49 run (bottom panel of Figure 2). In the sec-tions that follow, we focus our analysis on the LTS-t5-50model due to its more prominent structure.We note the presence of compact, high-density con-centrations (the small white points in Figures 2 and 3)in both runs. These points formed early on and slowlyincreased in number over the course of the simulations.These clumps originate because here, local regions col-lapse under self-gravity. Particles pile on top of eachother, and due to the nature of the kernel used by Gad-get2, once the distance between particles approachesthe smoothing length, the pressure gradient is no longercorrect and the particles stick together, creating high-density concentrations. These clumps are, in effect, sinkparticles that might have registered as stars if the sim-ulations had included a prescription for star formation. Many of these clumps originate in the inner rim of thedisk, which is a consequence of the piling up of parti-cles due to the artificial radial outflow. Thus, we ignorethese clumps and consider them as numerical artifacts.Their total mass is negligible, so they don’t affect theevolution of the disk.3.3.
Comparing with observations: clumps
In order to better compare our simulation with obser-vations of the CND, we must account for 1) the orienta-tion of the CND in the sky, and 2) the display functionof the observational data, be it linear or logarithmic.Thus, we orient the model disk to account for the ob-served orientation, with a disk inclination of 70 ◦ andwith the major axis at a position angle of 20 ◦ (e.g.,Mart´ın et al. 2012; Lau et al. 2013). We also display this Dinh et al.
Figure 4.
Column Number density map depicting the evolution of the Large Turbulent Scale model (LTS-t5-50) between twoenergy injections. Face-on view with the disk rotating counterclockwise. The imprint of the turbulence grid just after it isimposed is apparent at 1 .
71 Myrs (left panel). While the geometry of the turbulence injection is physically oversimplified, thisstructure is quickly dispersed by differential rotation, forming orbiting spiral streams that would undoubtedly be also producedby a more spatially continuous injection prescription. orientation-adjusted view with a linearly scaled colormap. An example is shown in Figure 5 (bottom leftpanel), which shows a snapshot of the LTS-t5-50 runat 1 .
83 Myrs. In addition, we convolve this columndensity map of the tilted LTS-t5-50 model with a 2DGaussian kernel to match the typical full width at half-maximum beam size of existing observations (3 . .
15 pc, e.g., Mart´ın et al. 2012; Lau et al. 2013; Hsiehet al. 2017). As shown in the lower right panel of Figure5, convolving the model gives rise to apparent clumpswhere the distorted streams become fortuitously pro-jected and smoothed. However, most of these apparentclumps are not identifiable as such in the face-on viewof the simulation (top row of Figure 5).The sizes of these apparent clumps ( ∼ . − . ∼ . .
54 pc), densities( ∼ to 10 cm − ), and radial distances from Sgr A*.In our model, however, apparent clumps emerge asa result of the disk projection, coupled with samplingwith finite spatial resolution. In addition, when consid-ering the two presentations of the orientation-adjusted model (lower left and right panels of Figure 5), we findthe clumps to coincide with the densest sections of thestreams. This result leads us to pose the question ofwhether at least some of the clumps seen in observa-tional studies might be, in reality, turbulence-induced,transient unresolved structures in the CND. Perhaps fur-ther investigation of the CND’s internal structure withhigher resolution instruments such as ALMA could set-tle this question.3.4. Comparing with observations: streams
The intrinsic structure in our LTS models consists oftidally stretched streams, which are evident when in-specting the column number density maps (see Figures3 and 5), especially the linearly scaled maps (e.g., topright and bottom left panels of Figure 5). However,we can use more sophisticated methods to investigatethese features in order to reduce biases when determin-ing which particles correspond to which stream.We employed the data clustering algorithm DBSCAN(part of the Python library scikit-learn , Pedregosa et al.2011) to categorize the streams in our LTS model. DB-SCAN views clusters as areas of high density (i.e., par-ticle concentrations) separated by areas of low density.As such, DBSCAN is able to find clusters of any shape,as opposed to the familiar “k-means” algorithm, whichworks better when clusters are approximately spheri-cal. There are two parameters to the DBSCAN al-gorithm: min samples , which represents the minimumnumber of points required to form a dense region, and eps ( (cid:15) ), which determines the minimum number of par-ticles in a neighborhood for that neighborhood to beconsidered a “point”. The algorithm starts with an ar-bitrary starting point. This point’s (cid:15) -neighborhood is urbulence in the CND Figure 5.
Column number density maps of the LTS-t5-50 model (domain size of 11x11 pc) at 1 .
83 Myrs. The upper two panelsshow the face-on log-scale column number density map (upper left) and its linear scale counterpart (upper right). The lowertwo panels compare the linear scale density map of the orientation adjusted model (lower left) with its convolved counterpart(lower right). The identified clumps are circled in yellow in both maps. There also appears to be a clump at the inner edge ofthe disk, however we do not consider it since it is likely a high density region due to the effects of the artificial radial outflow.The clumps in the convolved, orientation-adjusted view (bottom right) seem to be a confluence of distorted streams. Therefore,the clumps are mostly artifacts of both resolution and orientation. explored, and if it contains sufficiently many points, acluster is started. Otherwise, the point is labeled asnoise. Note that this point might later be found in asufficiently sized (cid:15) -environment of a different point andhence be made part of a cluster.We built a routine using DBSCAN and applied it tothe set of particle positions that have a density of atleast 10 cm − , in order to capture the densest parts ofthe streams. We set eps , which determines the neigh-borhood radius of a “point” in a cluster, to 0 .
05 pc, and min samples to 10. All particles within the inner 2 par- secs of the simulation were removed for this analysis,because the high-density region at the disk’s inner edgeformed due to the interaction of the gas with the artifi-cial radial outflow. Therefore this high density region isignored for our purposes.We applied our DBSCAN routine to the LTS-t5-50simulation at different times and were able to identifystreams in all cases. These streams are transient, withan approximate timescale dictated by the turbulence in-jection interval. The shape and quantity of the streamsundergo frequent changes throughout the simulation;0
Dinh et al.
Figure 6.
Results from employing DBSCAN on the LTS-t5-50 model (∆ E in = 10 ergs with ∆ t = 10 years) at 2.39 Myrs.The upper left panel displays seven identified streams (colored) imposed on the entire disk model (face-on view). The upperright panel accounts for the observed orientation. We do not consider the particles within a radius of 2 pc in the stream findinganalysis. We also note here that that the streams were first identified and colored in the face-on view before the model wasrotated, thus they are less evident in the rotated orientation. The bottom panel shows the plot of radial velocities vs positionangle of the particles in their corresponding streams (using the same colors). We limit the densities of the particles to be between n = 10 and n = 10 cm − . This trend closely follows the data points from Mart´ın et al. (2012), which are reproduced in theFigure: square and circular points correspond to measurements made in the northern and southern CND, respectively. urbulence in the CND − ∼ ◦ ,there is noticeable overlap between the purple and or-ange streams. It is possible that similar structures seenin CND observations may present a similar degeneracy,especially if observed structures have been categorizedby visual inspections. Thus, we note the importance of aclustering algorithm to find structures in the CND moreeffectively.However, the DBSCAN routine was less effectiveat identifying streams when it was applied to theorientation-adjusted view because the tilt of the diskdecreases particle separations and superposes streamsupon each other, as shown in the upper right panel ofFigure 6. In addition, streams that happened to con-nect at one or more locations are often unidentifiable asseparate entities.Consequently, our results indicate that there could bemore streams than might normally be observed due tothe disk’s orientation or the resolution. Using an elabo-rated version of DBSCAN or some equivalent algorithm,perhaps one that includes velocity information, could fa-cilitate the identification of streams and other structuresin both observations and models. DISCUSSIONIn this paper, we have presented 3D SPH simulationsthat explore the effects of turbulence on the dynamicalevolution of the CND. These simulations, which includeself-gravity and an adapted turbulence driving mech-anism, are based on two models representing differentdriving scales. The Small Turbulence Scale (STS) models did notshow sufficient structure to match the morphology re-vealed by observations (see Figure 1). This is not sur-prising since the scales at which turbulence is stirred inthis model are quite small (0 .
05 – 10 − pc) and theireffect on the much larger scale CND is only to endow itwith a grainy appearance. This implies that sources ofturbulence that operate on such small scales only havea minor impact on the overall dynamics of the CND.Sources that operate on these scales include not onlystellar winds, as previously mentioned, but also novae,pulsar wind nebulae, winds from AGB stars, winds fromWolf-Rayet and massive main-sequence stars, and mag-netorotational instabilities (MRIs) from the predomi-nantly torodial, few-milligauss magnetic field threadingthrough the CND (e.g., Werner et al. 1988; Hildebrandet al. 1990, 1993, C.D. Dowell et al. in preparation).In contrast, we found that turbulence driving in theLarge Turbulence Scales (STS) models, which includedthe effects of larger scale (2-0 .
03 pc) sources of turbu-lence, such as supernova blasts, had much greater effectson the disk’s evolution (see Figure 3). Here, the large-scale turbulence perturbations led to the formation ofstreams and clumps, some of which might resemble ob-servational features of the CND. However, these streamsare partially obscured when we account for both thedisk’s orientation and the general resolution of observa-tions. As such, our results suggest that not only canwe match the streams observed in the CND, but therecould be more streams than observed.In addition, convolving the orientation-adjustedmodel (see Figure 5) to match the resolution of obser-vations revealed clump-like structures similar to thosedetailed in a number of studies. Such “clumps” are diffi-cult to discern in the original unconvolved, orientation-adjusted model. They are also not evident in eitherthe linear or log-scale column-density maps of the face-on model. This suggests that at least some observedclumps could be artifacts of both the resolution and theorientation of the CND. We posit that the CND is sub-ject to episodic disturbances on large scales that, whensheared, leave a persistent pattern of spiral stream seg-ments. In the absence of continued episodic injectionsof turbulence, the disk returns to a smooth state on atime scale of a few million years.4.1.
Origin and Longevity of the CND
The results of our LTS model show that streams andapparent clumps both naturally arise in a long-lived disksubjected to large-scale perturbations. Thus, tidallyunstable clumps do not necessarily indicate that thedisk itself is a transient feature, as has often been ar-2
Dinh et al. gued (e.g., G¨usten et al. 1987; Requena-Torres et al.2012). Our LTS simulations show that these clumpsand streams are dynamic, constantly appearing and dis-sipating due to the injected turbulence. Simulations byBlank et al. (2016) also show that interactions betweenthe nuclear stellar cluster’s strong winds and the orbit-ing, annular disk cause instabilities at the inner edge ofthe CND. The existence of such mechanisms for continu-ously creating tidally unstable and therefore short-livedclumps and streams implies that such features cannot beused to conclude either that the CND is a transient fea-ture or that the clumps within it must be dense enough( n H > × cm − ) to be virialized and thereforetidally stable, as has often been argued (e.g., Vollmer &Duschl 2001a,b; Shukla et al. 2004; Christopher et al.2005; Montero-Casta˜no et al. 2009). Such high densitieswould imply extremely large CND masses of 10 − M (cid:12) , which are inconsistent with the optically thin far-IR and submillimeter fluxes from the CND (Genzel et al.2010; Etxaluze et al. 2011). The densities required fortidal stability are also much higher than have been in-ferred in many recent studies (Oka et al. 2011; Requena-Torres et al. 2012; Lau et al. 2013; Mills et al. 2013;Smith & Wardle 2014; Harada et al. 2015; Goicoecheaet al. 2018; Tsuboi et al. 2018). We conclude that clumpsand streams in the CND must truly be transient fea-tures.Based on the observation that the CND appears tobe disequilibrated because of the presence of internaltransient structures, many publications model the for-mation of the CND as a recent event that resulted fromthe infall of a dense cloud toward the central black hole(Sanders 1998; Bradford et al. 2005; Bonnell & Rice2008; Wardle & Yusef-Zadeh 2008; Hobbs & Nayak-shin 2009; Alig et al. 2011; Oka et al. 2011; Mapelliet al. 2012; Mapelli & Trani 2016; Trani et al. 2018;Goicoechea et al. 2018; Ballone et al. 2019). The idearests on the assumption that a relatively massive cloudwith a very small total angular momentum can somehowbe produced near the black hole, or can be producedfurther away, but be unimpeded in its accelerating tra-jectory toward the black hole. However, those assump-tions are questionable. Although there remains someuncertainty about the line-of-sight placement of cloudsin the central molecular zone (CMZ), gaseous structuresin the CMZ appear to move on orbits that conform tothe sense of Galactic rotation (e.g., Henshaw et al. 2016),even if they are influenced by a non-axisymmetric poten-tial or by radial accelerations caused by extreme accre-tion or starburst activity at the center. Therefore CMZclouds have considerable angular momentum. Creatinga cloud with near-zero angular momentum, or one with a retrograde orbit that could eliminate the angular mo-mentum of other clouds by colliding with them, wouldtherefore require it to be scattered at a large angle, butthere are no known massive perturbers that could dothat, including other clouds, because the sizes of cloudsare large enough that they cannot approach each otherclosely enough to cause scattering by more than about10 ◦ without colliding and merging, thereby preservingtheir angular momentum. Even with some undefinedscattering process at work, the phase-space volume intowhich clouds would need to be scattered to be broughtfrom a distance of more than several parsecs to the vicin-ity of the black hole is very small; the centrifugal barrierfor a CMZ cloud scattered at almost any angle would beencountered at a radius larger than that of the CND .Furthermore, any large cloud having a non-conformingvelocity would collide within a few dynamical times withother clouds and would be forced eventually into con-formity with the CMZ. Consequently, the production ofthe CND by a radially infalling cloud appears to be ex-tremely difficult. We therefore conclude that the CNDis not a transient feature, but is rather an enduringstructure that is fed piecemeal and relatively slowly andquasi-continuously by tidal streams from CMZ cloudsalready orbiting in its vicinity, as some authors have sug-gested (Ho et al. 1991; Wright et al. 2001; Hsieh et al.2017; Tsuboi et al. 2018).The evolution of the CND is also linked to its propen-sity to occasionally form stars and to its relationshipto the central black hole. The young nuclear cluster ofmassive stars occupying the central 0.5 pc presumablyformed in a starburst event from an earlier manifesta-tion of the CND, when the inner radius extended all theway in to the black hole and probably fed the black holeat a far higher rate than it is fed at present. This couldhappen again if the CND continues to be sustained byinwardly migrating material from the outside and if itundergoes viscous evolution that eventually causes itsinner edge to move inwards. Indeed, this could be arepetitive process in which the CND undergoes a limitcycle of activity punctuated in each cycle by a relativelybrief starburst event (Morris et al. 1999). To addresssuch a scenario, future simulations are needed to accountfor the effects of the magnetic field, CND feeding fromthe outside, and other physical processes not included inour models. Additionally, continued investigations withALMA and other high-resolution instruments will be in- Supernova blast waves can push ambient gaseous material frommoderately large distances into the zone of the CND (Palouˇset al. 2020), but the mass of material moved in that way is muchsmaller than the mass of the CND. urbulence in the CND
Software:
Figure 1 to 5 were done using theSPH visualization software
SPLASH (Price 2007).We used a modified version of the Gadget2 code(Springel 2005) which includes our turbulencemethod. The version of the code can be found athttps://github.com/jesusms007/CNDturbulence. Fi-nally, we used the unsupervised machine learning algo-rithm DBSCAN, which is part of the Python library scikit-learn (Pedregosa et al. 2011).APPENDIX A. ARTIFICIAL OUTFLOWWe showed in Salas et al. (2019) that our turbulence driving module enhances inward angular momentum transportof gas. To counteract the rapid filling of the CND cavity due to turbulence, we need to mimic the effects of the outflowfrom the NSC on the inner edge of the CND. Blank et al. (2016) modelled the NSC’s outflow by assuming a numberdensity of n = 100 cm − and a speed of v = 700 km/s, propagating radially outward from a starting radius of r = 0 . P outflow = ρ v from one side (where ρ and v are the massdensity and speed of the outflow, respectively), and a ram pressure from the opposite direction, P CND = ρ CND v CND (where ρ CND is the mass density of the neighboring gas, which is assumed to be identical to the parcel’s, and v CND is the parcel’s speed as seen from the CND’s frame of reference). The parcel’s speed will then be, v CND = (cid:114) ρ ρ CND v = (cid:114) µ n µ CND n CND v , (A1)where µ and µ CND are the mean mass per particle of the outflow’s ionized gas, and the molecular CND gas, respec-tively. The ratio (cid:112) µ /µ CND is subsumed in the factor f , described below.We then convert this treatment to work with Gadget2. However, here we consider the outflow’s density at r = 1 . /r , we thus adopt a value of n = 10 cm − . Furthermore, becausewe are not accounting for the additional hydrodynamical interaction between the outflow and the CND, a speed of v = 700 km/s, which is much larger than the escape velocity, would send SPH particles flying out to infinity, causingnumerical problems in the code. We simplify this problem by assuming that the largest speed this artificial outflowwill add to an SPH particle (which happens when the SPH particle number density, n i , is equal to the outflow’s, i.e., n i = n . For particles with n i < n , no radial speed is added. This limit is acceptable since we find in our simulationsthat they are very few particles with n i < n throughout the disk) equals a fraction of the particle’s escape speed, f v esc,i . This parameter f was adjusted in order to maintain a stable radius of the CND’s edge. In our testing, wefound f = 0 .
15 satisfied this condition.Thus, to mimic the NSC’s outflow, every timestep we add a radially outward speed to every SPH particle, i , withina radius r = 1 . v i ( r ) = f (cid:114) n n i v esc,i ( r ) , (A2)where v i ( r ), v esc,i ( r ) and n i are values corresponding to a given SPH particle.4 Dinh et al.
Figure 7.
Results from employing DBSCANS on the LTS-t5-50 model (∆ E in = 10 ergs with ∆ t = 10 years) before the rundisplayed in the results (1.75 Myrs). The upper left figure displays resulting four streams (colored) imposed on the entire diskmodel in the face on view. The upper right figure accounts for the rotation of the CND accounted for. The bottom figure mapsthe trend between the position angle and radial velocities of the particles in their corresponding streams (colored). This trendclosely follows that found in Mart´ın et al. 2012, whose values have been fitted above. The square points correspond to whatwas realized as the Northern CND, whereas the circular points belong to the Southern CND.B. FINDING STREAMS WITH DBSCANWe applied DBSCAN on the LTS-t5-50 model at multiple times after 1 . .
75 Myrs and 2 .
99 Myr. urbulence in the CND Figure 8.
Results from employing DBSCANS on the LTS-t5-50 model (∆ E in = 10 ergs with ∆ t = 10 years) after the rundisplayed in the results section (2.99 Myrs). The upper left figure displays resulting four streams (colored) imposed on theentire disk model in the face on view. The upper right figure is the same, but with the rotation of the CND accounted for. Thebottom figure maps the trend between the position angle and radial velocities of the particles in their corresponding streams(colored). This trend closely follows that found in Mart´ın et al. 2012, whose values have been fitted above. The square pointscorrespond to what was realized as the Northern CND, whereas the circular points belong to the Southern CND. Dinh et al.
REFERENCES
Alig, C., Burkert, A., Johansson, P. H., & Schartmann, M.2011, MNRAS, 412, 469,doi: 10.1111/j.1365-2966.2010.17915.xBallone, A., Mapelli, M., & Trani, A. A. 2019, MNRAS,488, 5802Balsara, D. S. 1995, Journal of Computational Physics, 121,357, doi: 10.1016/S0021-9991(95)90221-XBecklin, E. E., Gatley, I., & Werner, M. W. 1982, ApJ, 258,135, doi: 10.1086/160060Blank, M., Morris, M. R., Frank, A., Carroll-Nellenback,J. J., & Duschl, W. J. 2016, MNRAS, 459, 1721,doi: 10.1093/mnras/stw771Bonnell, I. A., & Rice, W. K. M. 2008, Science, 321, 1060,doi: 10.1126/science.1160653Bradford, C. M., Stacey, G. J., Nikola, T., et al. 2005, ApJ,623, 866, doi: 10.1086/428659Chevalier, R. A. 1974, ApJ, 188, 501, doi: 10.1086/152740Christopher, M. H., Scoville, N. Z., Stolovy, S. R., & Yun,M. S. 2005, ApJ, 622, 346, doi: 10.1086/427911Cioffi, D. F., McKee, C. F., & Bertschinger, E. 1988, ApJ,334, 252, doi: 10.1086/166834Clark, P. C., Glover, S. C. O., Klessen, R. S., & Bromm, V.2011, ApJ, 727, 110, doi: 10.1088/0004-637X/727/2/110Cowie, L. L., McKee, C. F., & Ostriker, J. P. 1981, ApJ,247, 908, doi: 10.1086/159100Cox, D. P. 1972, ApJ, 178, 159, doi: 10.1086/151775Do, T., Ghez, A. M., Morris, M. R., et al. 2009, ApJ, 703,1323, doi: 10.1088/0004-637X/703/2/1323Do, T., Lu, J. R., Ghez, A. M., et al. 2013, ApJ, 764, 154Dubinski, J., Narayan, R., & Phillips, T. G. 1995, ApJ,448, 226, doi: 10.1086/175954Dwarkadas, V. V., & Gruszko, J. 2012, MNRAS, 419, 1515,doi: 10.1111/j.1365-2966.2011.19808.xEtxaluze, M., Smith, H. A., Tolls, V., Stark, A. A., &Gonz´alez-Alfonso, E. 2011, AJ, 142, 134Federrath, C., Roman-Duval, J., Klessen, R. S., Schmidt,W., & Mac Low, M.-M. 2010, A&A, 512, A81,doi: 10.1051/0004-6361/200912437Genzel, R. 1989, in IAU Symposium, Vol. 136, The Centerof the Galaxy, ed. M. Morris, 393Genzel, R., Eisenhauer, F., & Gillessen, S. 2010, Reviews ofModern Physics, 82, 3121,doi: 10.1103/RevModPhys.82.3121Genzel, R., Watson, D. M., Crawford, M. K., & Townes,C. H. 1985, ApJ, 297, 766, doi: 10.1086/163574Ghez, A. M., Salim, S., Weinberg, N. N., et al. 2008, ApJ,689, 1044, doi: 10.1086/592738Goicoechea, J. R., Pety, J., Chapillon, E., et al. 2018,A&A, 618, A35 G¨usten, R., Genzel, R., Wright, M. C. H., et al. 1987, ApJ,318, 124, doi: 10.1086/165355Harada, N., Riquelme, D., Viti, S., et al. 2015, A&A, 584,A102Henshaw, J. D., Longmore, S. N., Kruijssen, J. M. D., et al.2016, MNRAS, 457, 2675Hildebrand, R. H., Davidson, J. A., Dotson, J., et al. 1993,ApJ, 417, 565, doi: 10.1086/173336Hildebrand, R. H., Gonatas, D. P., Platt, S. R., et al. 1990,ApJ, 362, 114Ho, P. T. P., Ho, L. C., Szczepanski, J. C., Jackson, J. M.,& Armstrong, J. T. 1991, Nature, 350, 309Hobbs, A., & Nayakshin, S. 2009, MNRAS, 394, 191Hsieh, P.-Y., Koch, P. M., Ho, P. T. P., et al. 2017, ApJ,847, 3, doi: 10.3847/1538-4357/aa8329Hsieh, P.-Y., Koch, P. M., Kim, W.-T., et al. 2018, ApJ,862, 150, doi: 10.3847/1538-4357/aacb27Jackson, J. M., Geis, N., Genzel, R., et al. 1993, ApJ, 402,173, doi: 10.1086/172120Lau, R. M., Herter, T. L., Morris, M. R., Becklin, E. E., &Adams, J. D. 2013, ApJ, 775, 37,doi: 10.1088/0004-637X/775/1/37Mac Low, M.-M. 1999, ApJ, 524, 169, doi: 10.1086/307784Mapelli, M., Hayfield, T., Mayer, L., & Wadsley, J. 2012,ApJ, 749, 168, doi: 10.1088/0004-637X/749/2/168Mapelli, M., & Trani, A. A. 2016, A&A, 585, A161Marr, J. M., Wright, M. C. H., & Backer, D. C. 1993, ApJ,411, 667Marshall, J., Lasenby, A. N., & Harris, A. I. 1995, MNRAS,277, 594Mart´ın, S., Mart´ın-Pintado, J., Montero-Casta˜no, M., Ho,P. T. P., & Blundell, R. 2012, A&A, 539, A29,doi: 10.1051/0004-6361/201117268McKee, C. F., & Ostriker, J. P. 1977, ApJ, 218, 148,doi: 10.1086/155667Mezger, P. G., Zylka, R., Salter, C. J., et al. 1989, A&A,209, 337Mills, E. A. C., G¨usten, R., Requena-Torres, M. A., &Morris, M. R. 2013, ApJ, 779, 47,doi: 10.1088/0004-637X/779/1/47Mills, E. A. C., Togi, A., & Kaufman, M. 2017, ApJ, 850,192Monaghan, J. J., & Gingold, R. A. 1983, Journal ofComputational Physics, 52, 374,doi: 10.1016/0021-9991(83)90036-0Montero-Casta˜no, M., Herrnstein, R. M., & Ho, P. T. P.2009, ApJ, 695, 1477,doi: 10.1088/0004-637X/695/2/1477 urbulence in the CND17