A critical look at the merger scenario to explain multiple populations and rotation in iron-complex globular clusters
MMNRAS , 1–13 (2016) Preprint 6 November 2018 Compiled using MNRAS L A TEX style file v3.0
A critical look at the merger scenario to explain multiplepopulations and rotation in iron-complex globular clusters
Elena Gavagnin (cid:63) , Michela Mapelli and George Lake Institute for Computational Science, Centre for Theoretical Astrophysics and Cosmology, Universit¨at Z¨urich, Winterthurerstrasse 190,CH-8057 Z¨urich, Switzerland INAF-Osservatorio Astronomico di Padova, Vicolo dell’Osservatorio 5, I–35122, Padova, Italy
ABSTRACT
Merging has been proposed to explain multiple populations in globular clusters(GCs) where there is a spread in iron abundance (hereafter, iron-complex GCs). Bymeans of N-body simulations, we investigate if merging is consistent with the obser-vations of sub-populations and rotation in iron-complex GCs. The key parameters arethe initial mass and density ratios of the progenitors. When densities are similar, themore massive progenitor dominates the central part of the merger remnant and the lessmassive progenitor forms an extended rotating population. The low-mass progenitorcan become the majority population in the central regions of the merger remnant onlyif its initial density is higher by roughly the mass ratio. To match the radial distribu-tion of multiple populations in two iron-complex GCs ( ω Cen and NGC 1851), the lessmassive progenitor needs to be four times as dense as the larger one. Our merger rem-nants show solid-body rotation in the inner parts, becoming differential in the outerparts. Rotation velocity V and ellipticity ε are in agreement with models for oblaterotators with isotropic dispersion. We discuss several kinematic signatures of a mergerwith a denser lower mass progenitor that can be tested with future observations. Key words:
Galaxy: globular clusters - stars: kinematics and dynamics - methods:numerical - galaxies: star clusters
Over a quarter of the objects in Messier’s catalog are glob-ular clusters (GCs), yet we still do not know how they wereformed. For many decades GCs were described as stellar sys-tems with homogeneous chemical composition and no agespread, despite early data showing multiple populations inM22 and ω Cen (Geyer 1967; Canon & Stobie 1973; Harris1974; Freeman & Rodgers 1975; Hesser & Bell 1980).
Hubble Space Telescope data show a clear bifurcation ofcolour in the main sequence (MS) of ω Cen (Anderson 1997),with more recent data showing at least four distinct redgiant branches (RGBs, Lee et al. 1999; Pancino et al. 2000).Currently, most observed GCs show signatures of multiplepopulations, both in the Milky Way (Gratton et al. 2004;Carretta et al. 2007; Kayser et al. 2008; Anderson et al. 2009;Carretta et al. 2009b, 2010a; Carretta 2015; Pancino et al.2010; Milone et al. 2010, 2012, 2013) and in the MagellanicClouds (Milone et al. 2008).Most GCs contain stars with similar heavy-element (cid:63)
E-mail: [email protected] abundances (especially [Fe/H]), but large ( > . dex) star-to-star abundance variations for elements lighter than Si (e.g.Cohen 1978; Peterson 1980; Sneden et al. 1991; Kraft et al.1992; Gratton et al. 2001; Carretta et al. 2009b; Johnsonet al. 2015). Moreover, the variations of light-element abun-dances are anti-correlated with one another (e.g. the O − Naanti-correlation, Gratton et al. 2001). This phenomenology isgenerally considered to be due to internal enrichment by pro-ton capture H-burning reactions at high temperature (e.g.Gratton et al. 2004).A minority of GCs also show significant Fe abundancevariations. In particular, ω Cen (Norris & Da Costa 1995;Lee et al. 1999; Bellini et al. 2010; D’Orazi et al. 2011; Pan-cino et al. 2011), M22 (Hesser et al. 1977; Marino et al.2009; Lee 2015), M2 (Piotto et al. 2012; Lardo et al. 2013;Milone et al. 2015), M54 (Sarajedini & Layden 1995; Bel-lazzini et al. 2008; Carretta et al. 2010b), NGC 1851 (Yong& Grundahl 2008; Milone et al. 2009; Carretta et al. 2010c,2011), NGC 5286 (Nataf et al. 2013; Marino et al. 2015),NGC 5824 (Saviane et al. 2012; Da Costa et al. 2014), Terzan5 (Ferraro et al. 2009; Massari et al. 2014) and M19 (John-son et al. 2015) are labelled as ’iron-complex’ GCs, because c (cid:13) a r X i v : . [ a s t r o - ph . GA ] J un E. Gavagnin et al. they have (i) a spread in [Fe/H] exceeding ∼ . dex, (ii)multiple photometric sequences, and (iii) a significant abun-dance spread for both light and heavy elements (Johnsonet al. 2015). Iron-complex GCs differ from other GCs in sev-eral ways. In most GCs, the stellar population showing noenrichment by proton capture accounts for about one thirdof the total GC mass, with little spread among GCs (Car-retta et al. 2009a). In contrast, in the iron-complex GCs,the ratio between the metal-poor and the metal-rich popu-lation changes from cluster to cluster. For example, in M19the metal-poor component is ∼ % of the entire popula-tion (Johnson et al. 2015), whereas ∼ % of spectroscopi-cally studied stars in M2 belong to the metal-poor compo-nent (Milone et al. 2015). Moreover, in the vast majority ofGCs, the proton-capture enriched population is more radi-ally concentrated than the most numerous one. In the iron-complex GCs the metal-poor population can be either moreconcentrated ( ω Cen, Bellini et al. 2009) or less concentrated(NGC 1851, Carretta et al. 2011) than the metal-rich one.Several theoretical models have been proposed to ex-plain the multiple populations (Bastian et al. 2013; see Ren-zini 2008 for a review). A first class of models appeals tomultiple star-formation events. After first population starsform out of pristine, metal-poor gas, the second populationof stars might form from either the ejecta of asymptotic giantbranch (AGB) stars (D’Ercole et al. 2008) or fast rotatingmassive stars (FRMS, Decressin et al. 2007). In the ‘AGBscenario’, winds and supernovae (SNe) of the first populationevacuate the residual gas. After ∼ , low-velocity windsfrom AGBs enriched in He and s-process elements start ac-cumulating at the centre and form the second population.However, the predicted mass of the second population is anorder of magnitude lower than what is observed, requiringa top-heavy first population initial mass function and anunusually efficient second population star formation.A second model, called the ‘early disc accretion model’(Bastian et al. 2013), proposes that the two populationsformed during the same star formation episode, but under-went different chemical enrichment. This model requires veryfast mass segregation and gas evaporation. With rapid masssegregation, the most massive stars sink to the centre wherehigh-mass stars in interacting binaries eject the primary’sHe-enriched envelope. This material pollutes the circumstel-lar discs of low mass stars that are still accreting, so theywill grow in mass thanks to these ejecta from more mas-sive (but still same-generation) stars. The main drawbacksof this model are disc lifetime and uniformity of enrichment.Even if circumstellar discs survive for the required − (De Marchi et al. 2013), the “rainfall” of enriched materialonto them is unlikely to be uniform (Kruijssen 2014).All the aforementioned scenarios are aimed to explainmultiple populations with no or negligible iron spread, whilethey fail to reproduce the [Fe/H] variations observed in theiron-complex GCs. So far, the only proposed scenario thatcan naturally account for a metallicity spread is the mergerbetween GCs (Sugimoto & Makino 1989; Makino et al. 1991; Recent studies highlight the possibility that the [Fe/H] spreadis spurious, at least in some GCs (e.g. M22, Mucciarelli et al.2015), because spectroscopically derived Fe abundances might beinaccurate due to non-local thermodynamical equilibrium effects. van den Bergh 1996; Amaro-Seoane et al. 2013; Pasquato& Chung 2016). In this scenario, the different metallicitiesare signatures of the progenitors and can be used as a tagto make predictions about the distribution and dynamicswithin the final merger remnant. Iron abundance is, in thisrespect, a good tag to identify uniquely the different popu-lations.The merger scenario might be consistent with the old-est metal-rich stars in ω Cen being a few Gyrs older thanthe oldest metal-poor stars (Villanova et al. 2014), a circum-stance that is against the predictions of self-enrichment sce-narios. Furthermore, a merger can explain the kinematicaldifferences in the velocity dispersion of the calcium-weak andcalcium-strong RGB stars in M22 (Lee 2015). The mergerscenario has been proposed also for NGC 1851, where themost metal-rich population is less concentrated than themetal-poor one (van den Bergh 1996; Carretta et al. 2010c,2011; Bekki & Yong 2012).Another advantage of the merger scenario is that it canaccount for signatures of rotation in GCs, which have beenobserved in several GCs with multiple populations, bothwith (e.g. ω Cen, M2, M22, M54, NGC1851, Pancino et al.2007; Lee 2015; Pryor et al. 1986; Kimmig et al. 2015; Bel-lazzini et al. 2012; Lardo et al. 2015) and without (e.g. Fabri-cius et al. 2014) a metallicity spread. If the two progeni-tors have non-zero relative orbital angular momentum, themerger remnant will likely preserve a signature of rotationin the merger plane. However, there is no evidence that GCswith a metallicity spread (the best-candidate merger rem-nants) have systematically higher rotation than the otherGCs. Moreover, other physical mechanisms can account forrotation in GCs (e.g. Mastrobuono-Battisti & Perets 2013;Vesperini et al. 2014; Bianchini et al. 2015).The main problem for the merger scenario is that twoGCs are expected to merge only if their relative velocity issmaller than (or of the same order of magnitude as) theirvelocity dispersion. The largest GC in the galaxy, ω Cen,has a dispersion of ∼ km s − , with typical values being ∼ − km s − . The relative velocities of current GCs in theMilky Way halo are at least one order of magnitude largerthan these values. This means that a merger between twoGCs that are in the halo of our galaxy is extremely un-likely. Two GCs will have a sufficiently low relative velocityto merge only if they formed in a small dwarf galaxy or in thesame molecular cloud. However, if the two progenitor clus-ters formed in the same molecular cloud and merged slightlyafter their formation, it is difficult to explain why the twopopulations have a different proton-capture enrichment andeven a different metallicity. Therefore, GCs in small dwarfgalaxies represent the most likely scenario where GC merg-ers will produce clusters that have a spread in metallicity.We take a critical approach to the merging scenario byexamining how the initial mass and density ratios of the pro-genitors affect the distribution and concentration of the sub-populations in the remnant (Section 3.1). Moreover, we alsoexamine the rotation signature of the merger product andwe show that the profile of rotation is related to the initialdensity ratio of the progenitors (Section 3.2). In the event ofequal-mass mergers, we expect that the denser initial pro-genitor will be more centrally concentrated in the remnant.In the case of unequal-mass mergers, the more massive pro-genitor will be closer to the centre than the less massive MNRAS , 1–13 (2016) erging to explain multiple populations and rotation in GCs progenitor and hence be more concentrated. We examinehow the density ratio can counter the mass ratio.This paper is organised as follows. In Section 2, we de-scribe the numerical tools and the initial conditions adopted.In Section 3, we present the main results of this work. Sec-tion 4 is dedicated to the discussion and conclusions. We used the starlab public software environment (Porte-gies Zwart et al. 2001) ported to GPUs (Gaburov et al. 2009)to run the simulations. To investigate the role of the relativemasses and densities of the progenitors, we performed a gridof simulations varying the mass ratio, i.e. M / M (where M is the mass of GC1 and M is the mass of GC2) and thedensity ratio, i.e. ρ / ρ (where ρ and ρ are the densitiesmeasured at the virial radius of GC1 and GC2, respectively).We consider mass ratios of 1, 2, 4, with density ratios of0.25, 0.5, 1, 2, 4. The range is motivated by the ratio of thepopulations in GCs (Johnson et al. 2015; Milone et al. 2015)and the absence of strong correlations between luminosityand density in present-day GCs (Harris 1996). The GCs aremodelled as non-rotating spherical King profiles (King 1966)with central dimensionless potential W = (this sets thecore radius R c = . R V ). The second GC (GC2) is alwayscomposed of 20 000 particles of equal mass m ∗ = (cid:12) for atotal mass of M (cid:12) . Its virial radius R V = pc. The firstGC (GC1) is varied to set the mass ratio and density ratio.To double (or quadruple the mass) of GC1, we double (orquadruple) the number of particles.The density ratio is set by adjusting the virial radius ofGC1, e.g. in the run M ρ the GC1 has twice the mass asGC2 and R V of GC1 is ∼ . × larger than the one of GC2,so that the density ratio between the two clusters is 1. Wenote that, assuming a fixed value for W , the density ratio isthe same at every fiducial radius, i.e. the core radius ( R c ),the tidal radius ( R t ) and the virial radius ( R V ).To prevent strong encounters and binary formation, weadopt a gravitational softening ε = . R V of the progenitorwith the smallest radius. The initial binary fraction is zeroand binaries do not form with this softening. We omit stellarand binary evolution to minimise the amount of free param-eters in these models. Stellar and binary evolution might af-fect the structural properties of GCs (Chernoff & Weinberg1990; Mapelli et al. 2013; Mapelli & Bressan 2013; Traniet al. 2014; Gieles 2013; Sippel et al. 2012) and will be con-sidered in a follow up study. Stars initially belonging to eachof the two progenitors are “tagged” with a different metal-licity flag. Initial conditions (ICs) are summarised in Table1. The two GCs are initially set on a parabolic orbit. Todefine the parabolic orbit we fixed the minimum encounterdistance (in the point-mass assumption), i.e. the pericentre r peri , to be half the sum of the virial radii of the two progen-itors GC1, GC2 [ r peri = . ( R V + R V ) ]. The initial distance D between the progenitors is four times the maximum valuebetween R t , and R t , , where R t , and R t , indicate the tidalradius of GC1 and GC2, respectively. The initial relative ve-locity is then calculated as the escape velocity at the initialposition.We choose a parabolic orbit because it is a represen- tative case for mergers (Alladin 1965; White 1978). Hyper-bolic encounters (with relative velocity much larger than theGC velocity dispersion) are the most common, as the phasespace for encounters increases with the cube of the velocity ofencounter and the cube of the impact parameter. However,the probability of merging encounters is sharply truncated(by failure to merge) when the orbits become very weaklyhyperbolic. In contrast, the two GCs will merge on a shortertimescale if they are on a bound orbit, but bound orbits areassociated with smaller values of the velocity. As we men-tioned in the introduction, the main drawback of the mergerscenario is that the observed relative velocities between GCsare generally larger than the value needed for a merger tobe successful. Thus, we consider bound orbits very unlikely.In summary, a parabolic orbit is representative of the mostlikely orbits leading to a merger.The half-mass relaxation time is (Spitzer 1987) t rlx ∼ × yr (cid:18) M M (cid:12) (cid:19) / (cid:18) R hm (cid:19) / (cid:18) m (cid:12) (cid:19) − (cid:18) ln Λ (cid:19) − , (1)where R hm is the initial half-mass radius, M the total mass, m is the particle mass and ln Λ is the Coulomb logarithm (sethere by the system size and gravitational softening). For ourprogenitors, the relaxation timescales are between and . . The initial crossing time at the virial radius inthe equal mass, equal density progenitors is ∼ Myr andscales as ρ − / . We run our simulations for 550 Myr. Thisis less than one half-mass relaxation timescale characteristicof the merger product in all cases, but two-body encountershave likely contributed to isotropising the velocities in theremnant. We examine the relative concentration and rotation of thetwo different populations in the merger remnant.
We plot the relative concentration using normalised densityprofiles of the sub-populations (i.e. each density profile isdivided by its progenitor’s mass). Figure 1 shows the densityprofile of GC1 and GC2 in green and in magenta respectively(where M GC (cid:62) M GC ).We plot the profiles of nine selected runs. The profilesare at time ∼ Myr since the beginning of the simulation.We see that the final density profiles of the merger remnantsare consistent with a single King model profile, although thetwo populations have different densities in the central re-gions. Depending on the run, we note that at small radii thenormalised density of GC1 members can be higher than thatof GC2 members or viceversa . This suggests that the initialmass and density ratios affect the relative central densityof the two populations in the final merger remnant (Figure1). Despite the normalisation to the progenitor’s mass, atlarge radii one curve is below the other in several panels.For example, in several plots of Figure 1 the magenta curveis below the green one (see especially the bottom right panel:
MNRAS000
MNRAS000 , 1–13 (2016)
E. Gavagnin et al.
Table 1.
Initial conditions of the simulations. Run (column 1): identifying name of the run, ‘M’ stands for mass ratio ( M / M ), ‘ ρ ’for density ratio ( ρ / ρ ), both followed by the values assumed, e.g. M4 ρ N (column 3): number of particles; M (column 4): total mass of the progenitor; R V (column 6): virial radius; D (column 7): initial distance between the progenitors’ centres of mass; r peri (column 9): orbital pericentre; V (column 10): initial relativevelocity. Run N M R V D r peri V [ M (cid:12) ] [pc] [pc] [pc] [km s − ]GC2 20k ρ ρ ρ ρ ρ ρ · ρ · ρ · ρ · ρ · ρ ·
10 176 7 7.03M4 ρ · ρ · ρ · ρ · since the profiles are normalised to the mass of the progeni-tors, this is a clear signature of mass loss during the mergingprocess).Figure 2 is a colour map of the relative concentrationof the two progenitors, defined as ‘ log ( R hm2 / R hm1 ) ’, i.e. thelogarithm of the ratio between the half-mass radius of GC2and GC1, in the initial conditions and at the end of thesimulations, for the whole grid of runs. The plot on the left-hand side in Figure 2 shows the ratio between the half-massradius of GC2 and GC1 in the ICs, the plot on the right-hand side shows the same quantity after the merger. Fromthe right panel in Figure 2, we see that when the initialdensities are equal, the more massive progenitor dominatesthe central part of the merger remnant and the less massiveprogenitor is more extended in the merger remnant. If theprogenitors have equal masses, the denser progenitor is moreconcentrated in the remnant. To compensate for an unequalmass ratio, the less massive progenitor must have a densitylarger by roughly the factor by which its mass is lower. If thesmaller mass progenitor is / A as massive, its initial densitymust be A times greater or alternatively, its radius must be A − / as large as the more massive one.In Figure 3 we compare the number ratio of sub-populations ( N / N ) in our simulated GCs with the obser-vations. Specifically, we plot the ratio of the minority ( N )to the majority ( N ) population against the radial distance from the centre, normalised to the half-mass (or half-light)radius. Observational data of three GCs are compared withour simulations: in ω Cen the metal-rich population is themost centrally concentrated and is the minority population(Bellini et al. 2009), in NGC 1851 (Carretta et al. 2011)and M22 (Carretta et al. 2011) the metal-poor populationis the more centrally concentrated (note that crowding pre-vents observing the very central regions of NGC 1851). InM22 the metal-rich population is the minority, while in NGC1851 the metal-poor population is the minority.The two runs shown in Figure 3 (M2 ρ ρ ω Cen and metalpoor in NGC 1851. We adopt a different progenitors model(with equal density) for M 22, where the minority population(metal-rich) is less concentrated.
Rotation is observed in nearby GCs (Anderson & King 2003;van den Bosch et al. 2006; Bellazzini et al. 2012; Lardo et al.2015; Fabricius et al. 2014), which can arise from a variety
MNRAS , 1–13 (2016) erging to explain multiple populations and rotation in GCs R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M1 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M1 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M1 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M2 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M2 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M2 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M4 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M4 ρ R (pc)10 -8 -7 -6 -5 -4 -3 -2 -1 D e n s i t y / M ( / p c ) M4 ρ Figure 1.
Normalised density profiles of the two populations in the final merger remnant. Note that the profiles all look like smoothKing models. Solid green line refers to GC1, dashed magenta line to GC2. Each profile has been normalised by the mass of the associatedprogenitor. The codename on the top of each plot refers to the run considered: ‘M’ stands for mass ratio M / M followed by its value,‘ ρ ’ for density ratio ρ / ρ followed by its value. From top to bottom the mass ratio increases by a factor of 2 every row and from left toright the density ratio increases by a factor of 4 every column. of mechanisms (Bertin & Varri 2008; Varri & Bertin 2012;Bianchini et al. 2013; Vesperini et al. 2014). While there isno connection demonstrated between rotation and multiplepopulations, Amaro-Seoane et al. (2013) pointed out that ω Cen, M 22, and NGC 2419 are among the most flattenedclusters.Flattening has been detected in several galactic GCs(White & Shawl 1987; Chen & Chen 2010) and couldbe explained by several physical factors, such as pressureanisotropy or external tides (van den Bergh 2008). An-other possible justification for the flattening is the inter-nal rotation of GCs (Fabricius et al. 2014). A correlationbetween flattening and iron-complex multiple populationswould favour the merger scenario.In this section, we look at the detailed kinematics ofour merger remnants, as a function of mass and density ra-tios. We want to quantify their amount of rotation and seewhether their degree of flattening correlates with rotation.All of our merger remnants have rotation, as a conse-quence of the parabolic orbits of their progenitors. In Figs. 4, 5 and 6, we show velocity maps for the complete rangeof initial mass ratios and the limiting cases of density ratios ρ / ρ = . , . In all cases, we plot line-of-sight velocitiesfor an observer sitting on the mid-plane perpendicular to therotation axis. For comparison with the observations (Fabri-cius et al. 2014), we used a Gaussian filter to create thevelocity map, progressively zoomed from left to right. Eventhe largest spatial scales of the final merger state (left-handcolumns) show a clear flattening and the two populationshave similar properties in configuration and velocity space.These maps illustrate some important trends: the rota-tion within 5 pc is generally solid body, it becomes differen-tial at − pc, and the rotation is cylindrical everywhere.The similarity of the maps shows that these features arecommon to all our simulations. Solid-body rotation is themost probable distribution function (maximal entropy) forthe relaxed core of a rotating N-body system (Lynden-Bell1967; Lightman & Shapiro 1978). Observations of clustersalso show solid body rotation over most of the half-mass MNRAS000
Normalised density profiles of the two populations in the final merger remnant. Note that the profiles all look like smoothKing models. Solid green line refers to GC1, dashed magenta line to GC2. Each profile has been normalised by the mass of the associatedprogenitor. The codename on the top of each plot refers to the run considered: ‘M’ stands for mass ratio M / M followed by its value,‘ ρ ’ for density ratio ρ / ρ followed by its value. From top to bottom the mass ratio increases by a factor of 2 every row and from left toright the density ratio increases by a factor of 4 every column. of mechanisms (Bertin & Varri 2008; Varri & Bertin 2012;Bianchini et al. 2013; Vesperini et al. 2014). While there isno connection demonstrated between rotation and multiplepopulations, Amaro-Seoane et al. (2013) pointed out that ω Cen, M 22, and NGC 2419 are among the most flattenedclusters.Flattening has been detected in several galactic GCs(White & Shawl 1987; Chen & Chen 2010) and couldbe explained by several physical factors, such as pressureanisotropy or external tides (van den Bergh 2008). An-other possible justification for the flattening is the inter-nal rotation of GCs (Fabricius et al. 2014). A correlationbetween flattening and iron-complex multiple populationswould favour the merger scenario.In this section, we look at the detailed kinematics ofour merger remnants, as a function of mass and density ra-tios. We want to quantify their amount of rotation and seewhether their degree of flattening correlates with rotation.All of our merger remnants have rotation, as a conse-quence of the parabolic orbits of their progenitors. In Figs. 4, 5 and 6, we show velocity maps for the complete rangeof initial mass ratios and the limiting cases of density ratios ρ / ρ = . , . In all cases, we plot line-of-sight velocitiesfor an observer sitting on the mid-plane perpendicular to therotation axis. For comparison with the observations (Fabri-cius et al. 2014), we used a Gaussian filter to create thevelocity map, progressively zoomed from left to right. Eventhe largest spatial scales of the final merger state (left-handcolumns) show a clear flattening and the two populationshave similar properties in configuration and velocity space.These maps illustrate some important trends: the rota-tion within 5 pc is generally solid body, it becomes differen-tial at − pc, and the rotation is cylindrical everywhere.The similarity of the maps shows that these features arecommon to all our simulations. Solid-body rotation is themost probable distribution function (maximal entropy) forthe relaxed core of a rotating N-body system (Lynden-Bell1967; Lightman & Shapiro 1978). Observations of clustersalso show solid body rotation over most of the half-mass MNRAS000 , 1–13 (2016)
E. Gavagnin et al. /M ρ / ρ Initial conditions < -0.60-0.40-0.2000.200.40> 0.60 l o g ( R h m / R h m ) /M ρ / ρ After merger < -0.60-0.40-0.2000.200.40> 0.60 l o g ( R h m / R h m ) Figure 2.
Colour map of initial (left) and final (right) ratio between the half-mass radius of GC2 and GC1 ( log ( R hm2 / R hm1 ) ). The x − and y − axis are the initial mass ratio and the initial density ratio of the two progenitors. The colour map quantifies the relativeconcentration of the two populations (in logarithmic scale), meaning the ratio of the two half-mass radii i.e. R hm2 / R hm1 . If the logarithmof this value is negative (magenta) GC2 is more centrally concentrated in the merger product; if it is positive (green), GC1 is morecentrally concentrated. The blue line marks the boundary between where GC1 is more centrally concentrated (green) and the situationswhere GC2 is more centrally concentrated (magenta). In both plots, the actual data grid at which R hm2 / R hm1 is evaluated is marked withblack dots, the colour map is then smoothed via interpolation in order to highlight the general trend. Note that the x − and y − axis areeffectively logarithmic. radius and differential rotation outside (Meylan & Mayor1986; Bianchini et al. 2013; Fabricius et al. 2014).Figure 7 shows the line-of-sight rotation profiles of allthe simulations for an observer sitting on the rotation plane.As in the velocity maps, we see that the inner rotationis solid-body, then it becomes differential at 5-10 pc. Thesolid-body rotation region is more extended in the runs withhigher mass ratio; the angular momentum of the less mas-sive object is preferentially deposited in the outskirts of theremnant.At the half mass radius, the merger remnant exhibitssolid-body or differential rotation depending on the initialdensity ratio between the progenitors. In Figure 8, we exam-ine the ratio of the rotation velocity at the half-mass radiusto the maximum rotation velocity V Rhm / V max , as a function ofdensity ratio. For equal-mass ratios the quantity V Rhm / V max is almost constant with respect to the density ratio (toppanel of Figure 8). Therefore each of these model clustershave transitioned from solid-body to differential rotation bythe half mass radius. In contrast, the trend for unequal massratios provides an interesting test of the model. When theless massive progenitor is less dense, it deposits its angularmomentum in the outer parts. In contrast, small-mass pro-genitors with larger density burrow into the centre. Whenthe minority population is more concentrated, the rotationcurve will peak at roughly the half-mass radius, whereaswhen the minority population is less concentrated, the peakoccurs outside the half-mass radius. Therefore, for unequal- mass ratios, V Rhm / V max decreases for increasing values of thedensity ratio.In order to compare the outcomes of our simulationswith observations, we study now the ( V / σ , ε ) diagram,which relates the ratio of the rotation velocity V and randommotion σ to the ellipticity ε , which measures the flattening.The expectation for isotropic rotators are derived fromthe tensor virial theorem (Chandrasekhar 1969). The ro-tation velocity is the square root of the mass weightedstreaming velocity squared. The velocity dispersion is theunordered kinetic energy. If the mass is stratified on concen-tric similar ellipsoids, the density profile drops out (Roberts1962; Chandrasekhar & Lebovitz 1962) and the ratio of theordered kinetic energy to the unordered one (or its squareroot V / σ ) is a function only of the ellipticity ε (Binney1978). The application to elliptical galaxies is straightfor-ward since V , σ and ε are all nearly constant with radius(Emsellem et al. 2007).For GCs, ε has a greater variation with radius and V is rising with an asymptote at a radius beyond the obser-vations. Hence, we look at how well ‘proxy’ and ‘measured’rotations relate to one another in the simulated merger rem-nants. As always with proxies and dimensionless parametersthat vary with radius, the results will be mixed.In Figure 9 we plot the ( V / σ , ε ) diagram, following theprescription of Fabricius et al. (2014) as proxies for V and ε , including both data from our simulations and observedGCs. Fabricius et al. (2014) fit a plane V ( x , y ) = ax + by + V sys (where V sys is the systemic velocity) to the velocity fields to MNRAS , 1–13 (2016) erging to explain multiple populations and rotation in GCs -1 R/R h N / N M2 ρ NGC1851 (Carretta+11)Omega Cen (Bellini+09) -1 R/R h N / N M2 ρ M22 (Lee+15)
Figure 3.
Ratio of the minority ( N ) to the majority ( N ) popula-tion versus the radial distance from the centre. The blue solid lineindicates our simulated models M2 ρ ρ ω Cen, Carretta et al. 2011 for NGC 1851, and Lee2015 for M 22). N / N is normalised to the half-mass radius and tothe half-light radius for the simulations and for the observations,respectively. determine the central velocity gradient, || ∇ V || = √ a + b .They take velocity dispersions σ and half-light radii R hl from Harris (1996) to create a proxy for rotational veloc-ity ∇ V · R hl , and find V / σ increasing with ellipticity. In ourcase, we define V in similar way ( ∇ V · R ) leaving though R as free parameter, with the intent to explore how this proxyfor V depends on the radius used to define it. ∇ V is also cal-culated within the radius considered each time. Specifically,in Figure 9 we considered three different values of radius R ,that are R hm , R hm / , and R hm / . Our choice is justified byFigure 7, where the solid-body rotation shifts to differentialrotation at radii varying from ∼ . R hm to ∼ . R hm . σ in our case is just the line-of-sight velocity dispersion.As for the ellipticity, we follow the prescriptions found inFabricius et al. (2014) and calculate ellipticity values ( ε )from the eigenvalues ( λ , λ ) of the covariance matrix of stellar positions (within the relevant radial dimension), i.e. ε = (cid:112) − λ / λ .Figure 9 shows that the result strongly depends on thechoice of radius. The V / σ ratio increases with ellipticity,but ellipticities and V / σ both increase with radius. If welooked at Figure 7, we might guess that something close to ∇ V · R hm would be the best proxy and certainly would nothave guessed that the plot using ∇ V · R hm / would look mostlike the oblate rotator (dashed line in the plot) and wouldbe most in agreement with the data from Fabricius et al.(2014) . Thus, the choice of rotational velocity in a ( V / σ , ε )diagram is not unique .Another possibility is to take as the rotation velocitythe maximum rotational velocity. In Figure 10, we plot theellipticity ε of our simulated merger remnants versus V max / σ ,where V max is the maximum rotational velocity (see Fig-ure 7). The result (shown as star symbols in Figure 10) com-pares favourably with the oblate rotator curve and observa-tions. Having set the initial orbits to parabolic, the valuesof ( V max / σ , ε ) for the simulated GCs are all in the sameportion of the oblate rotators curve. With time, the mergerremnants will radiate away angular momentum through two-body encounters (Fall & Frenk 1985). This will make themslide down on the curve to lower ellipticity and V / σ val-ues, closer to the observational data, because rotation andellipticity will decrease significantly as soon as the systemrelaxes and the two populations mix completely (velocitieswill isotropise and angular momentum will diminish).In Figure 9 and Figure 10, we plot not only the obser-vational sample of Fabricius et al. (2014), whose 11 GCs donot show any significant spread in Fe abundance, but alsodata of some iron-complex GCs (M 22 and M 54 from Bel-lazzini et al. 2012, M 2 from Pryor et al. 1986, ω Cen fromBianchini et al. 2013, NGC 1851 from Lardo et al. 2015).To derive the value of V max / σ for the iron-complex GCs,we use the double mean velocity amplitude (i.e. A rot ) whichis considered a good representation of V max (Pancino et al.2007). From the kinematical point of view, the iron-complexGCs for which V / σ and ε are available do not stand out incomparison with the sample of Fabricius et al. (2014). In this section, we discuss the results of our simulations inlight of the observational properties of iron-complex GCs.We focus on GCs with multimodal iron-complex abundancesbecause they have unique tags that can be mapped to pos-sible progenitors.In the merger scenario, we find that the minority popu-lation is less centrally concentrated unless the initial densityof the less massive progenitor is greater by more than themass ratio. In M22, the minority is metal poor and extended.The distribution compares well to equal density progenitorswith a mass ratio of 2:1. In ω Cen and NGC1851 the lessmassive population is more centrally concentrated than the The proxy for ∇ V adopted by Fabricius et al. (2014) wouldalways be higher than the true rotation velocity at the radius R ,because it comes from the best linear fit to the velocities within R and the second derivative of V with respect to R is negative (therotation curve is flattening).MNRAS000
Ratio of the minority ( N ) to the majority ( N ) popula-tion versus the radial distance from the centre. The blue solid lineindicates our simulated models M2 ρ ρ ω Cen, Carretta et al. 2011 for NGC 1851, and Lee2015 for M 22). N / N is normalised to the half-mass radius and tothe half-light radius for the simulations and for the observations,respectively. determine the central velocity gradient, || ∇ V || = √ a + b .They take velocity dispersions σ and half-light radii R hl from Harris (1996) to create a proxy for rotational veloc-ity ∇ V · R hl , and find V / σ increasing with ellipticity. In ourcase, we define V in similar way ( ∇ V · R ) leaving though R as free parameter, with the intent to explore how this proxyfor V depends on the radius used to define it. ∇ V is also cal-culated within the radius considered each time. Specifically,in Figure 9 we considered three different values of radius R ,that are R hm , R hm / , and R hm / . Our choice is justified byFigure 7, where the solid-body rotation shifts to differentialrotation at radii varying from ∼ . R hm to ∼ . R hm . σ in our case is just the line-of-sight velocity dispersion.As for the ellipticity, we follow the prescriptions found inFabricius et al. (2014) and calculate ellipticity values ( ε )from the eigenvalues ( λ , λ ) of the covariance matrix of stellar positions (within the relevant radial dimension), i.e. ε = (cid:112) − λ / λ .Figure 9 shows that the result strongly depends on thechoice of radius. The V / σ ratio increases with ellipticity,but ellipticities and V / σ both increase with radius. If welooked at Figure 7, we might guess that something close to ∇ V · R hm would be the best proxy and certainly would nothave guessed that the plot using ∇ V · R hm / would look mostlike the oblate rotator (dashed line in the plot) and wouldbe most in agreement with the data from Fabricius et al.(2014) . Thus, the choice of rotational velocity in a ( V / σ , ε )diagram is not unique .Another possibility is to take as the rotation velocitythe maximum rotational velocity. In Figure 10, we plot theellipticity ε of our simulated merger remnants versus V max / σ ,where V max is the maximum rotational velocity (see Fig-ure 7). The result (shown as star symbols in Figure 10) com-pares favourably with the oblate rotator curve and observa-tions. Having set the initial orbits to parabolic, the valuesof ( V max / σ , ε ) for the simulated GCs are all in the sameportion of the oblate rotators curve. With time, the mergerremnants will radiate away angular momentum through two-body encounters (Fall & Frenk 1985). This will make themslide down on the curve to lower ellipticity and V / σ val-ues, closer to the observational data, because rotation andellipticity will decrease significantly as soon as the systemrelaxes and the two populations mix completely (velocitieswill isotropise and angular momentum will diminish).In Figure 9 and Figure 10, we plot not only the obser-vational sample of Fabricius et al. (2014), whose 11 GCs donot show any significant spread in Fe abundance, but alsodata of some iron-complex GCs (M 22 and M 54 from Bel-lazzini et al. 2012, M 2 from Pryor et al. 1986, ω Cen fromBianchini et al. 2013, NGC 1851 from Lardo et al. 2015).To derive the value of V max / σ for the iron-complex GCs,we use the double mean velocity amplitude (i.e. A rot ) whichis considered a good representation of V max (Pancino et al.2007). From the kinematical point of view, the iron-complexGCs for which V / σ and ε are available do not stand out incomparison with the sample of Fabricius et al. (2014). In this section, we discuss the results of our simulations inlight of the observational properties of iron-complex GCs.We focus on GCs with multimodal iron-complex abundancesbecause they have unique tags that can be mapped to pos-sible progenitors.In the merger scenario, we find that the minority popu-lation is less centrally concentrated unless the initial densityof the less massive progenitor is greater by more than themass ratio. In M22, the minority is metal poor and extended.The distribution compares well to equal density progenitorswith a mass ratio of 2:1. In ω Cen and NGC1851 the lessmassive population is more centrally concentrated than the The proxy for ∇ V adopted by Fabricius et al. (2014) wouldalways be higher than the true rotation velocity at the radius R ,because it comes from the best linear fit to the velocities within R and the second derivative of V with respect to R is negative (therotation curve is flattening).MNRAS000 , 1–13 (2016) E. Gavagnin et al.
Figure 4.
Line-of-sight velocity maps at different scales at t= 550
Myr for the case with equal mass ratio between the progenitorsand ρ / ρ =0.25 (top row), ρ / ρ =4 (bottom row). From left to right, we zoom in the central parts of the remnant. The largest scales(left-hand columns) show a clear flattening. Examining these colour maps, the rotation within 5 pc is generally solid body (colour ischanging), it becomes differential at 5-10 pc (the colour stays constant outside this radius in the rotation plane) and it is cylindricaleverywhere (weak or no colour trend vertically). The similarity of all the maps reveals that these are common features of mergers. Figure 5.
The same as Figure 4, but for the case with 2:1 mass ratio between the progenitors and ρ / ρ =0.25 (top row), ρ / ρ =4(bottom row). MNRAS , 1–13 (2016) erging to explain multiple populations and rotation in GCs Figure 6.
The same as Figure 4, but for the case with 4:1 mass ratio between the progenitors and ρ / ρ =0.25 (top row), ρ / ρ =4(bottom row). majority population (Bellini et al. 2009). Merging only caresabout metallicity if there is a correlation between metallicityand mass or density. The minority is metal rich in ω Cen,while it is metal poor in NGC 1851. In light of our results,this means these GCs can be the result of a merger only ifthe less massive progenitor was the denser one. These trendsare best fit when the mass ratio is 2:1 and the less massiveprogenitor is four times as dense as the more massive one.It would be instructive to check with observationswhether the less massive progenitors are denser than moremassive ones. Figure 11 shows the relation between lumi-nosity (as a proxy for mass) and half-light radius in present-day GCs, from the catalogue of Harris (1996). This figureshows that there is no correlation between luminosity (hence,mass) and size in present-day Milky Way GCs. From thisfact, we cannot conclude very dense but small mass pro-genitors would be common, if merger progenitors were likepresent-day GCs. However, we also warn that consideringpresent-day GCs as representative of merger progenitors israther hazardous.The kinematical signatures of the merger remnant aresimilar to those observed in GCs. In our simulated remnants:1) the velocity dispersion is isotropic, 2) the merger productrotates close to solid body in the inner parts, then becomesdifferential, 3) rotation is cylindrical, 4) at the half mass ra-dius, the merger remnant exhibits solid-body or differentialrotation depending on the initial density ratio between theprogenitors, 5) the flattening of the remnant is consistentwith rotation. Both ε and V vary over radius, so definingappropriate values for a ( V / σ , ε ) plot is difficult. Differentchoices move points around in that plot, but the correlationbetween flattening and rotation in the remnants is similarto the expectations from the tensor virial theorem (Binney1978). As we already anticipated in the introduction, the mostsevere drawback of the merger scenario is that the relativevelocity between two clusters must be sufficiently low tomerge. Here ‘sufficiently low’ means that their relative ve-locity cannot be much larger than their velocity dispersion.The velocity dispersion of GCs is ≈ per cent of the veloc-ity dispersion of stars in the field of our Galaxy. This meansthat GCs move too fast to merge in our present-day Galaxy.Several studies propose that a sub-population of GCswere the nuclei of dwarf galaxies, with ω Cen as prototype(Majewski et al. 2000; Carraro & Lia 2000). If one GC were anucleus, the inspiral of a second GC would create conditionssimilar to an unequal mass merger.GCs can sink toward the centre of the host galaxy bydynamical friction. The dynamical friction timescale scalesas the inverse of the mass of the GC. Thus, the smallerthe GC, the longer it takes for it to sink to the centre bydynamical friction. For example, an object that has a massof ∼ per cent of the total mass of the host galaxy will spiralinto the centre by dynamical friction in roughly a dynamicaltime (Binney & Tremaine 2008).The smallest dwarf galaxy in the Local Group with GCsis Fornax, with five clusters (Larsen et al. 2012). The mostmassive among these GCs has not yet sunk into the centreby dynamical friction (Read et al. 2006). Thus, even Fornaxfailed to promote mergers or create a nucleus from its mostmassive GC.The Sagittarius dwarf galaxy is more promising (Grat-ton et al. 2012). At least five Milky Way GCs are thoughtto have been part of Sagittarius (Law & Majewski 2010).The velocity dispersion of Sagittarius is ∼
20 km s − . Thus,parabolic encounters between GCs would be rare, but notimpossible. Sagittarius does have a nuclear cluster. With avelocity dispersion of ∼ km s − , Sagittarius has a mass MNRAS000
20 km s − . Thus,parabolic encounters between GCs would be rare, but notimpossible. Sagittarius does have a nuclear cluster. With avelocity dispersion of ∼ km s − , Sagittarius has a mass MNRAS000 , 1–13 (2016) E. Gavagnin et al. -4-2024 M1 ρ -4-2024 M1 ρ -4-2024 M1 ρ -4-2024 M1 ρ -4-2024 M1 ρ -4-2024 M2 ρ -4-2024 V l o s ( k m / s ) M2 ρ -4-2024 M2 ρ -4-2024 M2 ρ -4-2024 M2 ρ -4-2024 M4 ρ -4-2024 M4 ρ -4-2024 M4 ρ -4-2024 M4 ρ
20 15 10 5 0 5 10 15 20
X (pc) -4-2024 M4 ρ Figure 7.
Line-of-sight velocity profile for all the simulationsfor an observer in the mid plane ( V los ). The half-mass radius isdenoted by an arrow in each plot. The thicker red horizontal linesdivide the panels by mass-ratio (the top group has M / M = , thecentral M / M = , the bottom one M / M = ). Every group hasplots for the 5 density ratios considered. V R h m / V m a x M1 V R h m / V m a x M2 ρ /ρ V R h m / V m a x M4 Figure 8.
Ratio of the rotation velocity at the half-mass radiusto the maximum rotation velocity, V Rhm / V max , as a function ofdensity ratio. From top to bottom: each panel refers to GCs withmass ratio M / M = , x − axis is effectivelylogarithmic. of × M (cid:12) within one kpc, so another cluster could in-spiral. Most dwarf galaxies have likely dissolved in the oldstellar halo of our Galaxy. At z = , there were roughly threetimes as many dwarf satellites as today (Kazantzidis et al.2008). So, there is some chance that several GCs mergedwithin dwarf galaxies in the past. Quantifying the rate ofsuch mergers is beyond the aims of this paper.Finally, it is possible that two GCs merge slightly after MNRAS , 1–13 (2016) erging to explain multiple populations and rotation in GCs † ∇ V · R / σ R hm R hm /2R hm /4 Fabricius+14M22 (Bellazzini+12)M2 (Pryor+86) ω Cen (Bianchini+13)M54 (Bellazzini+12)NGC1851 (Lardo+15)
Figure 9. ( ∇ V · R / σ , ε ) for each model at different values of R .The observational points (grey) are indicated in legend in the top-left corner. The black dashed line indicates the behaviour of anisotropic oblate rotator. their formation, when they are still part of the same pro-genitor molecular cloud. In this case, their relative velocityshould be of the same order of magnitude as the turbulentmotions inside the cloud (approximately few km s − ), en-abling a successful merger. There are clusters younger than100 Myr that are believed to be “caught in the act” of merg-ing while they are still within the parent cloud (Sabbi et al.2012).In summary, we confirm that the merger scenario mayprovide a viable explanation for multiple populations in iron-complex GCs. Our simulations show that the relative con-centration in the merger remnant betrays the initial densityratio of the progenitors. Moreover, the density ratio of theprogenitors leaves a signature in the rotation curves thatshould be object of further observations. ACKNOWLEDGEMENTS
We thank the anonymous referee for their critical readingof the manuscript and for their comments, which helpedus improving this work significantly. We thank Kim Venn,Giacomo Beccari, Eugenio Carretta and Raffaele Grattonfor useful discussions. The simulations were performed withthe Tasna GPUs cluster of ZBOX4 at University Zurich,PLX and Eurora clusters at CINECA (through CINECAAward N. HP10B338N6 and HP10CZVZHA ). We acknowl-edge the CINECA Award N. HP10B338N6, HP10CZVZHA † V m a x / σ M4 ρ [all] M2 ρ [all]M1 ρ [all] This workFabricius+14M22 (Bellazzini+12)M2 (Pryor+86) ω Cen (Bianchini+13)M54 (Bellazzini+12)NGC1851 (Lardo+15)
Figure 10. ( V / σ , ε ) diagram. The purple stars refer to this work,using the maximum rotation velocity ( V max ) for V . The observa-tional points (grey) are indicated in the legend. The dashed blackline shows the behaviour of isotropic oblate rotators. -1 R hl (pc)10 L ( L fl ) Figure 11.
Total luminosity versus half-light radius of MilkyWay GCs from Harris (1996) . and the University of Zurich for the availability of highperformance computing resources. EG acknowledges finan-cial support through SNF grant and Foundation MERAC2014 Travel grant. MM acknowledges financial support fromthe Italian Ministry of Education, University and Research(MIUR) through grant FIRB 2012 RBFR12PM1F, from
MNRAS000
MNRAS000 , 1–13 (2016) E. Gavagnin et al.
INAF through grant PRIN-2014-14 (Star formation andevolution in galactic nuclei) and from Foundation MERACthrough grant ‘The physics of gas and protoplanetary discsin the Galactic centre’.
REFERENCES
Alladin S. M., 1965, ApJ, 141, 768Amaro-Seoane P., Konstantinidis S., Brem P., Catelan M., 2013,MNRAS, 435, 809Anderson A. J., 1997, PhD thesis, UNIVERSITY OF CALIFOR-NIA, BERKELEYAnderson J., King I. R., 2003, AJ, 126, 772Anderson J., Piotto G., King I. R., Bedin L. R., GuhathakurtaP., 2009, ApJ, 697, L58Bastian N., Lamers H. J. G. L. M., de Mink S. E., LongmoreS. N., Goodwin S. P., Gieles M., 2013, MNRAS, 436, 2398Bekki K., Yong D., 2012, MNRAS, 419, 2063Bellazzini M., et al., 2008, AJ, 136, 1147Bellazzini M., Bragaglia A., Carretta E., Gratton R. G., LucatelloS., Catanzaro G., Leone F., 2012, A&A, 538, A18Bellini A., Piotto G., Bedin L. R., King I. R., Anderson J., MiloneA. P., Momany Y., 2009, A&A, 507, 1393Bellini A., Bedin L. R., Piotto G., Milone A. P., Marino A. F.,Villanova S., 2010, AJ, 140, 631Bertin G., Varri A. L., 2008, ApJ, 689, 1005Bianchini P., Varri A. L., Bertin G., Zocchi A., 2013, ApJ, 772,67Bianchini P., Renaud F., Gieles M., Varri A. L., 2015, MNRAS,447, L40Binney J., 1978, MNRAS, 183, 501Binney J., Tremaine S., 2008, Galactic Dynamics: Second Edition.Princeton University PressCanon R.-D., Stobie R.-S., 1973, MNRAS, 162, 207Carraro G., Lia C., 2000, A&A, 357, 977Carretta E., 2015, preprint, ( arXiv:1510.00507 )Carretta E., et al., 2007, A&A, 464, 967Carretta E., et al., 2009a, A&A, 505, 117Carretta E., Bragaglia A., Gratton R., Lucatello S., 2009b, A&A,505, 139Carretta E., Bragaglia A., Gratton R. G., Recio-Blanco A., Lu-catello S., D’Orazi V., Cassisi S., 2010a, A&A, 516, A55Carretta E., et al., 2010b, A&A, 520, A95Carretta E., et al., 2010c, ApJ, 722, L1Carretta E., Lucatello S., Gratton R. G., Bragaglia A., D’OraziV., 2011, A&A, 533, A69Chandrasekhar S., 1969, Ellipsoidal figures of equilibriumChandrasekhar S., Lebovitz N. R., 1962, ApJ, 136, 1037Chen C. W., Chen W. P., 2010, ApJ, 721, 1790Chernoff D. F., Weinberg M. D., 1990, ApJ, 351, 121Cohen J. G., 1978, ApJ, 223, 487D’Ercole A., Vesperini E., D’Antona F., McMillan S. L. W., Rec-chi S., 2008, MNRAS, 391, 825D’Orazi V., Gratton R. G., Pancino E., Bragaglia A., CarrettaE., Lucatello S., Sneden C., 2011, A&A, 534, A29Da Costa G. S., Held E. V., Saviane I., 2014, MNRAS, 438, 3507De Marchi G., Panagia N., Guarcello M. G., Bonito R., 2013,MNRAS, 435, 3058Decressin T., Meynet G., Charbonnel C., Prantzos N., Ekstr¨omS., 2007, A&A, 464, 1029Emsellem E., et al., 2007, MNRAS, 379, 401Fabricius M. H., et al., 2014, ApJ, 787, L26Fall S. M., Frenk C. S., 1985, in Goodman J., Hut P., eds, IAUSymposium Vol. 113, Dynamics of Star Clusters. pp 285–296Ferraro F. R., et al., 2009, Nature, 462, 483Freeman K.-C., Rodgers A.-W., 1975, ApJ, 201, L71 Gaburov E., Harfst S., Portegies Zwart S., 2009, New Astron., 14,630Geyer E. H., 1967, Z. Astrophys., 66, 16Gieles M., 2013, in Pugliese G., de Koter A., Wijburg M., eds,Astronomical Society of the Pacific Conference Series Vol. 470,370 Years of Astronomy in Utrecht. p. 339 ( arXiv:1209.2071 )Gratton R. G., et al., 2001, A&A, 369, 87Gratton R., Sneden C., Carretta E., 2004, ARA&A, 42, 385Gratton R. G., Carretta E., Bragaglia A., 2012, A&ARv, 20, 50Harris W. E., 1974, ApJ, 192, L161Harris W. E., 1996, AJ, 112, 1487Hesser J. E., Bell R. A., 1980, ApJ, 238, L149Hesser J. E., Hartwick F. D. A., McClure R. D., 1977, ApJS, 33,471Johnson C. I., Rich R. M., Pilachowski C. A., Caldwell N., MateoM., Bailey III J. I., Crane J. D., 2015, AJ, 150, 63Kayser A., Hilker M., Grebel E. K., Willemsen P. G., 2008, A&A,486, 437Kazantzidis S., Bullock J. S., Zentner A. R., Kravtsov A. V.,Moustakas L. A., 2008, ApJ, 688, 254Kimmig B., Seth A., Ivans I. I., Strader J., Caldwell N., AndertonT., Gregersen D., 2015, AJ, 149, 53King I. R., 1966, AJ, 71, 64Kraft R. P., Sneden C., Langer G. E., Prosser C. F., 1992, AJ,104, 645Kruijssen J. M. D., 2014, Classical and Quantum Gravity, 31,244006Lardo C., et al., 2013, MNRAS, 433, 1941Lardo C., et al., 2015, A&A, 573, A115Larsen S. S., Brodie J. P., Strader J., 2012, A&A, 546, A53Law D. R., Majewski S. R., 2010, ApJ, 718, 1128Lee J.-W., 2015, ApJS, 219, 7Lee Y.-W., Joo J.-M., Sohn Y.-J., Rey S.-C., Lee H.-C., WalkerA. R., 1999, Nature, 402, 55Lightman A. P., Shapiro S. L., 1978, Reviews of Modern Physics,50, 437Lynden-Bell D., 1967, MNRAS, 136, 101Majewski S. R., Patterson R. J., Dinescu D. I., Johnson W. Y.,Ostheimer J. C., Kunkel W. E., Palma C., 2000, in Noels A.,Magain P., Caro D., Jehin E., Parmentier G., Thoul A. A.,eds, Liege International Astrophysical Colloquia Vol. 35, LiegeInternational Astrophysical Colloquia. p. 619 ( arXiv:astro-ph/9910278 )Makino J., Akiyama K., Sugimoto D., 1991, Ap&SS, 185, 63Mapelli M., Bressan A., 2013, MNRAS, 430, 3120Mapelli M., Zampieri L., Ripamonti E., Bressan A., 2013, MN-RAS, 429, 2298Marino A. F., Milone A. P., Piotto G., Villanova S., Bedin L. R.,Bellini A., Renzini A., 2009, A&A, 505, 1099Marino A. F., et al., 2015, MNRAS, 450, 815Massari D., et al., 2014, ApJ, 795, 22Mastrobuono-Battisti A., Perets H. B., 2013, ApJ, 779, 85Meylan G., Mayor M., 1986, A&A, 166, 122Milone A. P., Piotto G., Bedin L. R., Sarajedini A., 2008, Mem.Soc. Astron. Italiana, 79, 623Milone A. P., Stetson P. B., Piotto G., Bedin L. R., Anderson J.,Cassisi S., Salaris M., 2009, A&A, 503, 755Milone A. P., et al., 2010, ApJ, 709, 1183Milone A. P., Marino A. F., Piotto G., Bedin L. R., Anderson J.,Aparicio A., Cassisi S., Rich R. M., 2012, ApJ, 745, 27Milone A. P., et al., 2013, ApJ, 767, 120Milone A. P., et al., 2015, MNRAS, 447, 927Mucciarelli A., Lapenna E., Massari D., Pancino E., Stetson P. B.,Ferraro F. R., Lanzoni B., Lardo C., 2015, ApJ, 809, 128Nataf D. M., Gould A. P., Pinsonneault M. H., Udalski A., 2013,ApJ, 766, 77Norris J. E., Da Costa G. S., 1995, ApJ, 447, 680MNRAS , 1–13 (2016) erging to explain multiple populations and rotation in GCs Pancino E., Ferraro F. R., Bellazzini M., Piotto G., Zoccali M.,2000, ApJ, 534, L83Pancino E., Galfo A., Ferraro F. R., Bellazzini M., 2007, ApJ,661, L155Pancino E., Rejkuba M., Zoccali M., Carrera R., 2010, A&A, 524,A44Pancino E., Mucciarelli A., Sbordone L., Bellazzini M., PasquiniL., Monaco L., Ferraro F. R., 2011, A&A, 527, A18Pasquato M., Chung C., 2016, preprint, ( arXiv:1602.00993 )Peterson R. C., 1980, ApJ, 237, L87Piotto G., et al., 2012, ApJ, 760, 39Portegies Zwart S. F., McMillan S. L. W., Hut P., Makino J.,2001, MNRAS, 321, 199Pryor C., Hartwick F. D. A., McClure R. D., Fletcher J. M.,Kormendy J., 1986, AJ, 91, 546Read J. I., Goerdt T., Moore B., Pontzen A. P., Stadel J., LakeG., 2006, MNRAS, 373, 1451Renzini A., 2008, MNRAS, 391, 354Roberts P. H., 1962, ApJ, 136, 1108Sabbi E., et al., 2012, ApJ, 754, L37Sarajedini A., Layden A. C., 1995, AJ, 109, 1086Saviane I., da Costa G. S., Held E. V., Sommariva V., GullieuszikM., Barbuy B., Ortolani S., 2012, A&A, 540, A27Sippel A. C., Hurley J. R., Madrid J. P., Harris W. E., 2012,MNRAS, 427, 167Sneden C., Kraft R. P., Prosser C. F., Langer G. E., 1991, AJ,102, 2001Spitzer L., 1987, Dynamical evolution of globular clustersSugimoto D., Makino J., 1989, PASJ, 41, 1117Trani A. A., Mapelli M., Bressan A., 2014, MNRAS, 445, 1967Varri A. L., Bertin G., 2012, A&A, 540, A94Vesperini E., Varri A. L., McMillan S. L. W., Zepf S. E., 2014,MNRAS, 443, L79Villanova S., Geisler D., Gratton R. G., Cassisi S., 2014, ApJ,791, 107White S. D. M., 1978, MNRAS, 184, 185White R. E., Shawl S. J., 1987, ApJ, 317, 246Yong D., Grundahl F., 2008, ApJ, 672, L29van den Bergh S., 1996, ApJ, 471, L31van den Bergh S., 2008, AJ, 135, 1731van den Bosch R., de Zeeuw T., Gebhardt K., Noyola E., van deVen G., 2006, ApJ, 641, 852This paper has been typeset from a TEX/L A TEX file prepared bythe author.MNRAS000