# Properties of globular clusters formed in dark matter mini-halos

Eduardo Vitral, Pierre Boldrini

► To cite this version:

Eduardo Vitral, Pierre Boldrini. Properties of globular clusters formed in dark matter mini-halos. Astronomy and Astrophysics - A&A, 2022, 667, pp.A112. 10.1051/0004-6361/202244530. hal-03482012

HAL Id: hal-03482012

https://hal.science/hal-03482012

Submitted on 20 Mar 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Properties of globular clusters formed in dark matter mini-halos*

Eduardo Vitral<fnref n="1" />● and Pierre Boldrini<fnref n="1" />,<fnref n="2" />

¹ Sorbonne Université, CNRS, UMR 7095, Institut d'Astrophysique de Paris, 98 bis Bd Arago, 75014 Paris, France
e-mail: vitral@iap.fr

² Université de Lorraine, CNRS, Inria, LORIA, 54000 Nancy, France

Received 18 July 2022 / Accepted 18 August 2022

## ABSTRACT

We seek to differentiate dynamical and morphological attributes between globular clusters (GCs) that were formed inside their own dark matter (DM) mini-halo from those who were not. We employed high-resolution full *N*-body simulations on a Graphics Processing Unit (GPU) of the GCs with and without a DM mini-halo, orbiting a Fornax-like dwarf galaxy. For GCs with DM, we observed that this dark extra mass triggers a tidal radius growth that allows the mini-halo to act as a protective shield against tidal stripping, being itself stripped beforehand. We demonstrate that this shielding effect becomes negligible when the tidal radius is smaller than the half-mass radius of the mini-halo. Contrary to previous predictions, we found that the inflation of outer stellar velocity dispersion profiles is expected for GCs with and without a mini-halo, as a result of the host's tidal field. Moreover, we observed that GCs with a DM mini-halo should have, in general, relatively more radial outer velocity anisotropy profiles throughout all of their orbits, smaller degrees of internal rotation, and as a consequence of the latter, smaller ellipticities for their stellar distribution. Due to dynamical friction, we observed a clear bimodal evolutionary distribution of GCs with and without DM in the integrals of motion space and show that for GCs originally embedded in DM, this method is not reliable for association with previous accretion events. Finally, we provide parametric mass profiles of disrupted DM mini-halos from GCs that are to be used in Jeans modelling and orbital integration studies.

**Key words.** globular clusters: general – dark matter – galaxies: kinematics and dynamics – galaxies: formation – galaxies: dwarf – methods: numerical

## 1. Introduction

There is no doubt that the nature of dark matter (DM) is one of the most elusive concepts in modern-day physics. However, the existence of this astrophysical component has been used and requested to explain a vast range of phenomena for a considerable amount of time, including the following: Back when Zwicky (1933, 1937) proposed that some sort of non-luminous matter could compose the amount of mass needed to explain the discrepancy between mass measurements of the Coma cluster based on the Virial theorem (e.g. Binney & Tremaine 2011) and based on the brightness and number of galaxies, up to recent measurements from the Planck Collaboration XIII (2016), Planck Collaboration VI (2020) which yield very robust fits of the cosmic microwave background (CMB) using a Λ-cold dark matter (ΛCDM) model accounting for the existence of DM. For the period in between, other important confirmations of this mysterious dark component were provided: Rubin & Ford (1970) and Rubin et al. (1980) showed that the rotation curves of outer stars in nearby galaxies needed an extra amount of mass (compared to observable, luminous matter) to explain their high velocity values and finally, gravitational lensing studies (e.g. Taylor et al. 1998) have also confirmed that the total amount of mass in many galaxy clusters corresponds to the dynamical measurements accounting for DM.

Such findings point to this DM as a fundamental component in galaxy formation, which is present in most galaxies as an enveloping halo, from the smallest to the highest

scales. Thus, in principle, it seems curious that dense collections of stars such as globular star clusters (GCs), that span from ~$10^5$ to $10^7 M_\odot$, with some of them thought to be accreted dwarf galaxies (e.g. Majewski et al. 2000; Bekki & Freeman 2003; Pechetti et al. 2022), do not seem to require any significant amount of DM to explain their dynamical mass (e.g. Shin et al. 2013; Conroy et al. 2011; Ibata et al. 2013; Moore 1996; Baumgardt et al. 2009; Lane et al. 2010; Feng et al. 2012; Hurst et al. 2015), although some recent works seem to suggest otherwise for a few particular sources (Carlberg & Grillmair 2022; Errani et al. 2022). In fact, Peebles (1984) proposed a formation scenario where GCs are formed inside their own DM mini-halo<fnref n="1" />, and further studies defended that if formed before the end of re-ionisation<fnref n="2" />, GCs could be smaller counterparts of galaxies (e.g. Bromm & Clarke 2002, Fig. 2 from Mamon et al. 2012; Silk & Mamon 2012 for a review on galaxy formation).

On the other hand, different formation scenarios, where GCs are not necessarily embedded in DM mini-halos also exist. For instance, GCs could be formed as bound gas clouds (Peebles & Dicke 1968), as galaxy DM-free fragments (e.g. Searle & Zinn 1978; Abadi et al. 2006), as relics of young massive clusters (YMCs, Portegies Zwart et al. 2010; Longmore et al. 2014) that formed in the high-redshift Universe (Kruijssen 2014, 2015), or simply that formed in situ

<footnote n="1">¹ Those DM mini-halos could have between $10^6$ and $10^8 M_\odot$, such as general DM sub-structures or sub-halos (Zavala & Frenk 2019).</footnote>

<footnote n="2">² Such early formation is consistent with GCs having typical ages up to ~13 Gyr (Marín-Franch et al. 2009), as well as with recent findings of high-redshift GCs by the James Webb Space Telescope (Mowla et al. 2022).</footnote>

<footnote>* Movies associated to Fig. 6 are available at
https://www.aanda.org</footnote>
along with its host galaxy (e.g. Harris 1991; Forbes et al. 1997). Moreover, recent cosmological simulations indicate realistic mechanisms through which all these scenarios (i.e. DM and DM-free) can actually happen (Trenti et al. 2015; Kimm et al. 2016; Ricotti et al. 2016; Keller et al. 2020; Ma et al. 2020), making it reasonable to argue that GCs likely originate from more than a single formation channel<fnref n="3" />.

Another reason that makes the understanding and further confirmation of different formation channels difficult is that much of the main consequences of these channels are better observable in the outskirts of GCs, where one often lacks good quality data. For example, the detection of tidal tails or stellar streams that might relate to an accretion event can be complicated by an observational bias where the stars in the stream are less luminous than the ones in the central and dense regions of the cluster, and thus more difficult to observe, comparatively (Balbinot & Gieles 2018). Similarly, a potential DM mini-halo could present a much more diffuse structure than the stellar component, so that its dynamical detectability might only be possible beyond several scale radii (Peñarrubia et al. 2017), where GC stars are usually confused with galactic field stars.

With the astrometric revolution brought by the Gaia mission (Gaia Collaboration 2018a,b, 2021), and the promising future discoveries of the James Webb Space Telescope (JWST, Gardner et al. 2006) and the Euclid mission (Laureijs et al. 2011; Lançon et al. 2021), the need to better constrain the expected differences between the multiple GC formation scenarios is a priority, so that these rich data sets can be fully exploited to better understand the many long-sought questions regarding GC formation. As a matter of fact, although many robust attempts to better model the observational implications of the DM mini-halo scenario have been made, the high computational cost of simulating a GC+DM system in a Milky Way (MW) type of galaxy forced these attempts to be placed in idealised scenarios: for instance, isolated GCs not experiencing tidal forces (Peñarrubia et al. 2017), or orbiting GCs in a static potential (Mashchenko & Sills 2005a). Preferably, simulations with clusters experiencing tidal forces, in a host galaxy fully composed of particles (and not just a static potential) would allow to take into account more correctly dynamical friction, tidal effects between globular clusters and their host galaxy, and dynamical shocks with larger structures of the host.

In this work, we aim to clearly separate the observational behaviours of GCs that are formed inside a DM mini-halo and which are devoid of one, both of them orbiting a host galaxy and thus permanently experiencing a tidal field. We do it by performing N-body simulations of a GC system with and without a DM embedding mini-halo, alongside with a host galaxy. In order to bypass the high computational cost mentioned above while still keeping a high resolution, we place our GCs in a dwarf spheroidal (dSph) galaxy, following the same prescriptions as the Fornax dSph, in a similar manner than done in Boldrini et al. (2020a)<fnref n="4" />. Such setup allows one to consider much less stars than it would be needed if the satellites were orbiting a MW-like galaxy and keep a satisfactory resolution, and also to avoid using a potential instead of particles. Besides, by taking a real

<table><caption>**Table 1.** Initial conditions of our simulations.</caption><thead><tr><th>Object</th><th><i>M</i><br>[10<sup>6</sup> <i>M</i><sub>☉</sub>]</th><th><i>N</i><br>[10<sup>4</sup>]</th></tr></thead><tbody><tr><td>Fornax dwarf</td><td></td><td></td></tr><tr><td>DM</td><td>1000</td><td>2000</td></tr><tr><td>Stars</td><td>38.2</td><td>76.4</td></tr><tr><td>Globular clusters</td><td></td><td></td></tr><tr><td>DM</td><td>20</td><td>40</td></tr><tr><td>Stars</td><td>1</td><td>2</td></tr><tr><td></td><td><i>r</i><br>[kpc]</td><td><i>v</i><sub>x</sub><br>[km s<sup>-1</sup>]</td><td><i>v</i><sub>y</sub><br>[km s<sup>-1</sup>]</td><td><i>v</i><sub>z</sub><br>[km s<sup>-1</sup>]</td></tr><tr><td>GC 1</td><td>5.32</td><td>13.90</td><td>1.31</td><td>14.38</td></tr><tr><td>GC 2</td><td>2.07</td><td>9.43</td><td>19.05</td><td>21.42</td></tr><tr><td>GC 3</td><td>1.95</td><td>37.49</td><td>4.55</td><td>6.21</td></tr><tr><td>GC 4</td><td>2.19</td><td>3.51</td><td>3.87</td><td>34.00</td></tr><tr><td>GC 5</td><td>2.05</td><td>19.62</td><td>29.80</td><td>15.10</td></tr></tbody></table>

**Notes.** From left to right, the upper columns give, for each component (dark matter and stars), the mass (in $10^6 M_☉$) and number of particles (times $10^3$). We also provide initial positions and velocities of our globular clusters (with respect to the Fornax dwarf) in the lower columns.

galaxy (i.e. Fornax) as our model, we are likely better exploring the dynamics and orbital evolution of different parameters, and thus reaching more realistic conclusions.

We describe our methods in Sect. 2, while our results and main DM signatures are presented in Sect. 3. We discuss the implications of our work and conclude it in Sects. 4 and 5. Throughout the rest of the paper, we reference the stellar component of the system of GC stars without DM as ⋆, and the stellar component of the system of GCs formed inside DM mini-halos as ⋆•. The dark matter is labelled as •.

## 2. Methods

### 2.1. Initial conditions

The initial conditions for the Fornax system were taken from Boldrini et al. (2020a), and although a complete description can be seen in their Sect. 2, we remind readers here of the most salient points, which are summarised in Table 1. We consider the scenario from that paper where the GCs were accreted recently, at $z = 0.36$ (i.e. 4 Gyr ago) by Fornax. It ensures that at $z = 0$, the GCs embedded in DM are still orbiting and no star clusters form in the centre of Fornax (in accordance with observations), and also that at 3 Gyr, the cluster's positions relative to Fornax are consistent with their observed projected distances.

The Fornax dSph is modelled with a stellar component following the Plummer (Plummer 1911) profile, with total mass of $3.82 \times 10^7 M_☉$ (de Boer & Fraser 2016) and Plummer scale radius of 1.038 kpc (Strigari et al. 2010). Its DM halo follows a NFW (Navarro et al. 1996) profile with total mass $10^9 M_☉$. Given the halo mass and redshift, the scale radius of the Fornax DM halo was estimated from cosmological N-body simulations (Prada et al. 2012).

The GCs on the other hand, started with a King (1962) stellar distribution, with total mass of $10^6 M_☉$ and a King radius of 1 pc, lower than the observed radius, since it is susceptible to increase through dynamical processes such as mass loss (Mackey & Gilmore 2003). When assigning their DM mini-halo, we prescribed a mass of $2 \times 10^7 M_☉$, also following a NFW profile, and computed the respective scale radius in a similar

<footnote n="3">³ Observational evidence for different GC populations has also been provided, for instance, by noticing significant colour gradients on the GC population from different galaxies (e.g. Cohen et al. 1998; Harris 2009; Wu et al. 2022).</footnote>

<footnote n="4">⁴ Indeed, Boldrini et al. (2020a) have shown that if the GCs from the Fornax dSph were formed in DM mini-halos, this could account for a natural explanation of the cusp-core and timing problems in this dwarf, adding robustness to this scenario.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Tidal radius growth due to an extra dark mass: Scheme illustrating the increase of tidal radius (blue) in the globular cluster system (orange) embedded in a dark matter mini-halo (grey), in the *left*, in comparison with the globular cluster system without dark matter, in the *right*, at a same orbital radius. This tidal radial expansion is mostly explained by its dependence on the satellite's mass, which is greater for the case with a DM mini-halo.</figcaption></figure>

fashion than for the dSph. The GC initial positions in phase space were set after analysing 7395 isolated dwarf galaxies with mass ~10⁹ $M_\odot$ in the Illustris TNG-100 simulation (Pillepich et al. 2018), whose mass resolution is $7.5 \times 10^6 M_\odot$, hence similar to the DM mini-halo mass. The positions were then selected by picking the cases where the sub-halos projected distances to the dwarf, at $z=0$, were in the same range as the ones observed for the Fornax dSph GCs. From this subset, a new set of positions was obtained by going back 4 Gyr in time, and finally, the final values were taken by computing the maximum weights in each dimension.

Our choice of ratio for the mass of the GC mini-halo and stellar component (i.e. $M_{DM}/M_{stars} = 20$) is consistent with recent works targeting this scenario. For instance, Mashchenko & Sills (2005a) simulated GCs with a stellar mass of $8.8 \times 10^4 M_\odot$ and a mini-halo virial mass of $10^7 M_\odot$, while Peñarrubia et al. (2017) tested models with mass ratios $M_{DM}/M_{stars} = \{0, 10, 100\}$. Indeed, GC systems with a DM mini-halo are very similar to orbiting (or accreted) DM sub-halos, which sit in a specific mass range for dwarf galaxies, predicted by the extended Press-Schechter (EPS) formalism (Bond et al. 1991; Lacey & Cole 1993). A Fornax-like galaxy is considered to have accreted roughly twelve DM sub-halos (see Fig. 1 from Boldrini et al. 2020b) with a mass ratio $M_{Fornax}/M_{sub-halo} \sim 10-10^2$ (Zavala & Frenk 2019), consistent with our values.

## 2.2. Simulations

To generate our N-body objects, we use the initial condition code MAGI (Miki & Umemura 2018). Adopting a distribution-function-based method, it ensures that the final realisation of the galaxy is in dynamical equilibrium (Miki & Umemura 2018). We perform our simulations with the high performance collisionless N-body code GOTHIC (Miki & Umemura 2017). This gravitational octree code runs entirely on a Graphics Processing Unit (GPU) and is accelerated by the use of hierarchical time steps in which a group of particles has the same time step (Miki & Umemura 2017). In order for the simulated GCs to fully relax prior to introducing the Fornax tidal field, they were previously left to evolve in isolation for 2 Gyr.

We evolve the Fornax-GC system over 4 Gyr in each scenario, i.e. for GCs with and without a DM mini-halo. We set

the particle resolution of all the live objects to $50 M_\odot$<fnref n="5" /> and the gravitational softening length to 0.1 pc. Numerical convergence tests have been performed in Boldrini et al. (2020a). The use of a host composed of live particles rather than just a static potential provides a better handling of dynamical friction, which is the deceleration of massive particles (here, the GC) after interacting with less massive particles (here, the stellar and DM particles from the host). This is necessary to interpret the evolution of orbital radii of GCs, both with and without a DM mini-halo, and to robustly quantify their differences. Furthermore, the pruning of tidal structures (Binney & Tremaine 2011, Sect. 7.3) due to the shocking of such structures with the host field tracers is an effect completely neglected when using static potentials. Hence, the study of stellar stream's properties and projected ellipticity are also better addressed when using live particles, such as we do.

## 2.3. Data analysis

Using the instantaneous orbital radius $r_{orbit}$ and the instantaneous satellite mass $M_{sat}$, we calculate the theoretical tidal radius of each GC (or GC+DM system) at each snapshot from the simulation as derived by Bertin & Varri (2008), and used in many GC-related studies (e.g. Khalaj & Baumgardt 2016; Daniel et al. 2017; Webb et al. 2019). For the system composed of GC stars plus DM particles (i.e. GC+DM), we measure the total tidal radius in our simulations by taking into account the mass of both DM and stars in our calculation. In order to determine the bound particles of the DM and stellar components of GCs, we follow the procedure described in van den Bosch & Ogiya (2018, see their Sect. 2.3). The bound mass is determined by computing iteratively the binding energy of each particle. Initially, the method assumes that each particle is bound. We use the publicly available code galpy<fnref n="6" /> (Bovy 2015) to compute the quantities $E$ and $L_z$, for the integrals of motion space.

<footnote n="5">⁵ We note that Banik & Bovy (2021) state that mass resolutions of 100 $M_\odot$ or less would be needed to avoid spurious numerical density variations in stellar streams, thus well in agreement with our mass resolution.</footnote>

<footnote n="6">⁶ Available at https://github.com/jobovy/galpy</footnote>
<figure><img src="image_3.png" /></figure>

**Fig. 2.** Dark matter shield: Velocity dispersion map of dark matter particles for GC5, projected in the X vs. Y plane and centred in the centre of mass of the globular cluster system. We display the last six pericentres of its orbit, where the tidal effects are stronger. The extension of bound globular cluster stars and bound dark matter particles are highlighted as dotted and dashed green lines, respectively, while the theoretical tidal radius, calculated according to Sect. 2.3, is displayed as a solid green circle. The maps are colour-coded logarithmically from blue (lower dispersion) to red (higher dispersion). The centres of Fornax and of the globular cluster are represented as a thick green cross and a plus sign, respectively. For this cluster, we notice that the empirical tidal radius, well traced by the blue region, remains always larger than the bound stars radii. This argues in favour of the dynamical presence of a dark matter shield, as illustrated in Fig. 1.

We performed most of the data analysis with internal routines of the astrometry modelling code BALROGO<fnref n="7" /> (Vitral 2021). This includes fits of the surface density (and ellipticity), construction of velocity anisotropy and dispersion profiles, and sky projection. Whenever targeting these themes, we mention in detail how we proceeded to perform the analysis.

# 3. Results: Dark matter signatures

In this section, we analyse the main implications of the presence of a DM mini-halo in the overall GC dynamics and morphology, according to our simulations. We remind that our analysis is modelled on the specific case of the Fornax dSph, in order to probe a realistic scenario.

## 3.1. Dark matter shield

In the beginning of our simulations with DM mini-halos, all GCs presented a DM envelope massive and concentrated enough, so their total tidal radius was considerably greater than the case without a DM mini-halo (see Fig. 1). Whenever we had such tidal radius increase, the DM particles are stripped beforehand GC stars, so the DM envelope works effectively as a shield against tidal stripping of the stellar component, being gradually

removed as the system experiences stronger tidal forces from the host galaxy.

In order to visualise this effect, we created a velocity dispersion map (see Appendix B.1 for details) of DM particles, which are much more spread than GC stars, and thus provide better spatial completeness to the map. Such map helps us to spot the transitory region where the DM particles start to effectively feel tidal effects of the host galaxy, which is characterised by a steep increase of the velocity dispersion produced by tidal heating of the system.

In Fig. 2, we show this velocity dispersion map for GC5, one of the clusters where this DM shield seemed most effective. We display the last six pericentres of its orbit, when the tidal effects are stronger. The extension of bound GC stars and bound DM particles are highlighted as dotted and dashed green lines, respectively, while the theoretical tidal radius, calculated according to Sect. 2.3, is displayed as a solid green circle. The centres of Fornax and of the GC are represented as a thick green cross and a plus sign, respectively.

We also observe, still in Fig. 2, a gradual decrease in the extension of the blue region, characterised by a low velocity dispersion, where the DM shield is effective. This means that at first, the DM shield is highly protective, and with time, as DM particles are stripped from the cluster, the system mass decreases, and the shield becomes weaker. The red colours, on the other hand, point to a region of high dynamical heating, which becomes more intense in the centre of Fornax and favours tidal stripping of the DM shield.

<footnote n="7">⁷ Code repository: https://gitlab.com/eduardo-vitral/balrogo</footnote>
<figure><img src="image_4.png" /><figcaption>Fig. 3. Velocity anisotropy: Evolution of the orbital radius (i.e. $r_{orbit}$, upper plots) and of the outer velocity anisotropy (i.e. $\beta_{\infty}$, lower plots), for the five simulated globular clusters. The clusters originally embedded in a dark matter mini-halo are depicted in solid green, and those without it are displayed in solid red. For $\beta_{\infty}$, we display the running mean over the five closest snapshots (time-wise), for better visualisation, which explains why the initial $\beta_{\infty}$ values seem to be different for the two scenarios. We observe typically lower $r_{orbit}$ and more radial (i.e. higher) $\beta_{\infty}$ for clusters with a dark matter mini-halo.</figcaption></figure>

This protective DM shield is a key mechanism to explain most of the DM signatures we observed in GCs embedded in DM. Although this phenomenon is ubiquitous in our GCs formed in DM mini-halos, the survival of such shield depends mostly on the orbital parameters of each cluster, and a specific discussion on this matter is addressed in Sect. 4.1. In the following, we focus on describing the main impacts of such shield with respect to the case where GCs are not embedded in DM mini-halos.

## 3.2. Orbital decay

The exchange of energy between the satellite and the field stars from its host galaxy will lead to the drag of the satellite, also referred to as dynamical friction. The use of live particles for both the satellite and its host, such as we do, allows to provide a better handling of this effect. The dynamical friction time of the cluster is the timescale needed for the satellite to reach the centre of mass of the host galaxy. It has been defined in Binney & Tremaine (2011, Eqs. (7)–(26)), and follows the relation $t_{fric} \propto 1/M$, where M the satellite’s total mass.

Applying this relation to our simulations, where in one case we have a GC system alone of $10^6 M_{\odot}$ and in the other case where the total system mass is that of the GC plus the DM mini-halo (i.e. $2 \times 10^7 M_{\odot}$), we find that the system with DM is supposed to sink to the centre twenty times faster than the system without the DM mini-halo. Indeed, when looking at the upper panels of Fig. 3, one can notice that the systems with DM (solid green) occupy much shorter orbital radii than the systems without DM (dashed red) throughout their orbits.

As a consequence, systems with a DM mini-halo tend to be located closer to their host centre sooner, and thus feel a stronger dynamical heating imposed by the host galaxy. This dynamical heating, on its turn, can potentially work to remove more GC stars (along with their DM envelope) than in the case without

DM. Hence, our simulations help to answer whether the presence of a DM mini-halo is rather protective or disruptive with respect to the cluster’s stellar component, a point that we discuss in detail in Sect. 4.1.

## 3.3. Velocity anisotropy

One of the points not addressed in previous works targeting a DM mini-halo in GCs (e.g. Mashchenko & Sills 2005a; Peñarrubia et al. 2017) is the evolution of the velocity anisotropy of the stellar component, in both cases with and without a DM mini-halo. The velocity anisotropy (‘anisotropy’ for short) is defined as in (Binney 1980):

$$ \beta(r) = 1 - \frac{\sigma_{\theta}^{2}(r) + \sigma_{\phi}^{2}(r)}{2 \sigma_{r}^{2}(r)}, \qquad (1) $$

where $\theta$ and $\phi$ are the tangential components of the coordinate system, while $\sigma_i^2$ stands for the velocity dispersion of the component $i$ of the coordinate system. In spherical symmetry, $\sigma_{\theta} = \sigma_{\phi}$.

We fit the anisotropy for each cluster and snapshot of our simulations in a Bayesian frame detailed in Appendix A. For that, we consider the generalised parametrisation from Osipkov (1979), Merritt (1985), which has been applied in many works modelling the dynamics of spherical stellar systems (e.g. Walker & Peñarrubia 2011; Mamon et al. 2013; Read et al. 2021; Vitral & Mamon 2021; Vitral et al. 2022)

$$ \beta_{gOM}(r) = \beta_0 + (\beta_{\infty} - \beta_0) \frac{r^2}{r^2 + r_{\beta}^2}, \qquad (2) $$

where $r_{\beta}$ is the anisotropy radius, and $(\beta_0, \beta_{\infty})$ are the anisotropy values at the centre of the system, and at infinity, respectively. We keep in mind that the initial $\beta(r)$ in our simulations is not necessarily isotropic as a result of the prior 2 Gyr relaxation run
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Stellar velocity dispersion: Series of plots for GC5. The *upper* plots relate to the simulations devoid of dark matter mini-halos (★), while the *lower* plots indicate the results for the globular clusters formed in such mini-halos (★●). The two columns on the *left* display hand-picked snapshots where the radial velocity dispersion profile (i.e. $\sigma_r(r)$) resembled better to an isolated case without a dark matter mini-halo (left), and with a massive dark matter mini-halo (right), according to Fig. 2 from Peñarrubia et al. (2017). The column on the *right* presents the evolution of the maximum value of the radial velocity dispersion profile (i.e. $\sigma_{r\text{max}}$), for each scenario concerning the dark matter mini-halo, colour-coded according to the distance of the cluster to the centre of the host galaxy (i.e. $r_{\text{orbit}}$), with two vertical dashed lines corresponding to the instants from the two columns on the *left*. These plots argue that the tidal field from the host galaxy tends to have a much greater impact on inflating the velocity dispersion than the presence of a dark matter mini-halo. In fact, such mini-halos help to protect the cluster from tidal effects, rather than contributing to it.</figcaption></figure>

we mention in Sect. 2.2, which can change the anisotropy shape, specially in the cluster outskirts ($\beta_0$ still lies closer to zero).

As we are most interested in the differences pertaining to the external dynamics of the GC systems, we chose to plot in Fig. 3 (lower plots) the evolution<fnref n="8" /> of $\beta_\infty$ for the GCs with a DM mini-halo (solid green) and without it (solid red). This allows us to quantitatively differentiate the impact of the Fornax tidal field on the cluster outer regions. Indeed, the external tidal field has been shown to play a major role on the outer anisotropy of GCs, by removing stars in radial orbits in the cluster's outskirts, which are more easily pruned by tides (e.g. Takahashi 1995; Tiongco et al. 2016a; Zocchi et al. 2016; Bianchini et al. 2017). As a result, one should expect more tangential (or equivalently, less radial) anisotropy profiles in systems that undergo more severe tidal interactions.

The variations of $\beta_\infty$ are more intense for clusters without the DM mini-halo, since they cover a broader range of orbital radii. Furthermore, an important diagnosis is that the outer anisotropy of clusters with the DM mini-halo is, in general, more radial, even though those clusters occupy regions closer to the host galaxy. We thus associate this tidal resilience with the dynamical presence of the DM shield discussed in Sect. 3.1.

We also notice that the initial value of $\beta_\infty$ (not necessarily zero) does not really impact these conclusions. This is because the values of $\beta_\infty$ not only depend considerably on the orbital radius, but their variation speed is almost identical to the respective orbital radius variation. This means that regardless of the initial anisotropy conditions, the values of $\beta_\infty$ adapt rapidly to the

ongoing tidal field and its further interpretation is more straightforward.

Hence, the outer anisotropy is an important parameter that allows to differentiate both scenarios. While clusters with more radial orbits (i.e. $\beta_\infty \ge 0.4$) could belong to both scenarios, clusters with tangential anisotropy, or less radial orbits (i.e. $\beta_\infty \le 0.4$) were mostly obtained in our simulations without a DM mini-halo. This does not necessarily mean that clusters with very tangential orbits did not have a DM mini-halo in the past, but rather that searches for a current mini-halo surrounding them would likely be in vain.

### 3.4. Velocity dispersion profile

Peñarrubia et al. (2017) simulated isolated GCs embedded in DM mini-halos and showed that due to the extra DM mass, an inflation of the radial velocity dispersion profile towards outer radii<fnref n="9" /> is to be expected. In the current work, we are able to test such predictions for a more realistic scenario, where the GC+DM subsystem is orbiting a host galaxy, and therefore experiencing tidal effects. For that, we measured the radial velocity dispersion profiles (see Appendix B.2 for details) of the bound stars in our clusters in both the cases with and without a DM mini-halo, through the course of their evolution in the Fornax tidal field.

In Fig. 4 (right plots), we show the values of the maximum radial velocity dispersion as a function of time, colour-coded by the distance to the centre of Fornax (i.e. $r_{\text{orbit}}$), in kpc. One of

<footnote n="8">⁸ We display the running mean over the five closest snapshots (time- wise), for better visualisation, which explains why the initial $\beta_\infty$ values seem to be different for the two scenarios.</footnote>

<footnote n="9">⁹ Their models presented an inflated structure at roughly $r \approx 20 r_{1/2}$, where $r_{1/2}$ is the half-mass radius.</footnote>
the most important realisations of our analysis is the fact that the overall shape of the velocity dispersion is much more impacted by the host galaxy’s tidal field than by an eventual DM mini-halo: The amplitude of the dispersion, traced by its value at the peak follows a periodic variation, with same period than the GC orbit, and has values almost uniquely dependent on the ongoing tidal forces.

In the plot, we can clearly observe that at pericentres (blueish), the velocity dispersion inflates as a whole: the tidal heating from the host galaxy is effectively felt more intensively, leading to a higher velocity dispersion. In contrast, at apocentres (reddish or blackish), the cluster is closer from an ideal isolated scenario, and tidal heating is less effective, leading to low velocity dispersion peaks.

To illustrate the much stronger dependence on tidal forces than on possible DM mini-halos, we selected two snapshots for GC5 for the case with and without a DM mini-halo (lower and upper plots, respectively). In these snapshots (Fig. 4, four left plots), we can verify that the radial velocity dispersion profile of both the DM and DM-free scenarios assume forms similar to both the isolated cases with no DM mini-halo (left), and with a massive DM mini-halo (right), as presented in Fig. 2 of Peñarrubia et al. (2017).

As a general trend, all the clusters had an increasing velocity dispersion close to pericentre, with multiple points of velocity dispersion inflation throughout the GC radial extension. In apocentres, as mentioned above, the clusters resembled better to an isolated case (with one or two inflation points), specially for the case with DM mini-halo, where the shapes retrieved by Peñarrubia et al. (2017) could be better spotted. The reason behind the best resemblance in apocentre for the case with DM mini-halos is directly connected to the DM shield, which protects the GC stars and approaches better the ideal isolated framework.

### 3.5. Rotation and ellipticity (flattening)

Satellites that are more affected by tidal effects should exchange more internal angular momentum from energy transfer with their host, and this will directly impact the degrees of internal rotation of the analysed system (e.g. Tiongco et al. 2016b, 2018). Hence, to quantify the differences pertaining to rotation between GCs formed in DM mini-halos and those who were not, we plot in the upper panel of Fig. 5, the evolution of their mean line-of-sight (LOS) velocity dispersion<fnref n="10" /> (with $1-\sigma$ Poisson error bars), since this quantity is used in many observational studies to constrain the degree of rotation of GCs (e.g. Bianchini et al. 2013; Boberg et al. 2017; Lanzoni et al. 2018). In order to recreate a two-dimensional projected distribution of stars, such as we observe in true data, and extract the respective LOS velocity, we placed, for each snapshot, the GC centre in the position $(\alpha, \delta) = (0, 0)$, at five kpc away from the observer<fnref n="11" />. The dispersion in

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Rotation: On the top, we show the mean evolution of the line-of-sight (LOS) velocity dispersion of bound globular cluster particles, for each simulated cluster, while the bottom plot displays the mean ellipticity (or flattening) of the stellar distribution. We display the running mean over the ten closest snapshots (time-wise), for better visualisation. Error bars are represent $1-\sigma$ Poisson uncertainties. Clusters originally embedded in a dark matter mini-halo (i.e. ★●) are displayed in green, while those devoid of the mini-halo (i.e. ★) are represented in red. As a general trend, globular clusters embedded in dark matter mini-halos exchange less internal angular momentum from tidal interactions, thus presenting smaller degrees of rotation, which in turn is associated with smaller ellipticities. This again attests the efficacy of the dark matter shield illustrated in Fig. 2.</figcaption></figure>

LOS velocity of all bound GC stars is neatly higher for clusters without a DM mini-halo, meaning that they suffered more prominent tidal interactions. In other words, GCs embedded in dark matter mini-halos exchange less internal angular momentum from tidal interactions, thus presenting smaller degrees of rotation.

An increased amount of rotation can also contribute to larger values of ellipticities (or flattening) for the stellar distribution from GCs (Fabricius et al. 2014; Kamann et al. 2018). We thus fitted the distribution of bound stars with a Sérsic (Sérsic 1963; Sersic 1968) asymmetric model, using the recipe described in Appendix C. The fits yielded a semi-major (a) and semi-minor axis (b), from which the ellipticity (or flattening) could be calculated as

$$e = \sqrt{1 - \left(\frac{b}{a}\right)^2} \qquad (3)$$

The lower panel from Fig. 5 was then constructed by computing, at each instant, the mean ellipticity<fnref n="12" /> of the five simulated clusters with and without DM. The error bars were calculated as $\sqrt{\langle\epsilon\rangle^2 + \sigma_e^2/n}$, where $\langle\epsilon\rangle$ is the mean uncertainty of the ellipticity fit, $\sigma_e$ is the standard deviation of the values from the $n=5$ clusters.

First, we observe a decreasing pattern in the mean ellipticities of both scenarios, which is expected from orbiting satellites,

<footnote n="10">¹⁰ We display the running mean over the closest ten snapshots, for better visualisation.</footnote>

<footnote n="11">¹¹ Very high distances would make the stars look closer, which can decrease the precision of our fits, and also erase sky-projection signatures likely to be observed in real data. Very short distances, on the other hand, could have unrealistic sky-projection effects. We argue that 5 kpc is a good compromise, as there is a fair amount of galactic GCs with well measured data, at roughly this distance.</footnote>

<footnote n="12">¹² Again, we display the running mean over the closest ten snapshots, for better visualisation.</footnote>
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Dark matter influence on tidal tails: Stellar distribution (and respective zoom in) of GC5 for the case where it was originally embedded in a dark matter mini-halo (right) and where it was devoid of it (left), both centred in the GC centre of mass. The chosen snapshot was such that GC5 was located at roughly the same orbital radius for both scenarios, for a fairer comparison. This plot highlights more prominent tails in the dark matter-free case, while clusters formed inside dark matter mini-halos present a more compact stellar envelope. Hence, the dark matter shield (Sect. 3.1) has the effect of delaying the formation of tidal tails. Movies are available [online](#).</figcaption></figure>

due to the stripping of outer stars by the tidal field (Akiyama 1991)<fnref n="13" />. Nonetheless, one notices that clusters embedded in a DM mini-halo show considerably smaller ellipticities (i.e. *e* ≤ 0.35) than the clusters devoid of DM, which presented rather *e* ≥ 0.35. Such relatively small ellipticities are indeed predictable for systems less affected by tidal forces (van den Bergh 2008), attesting the efficacy of the DM shield illustrated in Fig. 2.

### 3.6. Tidal tails

The study of tidal tails in GCs has been revolutionised by Gaia data and simulations. For instance, a troubling question that arises when simulating GCs on tidal fields is to understand why these simulations usually predict much more prominent tidal tails (e.g. Boldrini & Vitral 2021; Montuori et al. 2007) than what is observed for the majority of MW clusters. The answer to this timing problem is partially given by Balbinot & Gieles (2018), who showed with simulations that there is a preferential bias towards the escaping of low-mass stars, specially in denser clusters. Such trend reduces considerably the visibility of the tails. In addition, Gieles et al. (2021) recently defended this trend by showing that the visible and extended tails of Palomar 5 are well explained by a supra-massive population of stellar-mass black holes, which is a characteristic associated with less dense GCs (Kremer et al. 2020). On the other hand, the presence of a DM mini-halo could also reduce the prominence of tidal tails, since the mini-halo is expected to be stripped beforehand GC stars, thus delaying tail formation (e.g. Bromm & Clarke 2002;

Mashchenko & Sills 2005a; Saitoh et al. 2006; Bekki & Yong 2012; Boldrini & Vitral 2021).

In any case, the increase of better quality data such as Gaia EDR3 has allowed to go deeper into this question: Although some GCs, such as NGC 1851 and NGC 7089 (M2), were thought not to have tails based on ground-based imaging (Kuzma et al. 2016, 2018), further Gaia studies revealed long tails associated to them (Ibata et al. 2021). However, many clusters with no tidal features, or with only extended envelopes (without tails, see Piatti & Carballo-Bello 2020) are still present<fnref n="14" />, and there is still no consensus on which mechanisms are behind the lack of very extended tails in these GCs. On a similar note, Martin et al. (2022) recently used Gaia EDR3 and CFHT data to constrain a stellar stream, C-19, whose metallicity is consistent with it being a remnant of the oldest GC known in our Galaxy, even though the stream still presents a coherent structure to this day, thus requiring some sort of shielding mechanism to protect it from tidal stripping.

We decide to explore the consequences of a DM mini-halo on tidal tails, and we illustrate it in Fig. 6, with the stellar distribution of GC5 (one of the cases where the DM shield was most effective), for the case with (right) and without (left) an embedding DM mini-halo. The chosen snapshot was such that GC5 was located at roughly the same orbital radius for both scenarios, for a fairer comparison. We also provide [online](#) a video of the evolution of the tails in the two kind of clusters, with a respective zoomed-in version too. We observe much more prominent and obvious tails in the case without DM, with thick and well-defined streams measuring up to ≥10 kpc long. However, in the

<footnote n="13">¹³ Such an effect is better modelled by the use of live particles for both the host and the satellite, as in our simulations.</footnote>

<footnote n="14">¹⁴ See Table 3 from Zhang et al. (2022).</footnote>
<figure><img src="image_8.png" /></figure>

**Fig. 7.** Integrals of motion: Evolution of our simulated clusters in the integrals of motion (IOM) space, in particular, the $E \times L_z$ diagram. Their positions in this diagram are scattered in a sequential colour-map, starting from the same point in light grey at $t = 0$, and ending in darker tones (green for globular clusters originally formed in dark matter mini-halos, and red for clusters devoid of dark matter). The last snapshot is marked by a cross, for each cluster. This plot highlights a bimodal evolutionary distribution of clusters originally embedded in dark matter (moving towards lower energies) and those who were not (moving only slightly towards higher energies), and more importantly, it shows that clusters originally embedded in dark matter move significantly in IOM space, such that their association with past merger events through this diagram is not reliable.

case where the GCs are formed inside DM mini-halos, the stellar distribution remains roughly spherical. We connect this difference between two kind of GCs to the protection of the DM shield (as explained in Sect. 3.1, and in Figs. 1, 2), which reduces and delays tidal effects on the GCs where it is present.

Hence, although GCs embedded in DM mini-halos can still form stellar streams eventually, the ones devoid of DM seem to develop much longer, thicker and well-defined tails. Evidently, once the DM shield is destroyed, which can happen much before a Hubble time for clusters having specific orbital parameters (see Sect. 4.1), the development of tidal tails can be considered similar to the case without DM, and thus even clusters originally embedded in DM could present extended tails by present time, in those conditions. Having said that, the DM shield has the effect of delaying the formation of tidal tails, by making them much milder while the shield is present, and thus providing a possible answer to the timing problem mentioned in the beginning of this section.

## 3.7. Integrals of motion

By analysing the orbital evolution of simulated satellites in a static MW potential, Helmi & de Zeeuw (2000) proposed that the space of integrals of motion (IOM) was an optimal tool to recover the orbital history of accreted and disrupted satellites. In particular, by analysing the diagram of energy versus angular momentum in the Z galactic direction, they showed that GCs originating from an accreted satellite should remain in a similar position than its respective progenitor in IOM space, regardless of their phase and velocities being no longer the same. Even though energy is not conserved throughout orbital evolution of the satellite, the authors argued that those IOM positions should

not move significantly from their original positions. Indeed, by using such method, Massari et al. (2019) related the origin of the galactic GC population to different past major accretions events from the MW.

We decided to analyse the impact of a DM mini-halo in the evolution of a GC in IOM space by plotting their respective position in the energy-angular momentum (i.e. $E \times L_z$) diagram for our five analysed clusters, at each snapshot. Figure 7 displays this evolution for both scenarios: with a DM mini-halo in green, and without it in red, both beginning in the same position, marked by a light grey colour, and ending at their respective darker crosses.

One observes that clusters devoid of DM (red) do not significantly move in the diagram throughout their evolution, in agreement with Helmi & de Zeeuw (2000), which in turn serves as a validity check for our simulations. However, the same cannot be said for GCs originally embedded in DM, which move considerably during their orbital evolution. This is mainly due to the massive loss of energy from dynamical friction (see Sect. 3.2), and further changes on the host potential, which are better handled by our use of an evolving host composed of particles, instead of static potentials. Furthermore, this energy loss shifts the GC population formed in DM mini-halos to lower regions of the $E \times L_z$ diagram, creating a clear bimodal evolutionary distribution in Fig. 7, where green evolution tracks move downwards (lower energies), while the red ones tend to move only slightly upwards (higher energies). This does not imply though, that at a given time, a cluster with low energy necessarily had a DM mini-halo, since a DM-free GC could naturally form in such regions and not move significantly thereafter.

Given that there is no reason why these conclusions should change qualitatively when considering the case of the MW, Fig. 7 has important implications. Not only it shines a light on the bimodal evolutionary distribution commented above, which can be used to guide the search of GCs candidates that have been formed in a DM mini-halo, but most of all, it shows that GCs formed in DM mini-halos move significantly in IOM space due to their more intense energy loss from dynamical friction. Hence, the association of such clusters with past merger events by using the IOM method is not reliable. Indeed, Boldrini & Bovy (2022) recently showed, by performing orbital integrations, that there was no clear association between many MW GCs and other MW satellites, even though a large fraction of them are believed to be accreted from previous IOM analyses.

# 4. Discussion

In the previous section, we have consistently showed that, at least during the first orbits of a GC embedded in a DM mini-halo, the DM behaves as a protecting shield that prevents the development of prominent tidal features. Now, we discuss a different aspect, which is under which circumstances this DM shield is effective.

## 4.1. Survival of the dark matter mini-halo

The DM mini-halo of our simulated GCs was shown to effectively work as a shield against tidal effects. However, not only this mini-halo is gradually removed by being itself stripped beforehand GC stars, as the increase of mass in the GC+DM system (with respect to the GCs devoid of DM) also brings the system closer to the centre of the host galaxy<fnref n="15" />, where tidal

<footnote n="15">¹⁵ Through orbital decay from dynamical friction (see Sect. 3.2).</footnote>
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Survival of the dark matter mini-halo: Similarly to Fig. 2, we show the velocity dispersion map of dark matter particles for the five globular clusters we simulated, projected in the X vs. Y plane (with exception of GC4, projected on the Z vs. Y plane for better visualisation of its orbit) and centred in the centre of mass of the globular cluster system. We display the last pericentre passages of their respective orbits, where the tidal effects are stronger. One can notice that while GCs 3, 4 and 5 manage to keep the dark matter mini-halo (and thus the dark matter shield effect depicted in Fig. 2), GCs 1 and 2 do not. The extent of the red region for GCs 1 and 2 is related to their different orbital parameters. For GCs 3, 4 and 5, we again notice that the empirical tidal radius is well traced by the internal blue region.</figcaption></figure>

effects are stronger. As a result, the DM mini-halo is expected to become negligible with time, so that the question of interest becomes: When is a DM mini-halo effectively disrupted?

The answering of this question starts when comparing the five GCs we simulated, and their different orbital parameters. First, we can separate these five GCs into two categories: (i) The GCs where the DM shield is not effective by the end of our simulation and (ii) GCs whose shield still manages to protect them against strong tidal effects. Looking at Fig. 8, we can clearly assign GCs 1 and 2 to category (i), while GCs 3, 4 and 5 are better suited to category (ii), with GC 4 approaching to category (i). This assignment is due to the fact that the limit of bound stars in GCs 1 and 2 is comparable or greater than the respective tidal radius of the GC+DM system, meaning that the DM shield illustrated in Fig. 1 is no longer well observed. We also managed to observe that such disrupted mini-halos presented $r_t > r_{1/2,•}$ at some point of their orbital evolution (Fig. D.1), where $r_t$ is the tidal radius of the GC+DM system, and $r_{1/2,•}$ the half-mass radius of the bound DM particles. This criterion is similar to the relation presented in Hayashi et al. (2003) for the disruption of DM halos.

We argue that the DM loss in GCs 1 and 2 is accelerated by the many passages in pericentre (i.e. relatively small orbital periods), as well as the relatively small orbital radii ($r_{orbit}$) throughout the clusters' orbits. Those factors force these GCs to stay longer next to the centre of the host galaxy, where the tidal interactions are more severe. As a consequence, the DM mini-halos in GCs 1 and 2 suffer from stronger tidal effects and are stripped faster. In Table 2, we display the mean of the structural parameters of stars and DM particles, over the last ten snapshots of the simulation, and by analysing the difference between Cols. 5 and 8 (i.e. the half-mass radius of the DM mini-halo, and the tidal radius of the respective system), one can infer that the DM mini-halo survives better to worse, in the following order: GC5 and GC3 similarly, followed by GC4, GC1 and GC2, respectively.

<table><caption>**Table 2.** Mean of structural parameters from the last ten snapshots, considering only bound dark matter particles and bound stars.</caption><thead><tr><th>ID</th><th>M<sub>•</sub><br>[10<sup>5</sup> M<sub>☉</sub>]</th><th>M<sub>••</sub><br>[10<sup>5</sup> M<sub>☉</sub>]</th><th>M<sub>•</sub><br>[10<sup>5</sup> M<sub>☉</sub>]</th><th>r<sub>1/2,•</sub><br>[pc]</th><th>r<sub>1/2,••</sub><br>[pc]</th><th>r<sub>1/2,•</sub><br>[pc]</th><th>r<sub>LGC+DM</sub><br>[pc]</th><th>r<sub>LGC</sub><br>[pc]</th></tr><tr><th>(1)</th><th>(2)</th><th>(3)</th><th>(4)</th><th>(5)</th><th>(6)</th><th>(7)</th><th>(8)</th><th>(9)</th></tr></thead><tbody><tr><td>GC1</td><td>13.7</td><td>8.7</td><td>6.1</td><td>829</td><td>16.0</td><td>24.0</td><td>510</td><td>600</td></tr><tr><td>GC2</td><td>8.2</td><td>9.1</td><td>8.0</td><td>486</td><td>12.4</td><td>19.2</td><td>179</td><td>360</td></tr><tr><td>GC3</td><td>15.5</td><td>9.6</td><td>8.8</td><td>191</td><td>18.3</td><td>40.9</td><td>816</td><td>1223</td></tr><tr><td>GC4</td><td>10.9</td><td>9.7</td><td>9.0</td><td>354</td><td>12.9</td><td>23.0</td><td>483</td><td>586</td></tr><tr><td>GC5</td><td>9.3</td><td>9.7</td><td>8.6</td><td>181</td><td>17.7</td><td>41.8</td><td>935</td><td>1157</td></tr></tbody></table>

**Notes.** Columns are: (1) Globular cluster ID; (2) Total mass of the surviving dark matter mini-halo (in $10^5 M_☉$); (3) Total mass of the surviving globular cluster originally embedded in a dark matter mini-halo (in $10^5 M_☉$); (4) Total mass of the surviving globular cluster devoid of dark matter (in $10^5 M_☉$); (5) 3D half-number radius of the surviving dark matter mini-halo (in pc); (6) 3D half-number radius of the surviving globular cluster originally embedded in a dark matter mini-halo (in pc); (7) 3D half-number radius of the surviving globular cluster devoid of dark matter (in pc); (8) Tidal radius of the system composed of a globular cluster plus a dark matter mini-halo (in pc); (9) Tidal radius of the globular clustres simulated without a dark matter mini-halo (in pc). Dark matter mini-halos that survived better throughout the simulation have half-mass radii (Col. 5) smaller than the system's tidal radius (Col. 8).

## 4.2. Detectability of dark matter

### 4.2.1. On the amount of detectable dark matter

If some GCs manage to preserve part of their original DM mini-halo, it is important that we understand the observational limitations that may allow us (or not) to detect this DM component. For instance, as the mini-halo mass distribution is more diffuse than the GC stellar component, one might ask if the DM amount in the GC inner regions is significant. For that purpose, we plot
<figure><img src="image_10.png" /><figcaption>Fig. 9. Cumulative mass profiles: We display the cumulative mass profiles at the last snapshot of our simulations, for the dark matter particles (black), their embedded stellar component (green) and the stellar component from the simulations without dark matter (dashed red). The blue line depicts the MCMC fit of a Zhao (1996) αβγ model to the dark matter mini-halo (see Sect. 4.2.2), with a fainter blue region encompassing the [2.5, 97.5] percentiles. The dotted lines represent local slopes of the stellar mass distribution at the respective regions shown in the plot.</figcaption></figure>

in Fig. 9 the mass profile of different tracers (i.e. DM and GCs, for both scenarios) as a function of distance from the cluster's centre, at the last snapshot of our simulations.

Figure 9 reveals that up to ~10 pc<fnref n="16" />, the mass ratio between the DM component (solid black) and the stellar component of the GCs (solid green) is of the order of 1%, meaning that the influence of DM in the internal dynamics of these clusters is mostly negligible<fnref n="17" />. In the past, many studies (Shin et al. 2013; Conroy et al. 2011; Ibata et al. 2013; Moore 1996; Baumgardt et al. 2009; Lane et al. 2010; Feng et al. 2012; Hurst et al. 2015) measured very low mass-to-light ratios in GCs, and sometimes used that as an argument to rule out the presence of dark matter in them. We show, in agreement with what has been previously pointed out in Bromm & Clarke (2002), Mashchenko & Sills (2005a), Saitoh et al. (2006), Bekki & Yong (2012), that GCs originally embedded in DM lose much of their initial DM content due to tidal interactions. This naturally leads to mass ratios (stellar to total mass) very close to unity, up to distances as far as available data usually goes.

Figure 9 also compares the mass profile of the GCs formed in DM mini-halos and those devoid of it (dashed red). Clearly, clusters with DM have an inner mass much higher than GCs devoid of DM, with GCs 3 and 5, which managed to retain most of their initial DM mini-halo, having inner masses up to ten times higher than their DM-free counterparts. Hence, in addition to denser clusters producing an observational bias where tidal structures are more difficultly detected (Balbinot & Gieles 2018), we point that the ones that might be formed in DM mini-halos also tend to have at the same time higher inner densities and less prominent

tidal structures (see Sect. 3.6), which can render the search for stellar streams in these GCs particularly hard.

Finally, it is worth mentioning that the outer slopes of the stellar mass profile (dotted lines in Fig. 9) do not seem to show a necessary correlation with the analysed scenario (with or without DM), differently than what is shown in Peñarrubia et al. (2017, Fig. 2). We associate it once again, as in Sect. 3.4 for the velocity dispersion profile, with the presence of an intense tidal field (neglected in the former work) that shapes the outer regions of the stellar distribution in a more impactful way than an eventual DM mini-halo.

### 4.2.2. Analytical description

One of the main ways to look for a DM mini-halo in GCs is to perform mass modelling of observed data, which can be done by solving the Jeans equation (Binney 1980) for spherical systems with no streaming motions

$$ \frac{d(\rho\sigma_r^2)}{dr} + 2\frac{\beta(r)}{r}\rho(r)\sigma_r^2(r) = -\rho(r)\frac{G M_T(r)}{r^2}, \quad (4) $$

where we assume a total mass profile $M_T(r)$, a anisotropy profile $\beta(r)$, and a previously determined mass density profile $\rho(r)$ for the kinematic tracers (here, stars). The term $\rho\sigma_r^2$ is the dynamical pressure that counteracts gravity. Since the observable tracers are stars, we only need a parametrisation for the DM mass profile in order to model it. This parametric form is then incorporated in the total mass profile $M_T(r) = M_{\star}(r) + M_{\odot}(r)$. To obtain a suitable mass profile for the DM component, we tested two models, namely the Kazantzidis et al. (2004) model<fnref n="18" />,

$$ \rho_{\bullet}(r) = \rho_0 r^{-\gamma} e^{-r/a}, \quad (5) $$

and the Zhao (1996) $\alpha\beta\gamma$ model,

$$ \rho_{\bullet}(r) = \rho_0 r^{-\gamma} \left[ 1 + \left(\frac{r}{a}\right)^{\alpha} \right]^{(\gamma-\beta)/\alpha}, \quad (6) $$

<footnote n="18">¹⁸ This model is motivated from dynamical simulations of repeated tidal encounters (Kazantzidis et al. 2004).</footnote>

<footnote n="16">¹⁶ For our GCs embedded in DM, 10 pc corresponds to roughly 60% of the half mass radii of the stellar component, for the last snapshot of the simulation.</footnote>

<footnote n="17">¹⁷ The plot depicts this behaviour at the end of the simulations, but we stress that throughout the evolution of the cluster, the typical DM cumulative mass remains much lower than the inner GC cumulative mass, given the latter's particularly high values.</footnote>
where $\rho_0$ is a scaling factor, a the scale radius, $\gamma \in (0, 3)$ an inner slope, $\beta$ an outer slope, which we take to be 5, such as in the Plummer profile, and $\alpha$ a scaling exponent, taken as 2, according to the Plummer profile as well<fnref n="19" /> Then, we converted those densities to mass profiles by performing $M_\bullet(r) = \int_0^r 4\pi x^2 \rho_\bullet(x) dx$ with MATHEMATICA¹².

Similarly to what is shown in Mashchenko & Sills (2005a, Fig. 7), we noticed that the stripped DM particles from the original mini-halo had a tendency to form a second shell surrounding the remaining mini-halo. Such effect is only noticeable in the DM mass profile beyond a few hundred pc (see Fig. 9), but since we choose to analyse the mass profiles up to the outer regions of the clusters, we took that into account by assigning $M_\bullet(r) = M_{\bullet,1}(r) + M_{\bullet,2}(r)$, where $M_{\bullet,1}(r)$ is an inner component, most susceptible to be detected or to influence the dynamics of the embedded star cluster, and $M_{\bullet,2}(r)$ is an outer shell, modelled only for the sake of better constraining the outer slope of $M_{\bullet,1}(r)$. This yields, for the Kazantzidis et al. (2004) profile,

$$M_{\bullet,i}(r) = 4\pi a^{3-\gamma} \rho_0 \Gamma\left(3-\gamma, \frac{r}{a}\right) \quad (7)$$

and for the Zhao (1996) $\alpha\beta\gamma$ profile,

$$M_{\bullet,i}(r) = \frac{4\pi a^{3-\gamma} \rho_0}{3-\gamma} \left[ 1 + \left(\frac{a}{r}\right)^2 \right]^{(\gamma-3)/2} \quad (8)$$

where we note $\Gamma(b, x) = \int_0^x t^{b-1} e^{-t} dt$, which is the lower incomplete gamma function, and $i = \{1, 2\}$.

Since we fix the $\alpha\beta$ parameters when using the Zhao (1996) $\alpha\beta\gamma$ profile, we end up with the same amount of free parameters than the Kazantzidis et al. (2004) model, which allows us to compare the fits from the two models by simply choosing the one with highest likelihood as our preferred model. Both models performed well, and the preference of likelihoods depended not only on the simulated cluster, but on the considered radial range. In Fig. 9, we display in blue the fits of the Zhao (1996) $\alpha\beta\gamma$ model, which had an overall better performance for radii smaller than 1 kpc. Indeed, going beyond such distances is not of practical interest when modelling observed data, since at such lengths, the confusion with MW interlopers is considerably damaging for the data set. Hence, for data sets going up to a few tens of pc, Jeans modellers could simply use the Zhao (1996) $\alpha\beta\gamma$ model, with fixed $\alpha = 2$ and $\beta = 5$, and a single shell (i.e. $M_\bullet = M_{\bullet,1}$) to probe the existence of a DM mini-halo in GCs. We provide in Table D.1, the best fit parameters of our simulated clusters (at $z = 0$, for $M_{\bullet,1}$) in each mass model.

### 4.3. Comparison with the literature

We now compare our study with previous attempts to model GCs in DM mini-halos. Beginning with Mashchenko & Sills (2005a), the authors simulated GCs inside a tidal field, using a $M_{DM}/M_{stars}$ ratio of nearly a hundred. Their simulations used static potentials for the host, which tend to underestimate the impacts of dynamical friction and the trimming of outer structures due to passages in high density regions. Their work showed nevertheless that GCs embedded in DM mini-halos are more resilient to severe tidal stripping, in contrast to DM-free GCs, and also that higher mass to light ratios are found in the clusters outskirts, while the inner regions keep mass to light ratios similar to the purely baryonic case. They finally proposed that

presence of obvious stellar streams is probably the only observational evidence that can reliably rule out a presence of significant amounts of DM in GCs.

Our conclusions are in general very similar to those from Mashchenko & Sills (2005a), with some important caveats. Indeed, we find that the presence of a DM mini-halo tends to initially protect the GC from tidal effects, but our more suitable handling of dynamical friction from a time-changing potential with live particles suggests that clusters with a DM mini-halo sink faster to the clusters centre, where the tidal field is more intense. As a result, clusters whose orbital parameters tend to naturally place them closer to the host will have this effect strengthened, resulting in a faster disruption of their mini-halo, and further cluster dissolution. Thus, the question if GCs embedded in DM are always more resilient to tidal effects is not as simple, and some of them might actually be disrupted faster, if originally orbiting near the host's centre (see Sect. 4.1).

As in Mashchenko & Sills (2005a), we also find that the DM mass profile is negligible, with respect to the stellar component, up to at least roughly a half-mass radius from the clusters centre (Fig. 9). However, we do not associate obvious stellar streams as the only indicator of a DM-free GC. In fact, it is important to mention that although systems with DM do not form obvious stellar streams, they can still form milder and thinner ones before the DM mini-halo is entirely disrupted. Another important DM-free indicators could be a very tangential velocity anisotropy in the outskirts and much higher degrees of rotation, which indirectly leads to higher ellipticities as well. These indicators suggest that the GC stellar component is highly impacted by the host's tidal field, hinting that a DM mini-halo is likely not present at current time<fnref n="20" />. We also propose to use IOM space locations to evaluate if the existence of a DM mini-halo is likely (Fig. 7), since that clusters formed in this scenario tend move towards regions of lower energy in the $E \times L_z$ diagram, due to dynamical friction, while still keeping in mind that DM-free clusters can eventually form in these regions and remain there during their evolution. Our new analysis of the IOM evolution additionally shows that associations of GCs and accreted satellites through nearby locations in IOM space are not reliable for clusters that were embedded in a DM mini-halo.

Recently, Peñarrubia et al. (2017) analysed the effects of isolated GCs embedded in DM mini-halos, with DM to stellar mass ratios of ten and a hundred. Their analyses predicted a flattening or increase of the velocity dispersion profile at outer radii for GCs with DM, along with shallower surface density profiles. They also provided analytical forms for the DM mass profile, which can be used in Bayesian approaches. Although such predictions seem indeed valid for isolated systems, when placing them inside a realistic tidal field, such signatures are mostly erased in the cluster outskirts. We have showed in Sect. 3.4 that an eventual inflation in the velocity dispersion profile can also happen in DM-free clusters due to tidal heating, sometimes in even more pronounced proportions (see Fig. 4). When at apocentres though, where tidal interactions are milder, we retrieved similar velocity dispersion shapes than Peñarrubia et al. (2017) for GCs embedded in DM mini-halos, since those are better shielded from tides and resemble better the isolated case. We did not find

<footnote n="19">¹⁹ When fixing those parameters, the goodness of fit was improved.</footnote>

<footnote n="20">²⁰ It is important to notice however, that as pointed in Mashchenko & Sills (2005a), such indicators help to rule out a current DM mini-halo, but not necessarily a DM mini-halo formation scenario, since such mini-halos can be disrupted much before a Hubble time.</footnote>
a clear relation between the outer mass slopes of GC with and
without DM, likely due to the same arguments above.

As in Peñarrubia et al. (2017), we provide an analytical
shape for the mass profiles of DM mini-halos, that can be
used in Bayesian analyses, specially in Jeans modelling rou-
tines (Read et al. 2021). In our case though, our mass profiles
were derived also for partially disrupted mini-halos, such as it is
expected from systems evolving in a tidal field. We thus consider
our parametric forms better suited to the realistic cases of GCs
orbiting the MW or other Local Group galaxies.

Another recent simulation of GCs pertaining to the DM
formation scenario was done by Carlberg & Keating (2022).
However, the authors do not consider a single mini-halo per
cluster such as in our work (or such as in Mashchenko & Sills
2005a; Peñarrubia et al. 2017), but they consider the scenario
where many GCs are formed inside higher mass sub-halos<fnref n="21" />.
Although that analysis is also interesting, a direct comparison
with our work would be unfair, given such very different initial
conditions.

4.4. General caveats

As the simulations we ran required a substantial computer power
and running time, we needed to constrain some of our assump-
tions and exploratory analyses, in order to meet all the numeri-
cal convergence tests performed in Boldrini et al. (2020a). In this
section, we briefly discuss novel aspects that can be envisaged in
future works and how they could affect some of our results, by
changing a few of our initial assumptions.

4.4.1. Mass profiles

Aiming for a more homogeneous treatment and interpretation,
we assume all five GCs and respective mini-halos to have the
same mass distributions. Although such assumption is reason-
able, given the observed stellar masses of the Fornax dSph GCs
(de Boer & Fraser 2016), the more general scenario of Galactic
GCs for example encompasses a much broader mass and scale
radii range (Baumgardt & Hilker 2018), which can deliver more
different results. Another way to look into this variable is that all
of our GCs have different filling factors<fnref n="22" />, and one could even-
tually envisage a numerical simulation where different GCs are
formed with a similar filling factor.

In such scenario, clusters and respective mini-halos that lie
closer to the host would be more compact than in our simula-
tions, which does seem to happen to some extent in the Milky
Way (see GC half-number radii from Vitral 2021 and orbital
radii from Baumgardt et al. 2019, for example). By being more
compact, inner clusters could have their DM shield effective for
longer timescales, and thus present milder or no stellar streams
at present age, as well as smaller degrees of rotation, elliptic-
ity and higher velocity anisotropies, according to our findings.
Similarly, if innermost GCs and respective mini-halos happen
to be less massive, dynamical friction would be less prominent,
and the GC+DM system would sink more slowly to regions
of stronger tides. Again, the shielding consequences mentioned
above would then last for longer in such clusters. We thus pro-
pose, for instance, that future works could explore a greater vari-
ety of GC and mini-halo mass profiles for clusters starting at
similar positions with respect to their host.

4.4.2. Particle masses

In our simulation, all live objects (i.e. stars and DM particles)
are set up with the same mass of 50 $M_\odot$. This numerical recipe
helps to avoid possible heating and segregation effects that might
affect the long-term internal evolution of the GCs embedded in
DM. We thus assume the simplest case of equal particle mass in
order not to bias our results in the case of heavier or lighter DM
particles, and choose to study rather the external regions of the
GCs. Indeed, mass segregation effects are less severe in external
cluster regions due to a lower density of particles, which in turn
limits the amount of dynamical interactions that lead to energy
equipartition.

Nonetheless, if usual DM particles or sub-structures have
masses considerably lower or higher than stars, our simulations
could then fail to spot possible mass segregation effects spe-
cially related to the cluster internal dynamics. Meanwhile, the
outer dynamics should not be considerably affected (as explained
above), such that our results would remain valid up to a reason-
able extent.

Mashchenko & Sills (2005b) studied the aspects of clus-
ter relaxation, using a particle mass ratio of $m_{\star\bullet}/m_{\bullet} = \{0.044, 0.44, 0.88\}$. They found that because the inner regions of the
GC+DM system are dominated by stars (as we confirm in
this study), the dynamics in the cluster's core is not severely
changed, and the impact of DM only starts to be visible at the
radius where the amount of DM particles and stars are equiv-
alent. At this point, energy equipartition leads the more mas-
sive DM particles to scatter the surrounding stars, creating a
cutoff radius similar to what is observed in the density pro-
files of many Milky Way GCs (e.g. Trager et al. 1995). Equiv-
alently, one could expect that less massive DM particles would
themselves be scattered and deplete the mini-halo faster, thus
decreasing the effectiveness of the shielding mechanism we
demonstrated.

4.4.3. Different orbits

When comparing simulated GCs with and without a DM mini-
halo in a tidal field, the important mass contribution from the
mini-halo is such that the orbits of these two scenarios are
bound to be different, as a consequence of dynamical friction
(e.g. upper panels of Fig. 3). In such setup, it can be some-
times difficult to disentangle which dynamical effects are due
to the DM mini-halo, or to the ongoing tidal interactions at a
given time step. For that reason, we choose to analyse the evo-
lution of these two scenarios starting from the same distance
relative to the host, in an attempt to provide a fairer compari-
son. This forcibly leads clusters formed inside a DM mini-halo
to inner regions of the host's potential, where tidal forces are
stronger.

Hence, all indications of milder tidal effects in GCs embed-
ded in DM add robustness to our conclusion of a DM shielding
mechanism. This is particularly true for the DM consequences
we highlighted on stellar streams, velocity anisotropy, rotation
and ellipticity. However, populating our simulations with more
GCs would undeniably provide a more complete picture on the
distribution of GC positions in IOM space at a given time, for
instance. Also, by simulating clusters at higher distances from
their host, they would suffer less dynamical friction<fnref n="23" />, which
could allow us to compare velocity dispersion profiles of clus-
ters at more similar orbital radii. Due to the high cost of our

<footnote n="21">²¹ Their minimum sub-halo mass is of $5 \times 10^8 M_\odot$, much higher than the values considered here.</footnote>

<footnote n="22">²² The filling factor is the ratio of the Jacobi to half-mass radii: $r_J/r_h$.</footnote>

<footnote n="23">²³ This is because at higher distances, there are less host's particles for the clusters to interact with, and thus exchange energy.</footnote>
simulations, we restrained such study to the five observed clus-
ters in the Fornax dSph, so as to approach a realistic case, but
adding similar clusters at a broader range of orbital radii seems
to be a tantalising direction for future studies.

**5. Summary and conclusions**

We have used full *N*-body simulations on GPU of the globular
cluster system in the Fornax dwarf spheroidal galaxy to probe
the differences from globular clusters originally embedded in a
dark matter mini-halo, and those devoid of it. For that, we sim-
ulated the two cases, where globular clusters have or not an
enveloping dark matter mini-halo, using the initial conditions
from Boldrini et al. (2020a), where all tracers (stars and dark
matter) are assumed as *N*-body particles, for both the globular
clusters and the host galaxy. Although the extra mass of the sys-
tem containing dark matter accelerates orbital decay, sinking the
clusters to regions of stronger tidal fields, the same extra mass
triggers an increase of the tidal radius, which generally goes well
beyond the limiting region of globular cluster bound stars (see
Fig. 1). As a result, this tidal radius growth causes the dark mat-
ter mini-halo to work as a protective shield (see Fig. 2), which is
itself stripped beforehand globular cluster stars. Below, we sum-
marise some of the main diagnostics from our work.

- Globular clusters with a dark matter mini-halo should have, in general, more radial outer velocity anisotropy profiles, throughout all their orbit.

- Inflations on the outer velocity dispersion profile are expected in both clusters with and without a dark matter mini-halo, as a result of the host's tidal field.

- Globular clusters with a dark matter mini-halo should have smaller degrees of internal rotation, and as a consequence, smaller ellipticities for their stellar distribution.

- Globular clusters without a dark matter mini-halo develop more prominent stellar streams, and present a more diffuse stellar distribution.

- Due to stronger dynamical friction, clusters originally embedded in a dark matter mini-halo evolve towards considerably lower energy regions of the integrals of motion space (differently from their dark matter-free counterparts), and cannot be reliably associated to previous accretion events due to their fast displacements in such space.

- Even clusters that retain large amounts of dark matter tend to have inner densities vastly dominated by baryonic components (up to a factor of a hundred in mass), and dark matter searches should target rather the cluster's outskirts.

- We provide parametric mass profiles for dark matter mini-halos that evolved in a tidal field, which is of much usefulness for dark matter mini-halo searches employing Jeans mass-modelling (e.g. Carlberg & Grillmair 2022).

All in all, our work helps to confirm previous predictions (e.g. the shielding mechanism from dark matter mini-halos), but also challenges a few preceding points. For instance, we can mention the shape of stellar velocity dispersion profiles, which we show to be more dependent on the tidal field rather than the presence of a mini-halo. In addition, we also point to the long-term survival of clusters embedded in dark matter, which in some scenarios can be disrupted faster due to pronounced dynamical friction, since the latter leads them to regions of stronger tidal forces, where the shielding mechanism will vanish in shorter timescales, all while leaving the cluster closer to its host. Such new predictions are due mostly to our better handling of tidal interactions and dynamical friction, respectively. We also provide the first

direct analysis pertaining to the evolution of the stellar veloc-
ity anisotropy, rotation, ellipticity and position in integrals of
motion space for the dark matter mini-halo scenario. Finally, we
evaluate some of the caveats of our analyses, and point towards
novel aspects that can be further explored in more detail with the
help of numerical simulations.

**Acknowledgements.** We thank the anonymous referee for the constructive report, with many insightful comments that have helped us to improve the quality of our results and clarify some descriptions in the manuscript. We also thank Geraint Lewis and Amaury Micheli for useful tips concerning the 3D visualisation of stellar streams. Eduardo Vitral was funded by an AMX doctoral grant from École Polytechnique. This work was supported by the EXPLORAGRAM Inria AeX grant. We greatly benefited from the public software PYTHON (Van Rossum & Drake 2009) packages BALROGO (Vitral 2021), SCIPY (Jones et al. 2001), NUMPY (van der Walt et al. 2011) and MATPLOTLIB (Hunter 2007). We also used GEOGEBRA (Hohenwarter 2002), as well as the SPYDER Integrated Development Environment (Raybaut 2009).

**References**

Abadi, M. G., Navarro, J. F., & Steinmetz, M. 2006, MNRAS, 365, 747

Akiyama, K. 1991, Earth Moon Planets, 54, 203

Balbinot, E., & Gieles, M. 2018, MNRAS, 474, 2479

Banik, N., & Bovy, J. 2021, MNRAS, 504, 648

Baumgardt, H., & Hilker, M. 2018, MNRAS, 478, 1520

Baumgardt, H., Côté, P., Hilker, M., et al. 2009, MNRAS, 396, 2051

Baumgardt, H., Hilker, M., Sollima, A., & Bellini, A. 2019, MNRAS, 482, 5138

Bekki, K., & Freeman, K. C. 2003, MNRAS, 346, L11

Bekki, K., & Yong, D. 2012, MNRAS, 419, 2063

Bertin, G., & Varri, A. L. 2008, ApJ, 689, 1005

Bianchini, P., Varri, A. L., Bertin, G., & Zocchi, A. 2013, ApJ, 772, 67

Bianchini, P., Sills, A., & Miholics, M. 2017, MNRAS, 471, 1181

Binney, J. 1980, MNRAS, 190, 873

Binney, J., & Mamon, G. A. 1982, MNRAS, 200, 361

Binney, J., & Tremaine, S. 2011, Galactic Dynamics: Second Edition, Princeton Series in Astrophysics (Princeton: Princeton University Press)

Boberg, O. M., Vesperini, E., Friel, E. D., Tiongco, M. A., & Varri, A. L. 2017, ApJ, 841, 114

Boldrini, P., & Bovy, J. 2022, MNRAS, 516, 4560

Boldrini, P., & Vitral, E. 2021, MNRAS, 507, 1814

Boldrini, P., Mohayaee, R., & Silk, J. 2020a, MNRAS, 492, 3169

Boldrini, P., Mohayaee, R., & Silk, J. 2020b, MNRAS, 495, L12

Bond, J. R., Cole, S., Efstathiou, G., & Kaiser, N. 1991, ApJ, 379, 440

Bovy, J. 2015, ApJS, 216, 29

Bromm, V., & Clarke, C. J. 2002, ApJ, 566, L1

Carlberg, R. G., & Grillmair, C. J. 2022, ApJ, 935, 14

Carlberg, R. G., & Keating, L. C. 2022, ApJ, 924, 77

Ciotti, L., & Bertin, G. 1999, A&A, 352, 447

Cohen, J. G., Blakeslee, J. P., & Ryzhov, A. 1998, ApJ, 496, 808

Conroy, C., Loeb, A., & Spergel, D. N. 2011, ApJ, 741, 72

Daniel, K. J., Heggie, D. C., & Varri, A. L. 2017, MNRAS, 468, 1453

de Boer, T. J. L., & Fraser, M. 2016, A&A, 590, A35

Dejonghe, H. 1987, MNRAS, 224, 13

Errani, R., Navarro, J. F., Ibata, R., et al. 2022, MNRAS, 514, 3532

Fabricius, M. H., Noyola, E., Rukdee, S., et al. 2014, ApJ, 787, L26

Feng, L., Yuan, Q., Yin, P.-F., Bi, X.-J., & Li, M. 2012, JCAP, 04, 030

Forbes, D. A., Brodie, J. P., & Grillmair, C. J. 1997, AJ, 113, 1652

Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

Gaia Collaboration (Helmi, A., et al.) 2018a, A&A, 616, A12

Gaia Collaboration (Babusiaux, C., et al.) 2018b, A&A, 616, A10

Gaia Collaboration (Brown, A. G. A., et al.) 2021, A&A, 649, A1

Gardner, J. P., Mather, J. C., Clampin, M., et al. 2006, Space Sci. Rev., 123, 485

Gieles, M., Erkal, D., Antonini, F., Balbinot, E., & Peñarrubia, J. 2021, Nat. Astron., 5, 957

Harris, W. E. 1991, ARA&A, 29, 543

Harris, W. E. 2009, ApJ, 699, 254

Hayashi, E., Navarro, J. F., Taylor, J. E., Stadel, J., & Quinn, T. 2003, ApJ, 584, 541

Helmi, A., & de Zeeuw, P. T. 2000, MNRAS, 319, 657

Hohenwarter, M. 2002, Master's Thesis, Paris Lodron University, Salzburg, Austria
Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90

Hurst, T. J., Zentner, A. R., Natarajan, A., & Badenes, C. 2015, Phys. Rev. D,
91, 103514

Ibata, R., Nipoti, C., Sollima, A., et al. 2013, MNRAS, 428, 3648

Ibata, R., Malhan, K., Martin, N., et al. 2021, ApJ, 914, 123

Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open Source Scientific
Tools for Python, http://www.scipy.org

Kamann, S., Husser, T. O., Dreizler, S., et al. 2018, MNRAS, 473, 5591

Kazantzidis, S., Mayer, L., Mastropietro, C., et al. 2004, ApJ, 608, 663

Keller, B. W., Kruijssen, J. M. D., Pfeffer, J., et al. 2020, MNRAS, 495, 4248

Khalaj, P., & Baumgardt, H. 2016, MNRAS, 457, 479

Kimm, T., Cen, R., Rosdahl, J., & Yi, S. K. 2016, ApJ, 823, 52

King, I. 1962, AJ, 67, 471

Kremer, K., Ye, C. S., Rui, N. Z., et al. 2020, ApJS, 247, 48

Kruijssen, J. M. D. 2014, Class. Quant. Grav., 31, 244006

Kruijssen, J. M. D. 2015, MNRAS, 454, 1658

Kuzma, P. B., Da Costa, G. S., Mackey, A. D., & Roderick, T. A. 2016, MNRAS,
461, 3639

Kuzma, P. B., Da Costa, G. S., & Mackey, A. D. 2018, MNRAS, 473, 2881

Lacey, C., & Cole, S. 1993, MNRAS, 262, 627

Lançon, A., Larsen, S., Voggel, K., et al. 2021, ArXiv e-prints
[arXiv:2110.13783]

Lane, R. R., Kiss, L. L., Lewis, G. F., et al. 2010, MNRAS, 406, 2732

Lanzoni, B., Ferraro, F. R., Mucciarelli, A., et al. 2018, ApJ, 861, 16

Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, ArXiv e-prints
[arXiv:1110.3193]

Longmore, S. N., Kruijssen, J. M. D., Bastian, N., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 291

Ma, X., Grudić, M. Y., Quataert, E., et al. 2020, MNRAS, 493, 4315

Mackey, A. D., & Gilmore, G. F. 2003, MNRAS, 340, 175

Majewski, S. R., Patterson, R. J., Dinescu, D. I., et al. 2000, in Liege International Astrophysical Colloquia, eds. A. Noels, P. Magain, D. Caro, et al., 35, 619

Mamon, G. A., Tweed, D., Thuan, T. X., & Cattaneo, A. 2012, Astrophys. Space
Sci. Proc., 28, 39

Mamon, G. A., Biviano, A., & Boué, G. 2013, MNRAS, 429, 3079

Marín-Franch, A., Aparicio, A., Piotto, G., et al. 2009, ApJ, 694, 1498

Martin, N. F., Venn, K. A., Aguado, D. S., et al. 2022, Nature, 601, 45

Mashchenko, S., & Sills, A. 2005a, ApJ, 619, 258

Mashchenko, S., & Sills, A. 2005b, ApJ, 619, 243

Massari, D., Koppelman, H. H., & Helmi, A. 2019, A&A, 630, L4

Merritt, D. 1985, AJ, 90, 1027

Miki, Y., & Umemura, M. 2017, New Astron., 52, 65

Miki, Y., & Umemura, M. 2018, MNRAS, 475, 2269

Montuori, M., Capuzzo-Dolcetta, R., Di Matteo, P., Lepinette, A., & Miocchi, P. 2007, ApJ, 659, 1212

Moore, B. 1996, ApJ, 461, L13

Mowla, L. A., Iyer, K. G., Desprez, G., et al. 2022, ApJ, 937, L35

Navarro, J. F., Frenk, C. S., & White, S. D. M. 1996, ApJ, 462, 563

Osipkov, L. P. 1979, Sov. Astron. Lett., 5, 42

Pechetti, R., Seth, A., Kamann, S., et al. 2022, ApJ, 924, 48

Peebles, P. J. E. 1984, ApJ, 277, 470

Peebles, P. J. E., & Dicke, R. H. 1968, ApJ, 154, 891

Peñarrubia, J., Varri, A. L., Breen, P. G., Ferguson, A. M. N., & Sánchez-Janssen,
R. 2017, MNRAS, 471, L31

Piatti, A. E., & Carballo-Bello, J. A. 2020, A&A, 637, L2

Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077

Planck Collaboration XIII. 2016, A&A, 594, A13

Planck Collaboration VI. 2020, A&A, 641, A6

Plummer, H. C. 1911, MNRAS, 71, 460

Portegies Zwart, S. F., McMillan, S. L. W., & Gieles, M. 2010, ARA&A, 48, 431

Prada, F., Klypin, A. A., Cuesta, A. J., Betancort-Rijo, J. E., & Primack, J. 2012, MNRAS, 423, 3018

Raybaut, P. 2009, Spyder-documentation, pythonhosted.org

Read, J. I., Mamon, G. A., Vasiliev, E., et al. 2021, MNRAS, 501, 978

Ricotti, M., Parry, O. H., & Gnedin, N. Y. 2016, ApJ, 831, 204

Rubin, V. C., & Ford, W. K., Jr. 1970, ApJ, 159, 379

Rubin, V. C., Ford, W. K., Jr., & Thonnard, N. 1980, ApJ, 238, 471

Saitoh, T. R., Koda, J., Okamoto, T., Wada, K., & Habe, A. 2006, ApJ, 640, 22

Searle, L., & Zinn, R. 1978, ApJ, 225, 357

Sérsic, J. L. 1963, Bull. Assoc. Argent. Astron., 6, 41

Sersic, J. L. 1968, Atlas de Galaxias Australes (Cordoba, Argentina:
Observatorio Astronomico)

Shin, J., Kim, S. S., & Lee, Y.-W. 2013, J. Korean Astron. Soc., 46, 173

Silk, J., & Mamon, G. A. 2012, Res. Astron. Astrophys., 12, 917

Strigari, L. E., Frenk, C. S., & White, S. D. M. 2010, MNRAS, 408, 2364

Takahashi, K. 1995, PASJ, 47, 561

Taylor, A. N., Dye, S., Broadhurst, T. J., Benítez, N., & van Kampen, E. 1998, ApJ, 501, 539

Tiongco, M. A., Vesperini, E., & Varri, A. L. 2016a, MNRAS, 455, 3693

Tiongco, M. A., Vesperini, E., & Varri, A. L. 2016b, MNRAS, 461, 402

Tiongco, M. A., Vesperini, E., & Varri, A. L. 2018, MNRAS, 475, L86

Trager, S. C., King, I. R., & Djorgovski, S. 1995, AJ, 109, 218

Trenti, M., Padoan, P., & Jimenez, R. 2015, ApJ, 808, L35

van de Ven, G., & van der Wel, A. 2021, ApJ, 914, 45

van den Bergh, S. 2008, AJ, 135, 1731

van den Bosch, F. C., & Ogiya, G. 2018, MNRAS, 475, 4066

van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Comput. Sci. Eng., 13, 22

Van Rossum, G., & Drake, F. L. 2009, Python 3 Reference Manual (Scotts Valley,
CA: CreateSpace)

Vasiliev, E. 2019, MNRAS, 489, 623

Vitral, E. 2021, MNRAS, 504, 1355

Vitral, E., & Mamon, G. A. 2021, A&A, 646, A63

Vitral, E., Kremer, K., Libralato, M., Mamon, G. A., & Bellini, A. 2022, MNRAS, 514, 806

Walker, M. G., & Peñarrubia, J. 2011, ApJ, 742, 20

Webb, J. J., Bovy, J., Carlberg, R. G., & Gieles, M. 2019, MNRAS, 488, 5748

Wu, Y., Liu, C., Peng, E. W., et al. 2022, ApJ, 926, 149

Zavala, J., & Frenk, C. S. 2019, Galaxies, 7, 81

Zhang, S., Mackey, D., & Da Costa, G. S. 2022, MNRAS, 513, 3136

Zhao, H. 1996, MNRAS, 278, 488

Zocchi, A., Gieles, M., Hénault-Brunet, V., & Varri, A. L. 2016, MNRAS, 462, 696

Zwicky, F. 1933, Helv. Phys. Acta, 6, 110

Zwicky, F. 1937, ApJ, 86, 217
# Appendix A: Velocity anisotropy

The Bayesian computation of the velocity anisotropy in BALROGO follows the procedure described in the following. First, we consider the local velocity distribution $h(v, r)$ as a Gaussian product over all spherical coordinates (i.e. $r, \theta, \phi$), in a similar fashion than done in MAMPOSST-PM (Mamon et al. 2013; Read et al. 2021)

$$h(v, r) = \mathcal{G}[v_r, \mu_r, \sigma_r(r)] \mathcal{G}[v_\theta, \mu_\theta, \sigma_\theta(r)] \mathcal{G}[v_\phi, \mu_\phi, \sigma_\phi(r)]. \quad (A.1)$$

<figure><img src="image_11.png" /><figcaption>**Fig. A.1.** Goodness of fit: Bayesian fit of the velocity anisotropy with BALROGO. The upper plot shows the case without a dark matter mini-halo, while the lower plot shows the case with a dark matter mini-halo. The black dots are the binned values of the empirical velocity anisotropy, logarithmically spaced and with error bars representing the 3-$\sigma$ Poisson error. The solid coloured line (red or green) represent the best likelihood fit of the Osipkov (1979), Merritt (1985) model (Eq. [2]), with respective shaded regions encompassing the [2.5, 97.5] percentiles of the fit.</figcaption></figure>

Differently than in MAMPOSST-PM, which solves the Jeans equation (Binney & Mamon 1982) to derive $\sigma_r(r)$, we fit $\sigma_r(r)$ with a generalisation of the parametric Plummer form from Dejonghe (1987, similarly used in Vasiliev 2019).

$$\sigma_r(r) = \frac{\sigma_0}{[1 + (r/r_\sigma)^\alpha]^\gamma}, \quad (A.2)$$

where $\sigma_0$ is the radial velocity dispersion at the cluster's centre, $r_\sigma$ a characteristic radius, $\alpha$ and $\gamma$ characteristic exponents, which we allow to vary between {1, 8}, and {0.05, 2}, respectively. Then, $\sigma_\theta(r)$ and $\sigma_\phi(r)$ can be derived, by assuming spherical symmetry and using the definition of the velocity anisotropy (Eq. [1]), such that $\sigma_\theta^2(r) = \sigma_\phi^2(r) = \sigma_r^2(r)[1 - \beta(r)]$. We then use the gOM expression for the anisotropy (Eq. [2]) in order to write $\beta(r)$. Since we do not use the Jeans equation, we can also account for the cluster's systemic motions, and so we assign $\mu_i$ as the mean of the velocities in the $i$ coordinate (differently from MAMPOSST-PM, which fixes $\mu_i = 0$).

Finally, we minimise the logarithm of the likelihood, which is basically the expression from Eq. A.1, multiplied over each tracer. We do it by using a Markov chain Monte Carlo (MCMC) approach from the PYTHON package EMCEE (Foreman-Mackey et al. 2013), with 10^4 steps^24 and 15 chains. In Figure A.1, we provide a goodness of fit of the velocity anisotropy of GC5 for the last snapshot from our simulations, for the case with and without a DM mini-halo.

# Appendix B: Velocity dispersion

We describe here how we constructed velocity dispersion profiles for the output of our simulations.

## B.1. Dispersion map

The velocity dispersion maps displayed in Figures 2 and 8 was constructed by first binning the projected X vs. Y map with PYTHON'S HEXBIN routine, setting the argument `gridsize=100`. For each bin, we then computed the 3D velocity dispersion of the DM component by summing quadratically the velocity dispersion on X, Y and Z directions. Whenever there were more than hundred particles inside the bin, we used all the bin's particles, otherwise we completed the sample by picking the closest particles to the bin's centre of mass, until the threshold of a hundred particles was attained.

The process above assured that each bin had a statistically significant number of tracers, which helped us to reach a better spatial resolution, eventually. However, as the results still presented an important amount of statistical noise, we decided to smooth our maps with PYTHON'S GAUSSIAN_FILTER routine, with the argument `sigma=3`. We finally displayed the outcome of this procedure in a map colour-coded logarithmically from blue (lower dispersion) to red (higher dispersion).

## B.2. Dispersion radial profiles

To build the radial dispersion profiles from Figure 4, we used the routines `ANGLE.CART_TO_SPH` and `DYNAMICS.DISPERSION` from the BALROGO PYTHON package, to first convert the data into spherical coordinates and then compute the velocity dispersion of the radial component, as a function of the distance to the cluster's centre. The way we set BALROGO to compute this dispersion is by first dividing the radial extent into thirty equally spaced logarithmic bins and then calculating the respective dispersion and Poisson error associated to it.

Next, it smooths the profile in order to remove statistical noise, by fitting a 10-th order polynomial to it with the NUMPY.POLYFIT routine. 1-$\sigma$ regions were obtained by considering the covariance matrix returned by the same algorithm. We finally chose to consider the results inside a more restrict spatial extent in order to neglect the potentially bad fits of the polynomial on the borders of our data.

# Appendix C: Asymmetric surface density

Here we describe how the BALROGO's POSITION.ELLIPSE_LIKELIHOOD routine performs asymmetric Sérsic fits with $(\alpha, \delta)$ data. First, it projects the $(\alpha, \delta)$ data according to classical spherical trigonometry relations, translating it so it is centred at the origin:

$$x_p = \cos \delta \sin (\alpha - \alpha_0), \quad (C.1a)$$

$$y_p = \sin \delta \cos \delta_0 - \cos \delta \sin \delta_0 \cos (\alpha - \alpha_0), \quad (C.1b)$$

where $(\alpha_0, \delta_0)$ is the centre of the cluster. Next, it rotates the axis so the data can be easily handled:

$$x = x_p \cos \theta + y_p \sin \theta, \quad (C.2a)$$

$$y = -x_p \sin \theta + y_p \cos \theta, \quad (C.2b)$$

where $\theta$ is the angle between the original reference frame and the new one. With this new set, we are able to define a likelihood

<footnote n="24">24 The first 5000 were discarded as the initial burn-in phase.</footnote>
function of the stellar distribution as:

$$ \mathcal{L} = \prod_{i} \frac{\Sigma(m)}{N_{\text{tot}}}. \qquad (C.3) $$

with $m = \sqrt{(x/a)^2 + (y/b)^2}$, and where (a, b) are the semi-axis of the ellipse. The surface density $\Sigma(m)$ and the number of tracers at infinity $N_{\text{tot}}$ are defined such as in van de Ven & van der Wel (2021), so that we have:

$$ \frac{\Sigma(m)}{N_{\text{tot}}} = \frac{b_n^{2n} \exp[-b_n m^{1/n}]}{2\pi a b n \Gamma(2n)}. \qquad (C.4) $$

In equation C.4, n is the Sérsic index, $\Gamma(x)$ is the gamma function of the variable x and the $b_n$ is computed with the precise approximation from Ciotti & Bertin (1999). Thus, the fit finds the parameters which maximise $\mathcal{L}$.

Finally, to derive statistical errors of our Bayesian estimates, we use PYTHON'S NUMDIFFTOOLS.HESSIAN method to compute the Hessian matrix of the probability distribution function (i.e. Eq. [C.4]). After, we assign the uncertainties of each parameter as the square root of the respective diagonal position of the inverted Hessian matrix.

## Appendix D: Extra material

<figure><img src="image_12.png" /><figcaption>**Fig. D.1.** Evolution of scale radii: Evolution of the GC+DM system tidal radius ($r_t$, blue) and of the half-mass radius of the bound dark matter particles ($r_{1/2,\bullet}$, black), for each simulated cluster. We observe that clusters who lost most of their dark matter mini-halo (i.e. GC1 and GC2) presented $r_t > r_{1/2,\bullet}$ at some point of their orbital evolution, similarly to the relation for sub-halo disruption from Hayashi et al. (2003).</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. D.2.** Evolution of the dark mass: Evolution of the total mass of bound dark matter particles, for each simulated cluster. Bound particles were defined according to Section 2.3.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Fig. D.3.** Evolution of the stellar mass: Evolution of the total stellar mass of bound globular cluster particles, for each simulated cluster. Clusters originally embedded in a dark matter mini-halo (i.e. ★●) are displayed in green, while those devoid of the mini-halo (i.e. ★) are represented in red. Bound particles were defined according to Section 2.3.</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Fig. D.4.** Evolution of the half mass radius: Evolution of the half mass radius of bound globular cluster particles, for each simulated cluster. Clusters originally embedded in a dark matter mini-halo (i.e. ★●) are displayed in *green*, while those devoid of the mini-halo (i.e. ★) are represented in *red*. Bound particles were defined according to Section 2.3.</figcaption></figure>

<table><caption>**Table D.1.** Best fit parameters.</caption><thead><tr><th></th><th>a<sub>Z</sub><br>[pc]</th><th>γ<sub>Z</sub></th><th>f<sub>Z</sub></th><th>a<sub>K</sub><br>[pc]</th><th>γ<sub>K</sub></th><th>f<sub>K</sub></th></tr></thead><tbody><tr><td>GC 1</td><td>67.2</td><td>0.66</td><td>0.03</td><td>28.9</td><td>0.27</td><td>0.03</td></tr><tr><td>GC 2</td><td>58.1</td><td>0.93</td><td>0.03</td><td>27.3</td><td>0.57</td><td>0.02</td></tr><tr><td>GC 3</td><td>175.6</td><td>0.94</td><td>0.13</td><td>101.0</td><td>0.94</td><td>0.13</td></tr><tr><td>GC 4</td><td>92.0</td><td>0.92</td><td>0.05</td><td>46.1</td><td>0.82</td><td>0.04</td></tr><tr><td>GC 5</td><td>119.8</td><td>1.03</td><td>0.06</td><td>66.7</td><td>0.96</td><td>0.06</td></tr></tbody></table>

**Notes.** Best likelihood parameters from our mass profile fits for each simulated globular cluster at $z = 0$. The first three columns yield the fits when using a Zhao (1996) $\alpha\beta\gamma$ model with fixed $\alpha = 2$ and $\beta = 5$ (as in the classic Plummer model) for better convergence, while the last three columns yield the fits when adopting a Kazantzidis et al. (2004) model. Each three columns display the scale radius $a$, inner slope $\gamma$ and mass fraction $f$ relative to the first dark matter shell explained in Section 4.2.2, defined as $M_{\star,1}$. For details, see Eqs. 5 and 6.