# On the effects of photoionization feedback on second-generation star formation in globular clusters of different masses

A Yaghoobi, Joakim Rosdahl, F Calura, S Ataiee

► To cite this version:

A Yaghoobi, Joakim Rosdahl, F Calura, S Ataiee. On the effects of photoionization feedback on second-generation star formation in globular clusters of different masses. Monthly Notices of the Royal Astronomical Society, 2024, 528 (4), pp.5477-5486. 10.1093/mnras/stad3771. hal-04759210

HAL Id: hal-04759210

https://cnrs.hal.science/hal-04759210v1

Submitted on 30 Oct 2024

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
# On the effects of photoionization feedback on second-generation star formation in globular clusters of different masses

A. Yaghoobi,<sup>1,2*</sup> J. Rosdahl,<sup>1,2</sup> F. Calura<sup>1,3</sup> and S. Ataiee<sup>1,4</sup>

Accepted 2023 December 11. Received 2023 November 3; in original form 2023 July 28

## ABSTRACT

We simulate the formation of second-generation (SG) stars in young clusters with masses of $10^5$ and $10^6 M_\odot$ within 30–100 Myr after the formation of clusters. We assume the clusters move through a uniform interstellar medium with gas densities of $10^{-24}$ and $10^{-23} \text{ g cm}^{-3}$ and consider the stellar winds from asymptotic giant branch (AGB) stars, gas accretion on to the cluster, ram pressure, star formation, and photoionization feedback of our stellar systems including binary stars. We find that SG stars can be formed only within the $10^6 M_\odot$ cluster in the high-density simulation, where the cluster can accrete sufficient pristine gas from their surrounding medium, leading to efficient cooling required for the ignition of SG formation and sufficient dilution of the AGB ejecta. Hence, our results indicate that a denser environment is another requirement for the AGB scenario to explain the presence of multiple populations in globular clusters. On the other hand, the ionizing feedback becomes effective in heating the gas in our low-density simulations. As a result, the clusters cannot accumulate a considerable amount of pristine gas at their centre. The gas mass within the clusters in these simulations is similar to that in young massive clusters (YMCs). Hence, our studies can provide a possible reason for the lack of gas, star formation, and SG stars in YMCs. Our results indicate that the ionizing stellar feedback is not a severe problem for SG formation; rather, it can help the AGB scenario to account for some observables.

**Key words:** hydrodynamics – radiative transfer – methods: numerical – stars: formation – globular clusters: general.

## 1 INTRODUCTION

Globular clusters (GCs), containing low-metallicity ([Fe/H] $\le 1$), low-mass ($\le 0.75 M_\odot$), and old ($\gtrsim 1 \text{ Gyr}$) stars (Cadelano et al. 2020), can provide useful information as early star-forming sub-systems to study the formation history of their parent galaxies (Larson 1996). Recent observations have revealed that GC stars display anomalous variations in light elements (e.g. He, C, N, O, Na, Al, and sometimes Mg) that are not expected to be due to low-mass stars during stellar evolutionary processes (Prantzos, Charbonnel & Iliadis 2007; Carretta et al. 2009; Piotto et al. 2015). Together with another population characterized by a composition similar to normal field stars, such anomalous stars build multiple stellar populations (MSPs) observed in GCs. The anomalous stars are enriched in N and Na and depleted in C and O, whereas the sum of the C, N, and O abundances inside these stars are observed to be generally constant within the measurement errors (Dickens et al. 1991; Yong, Grundahl & Norris 2015). Additionally, except for some massive GCs such as $\omega$ Cen (Johnson & Pilachowski 2010), no considerable spread in Fe abundance and heavy elements has been found in GCs hosting MSPs (Gratton, Sneden & Carretta 2004; Gratton, Carretta & Braglia 2012). This specific abundance pattern makes it difficult to explain how

MSPs have been formed in GCs; hence, their formation remains a mystery.

Such an abundance pattern is expected from the yields of hot hydrogen burning within intermediate-mass and massive stars (Renzini et al. 2015). Hence, most scenarios proposed for the formation of MSPs are based on the ejecta of these stars, fuelling the formation of the SG stars, where MSPs are subsequent generations. The population with the same composition as the field stars at the same metallicity [Fe/H] is referred to as the first generation (FG) or first population. Stars within GCs showing variation in the light elements compared to the first population are commonly categorized as the second generation (SG) or second population (Bastian & Lardo 2018; Milone & Marino 2022). To date, the suggested FG candidates are intermediate-mass asymptotic giant branch (AGB) stars ($4 \le m/M_\odot \le 8$; D’Ercole et al. 2008; Conroy & Spergel 2011; D’Ercole, D’Antona & Vesperini 2016; Bekki 2017), fast-rotating massive stars (FRMS, $20 \le m/M_\odot \le 120$; Decressin et al. 2007a; Decressin, Charbonnel & Meynet 2007b; Krause et al. 2013), and supermassive stars ($m/M_\odot \ge 10^4$; Denissenkov & Hartwick 2014; Denissenkov et al. 2015). Alternative scenarios propose different mechanisms in which the MSPs may be coeval. For example, in the early disc accretion scenario, the FG low-mass stars accrete enriched stellar winds from more massive interacting binary stars, resulting in different abundances (Bastian et al. 2013a).

Each proposed scenario has its own advantages and drawbacks, but none of them has been able to explain all the observed properties

<footnote n="1"><sup>1</sup>Department of Physics, Faculty of Sciences, Ferdowsi University of Mashhad, Mashhad 91775-1436, Iran</footnote>

<footnote n="2"><sup>2</sup>Université Claude Bernard Lyon 1, CRAL UMR5574, ENS de Lyon, CNRS, Villeurbanne, F-69622, France</footnote>

<footnote n="3"><sup>3</sup>INAF – OAS, Osservatorio di Astrofisica e Scienza dello Spazio di Bologna, via Gobetti 93/3, I-40129 Bologna, Italy</footnote>

<footnote>*E-mail: asiyeh.yaghoobi@gmail.com</footnote>
of MSPs (see Renzini et al. 2015; Bastian & Lardo 2018), such as the observed chemical anomalies. To solve this issue, the ejecta need to be diluted by pristine gas, i.e. the gas with the same composition as the one that formed the FG stars. In addition, most models, such as AGB and FRMS scenarios, encounter a ‘mass budget problem’ to produce a cluster with an SG fraction similar to observed values, between ~30 and 90 per cent for Galactic GCs (Milone et al. 2017). One problem is that, assuming a standard initial mass function (IMF), the FG ejecta do not contain enough mass to make up the entirety of the observed SG mass. This problem might be alleviated in specific conditions, such as a top-heavy IMF for FG stars and a binary fraction close to unity (D’Ercole et al. 2008; Khalaj & Baumgardt 2015; Krause et al. 2016; Vesperini et al. 2021). While generally not clear from observations, the epoch of SG formation and age difference between populations depend on the proposed scenario. In models based on massive and very massive stars, SG formation should occur within a short time-scale after the formation of FG stars and before the start of Type II supernova (SN) explosions to avoid the incorporation of newly produced Fe and a significant metallicity spread (e. g. Krause et al. 2013). On the other hand, in the AGB scenario, the SG formation starts after the end of SNe (about 30 Myr) and lasts for roughly a few ten Myrs (D’Ercole et al. 2008; D’Ercole, D’Antona & Vesperini 2016; Calura et al. 2019).

The presence of MSPs with anomalous chemical compositions is found to be specific to only massive GCs older than 2 Gyr (Bastian & Lardo 2018). It is not clear yet whether young massive clusters (YMCs; with a mass range of $[10^4–10^5.5] M_\odot$ and ages $\lesssim 100$ Myr

host MSPs with such anomalous chemical compositions (Portegies Zwart, McMillan & Gieles 2010). Nevertheless, observationally these clusters do not show any significant amount of gas (Bastian et al. 2013b; Bastian, Hollyhead & Cabrera-Ziri 2014; Cabrera-Ziri et al. 2015; Hollyhead et al. 2015) and any evidence for ongoing star formation (SF; Bastian et al. 2013b). As their age can be compared with the proposed epoch for SG formation in scenarios, in principle these facts can constrain SG formation models (Bastian et al. 2013b; Bastian & Strader 2014; Bastian, Hollyhead & Cabrera-Ziri 2014), assuming that YMCs are equivalent to proto-GCs. But some rely on the idea that the SG formation is limited to only the early Universe, not the present-day one (Renzini et al. 2015; D’Ercole, D’Antona & Vesperini 2016). For example, D’Ercole, D’Antona & Vesperini (2016) show that the AGB scenario needs ambient surrounding properties of GCs consistent with those of star-forming disc galaxies observed at redshifts $> 2$. Moreover, the fact that YMCs are almost gas free (Bastian et al. 2013b; Bastian, Hollyhead & Cabrera-Ziri 2014; Cabrera-Ziri et al. 2015; Hollyhead et al. 2015) indicates that they efficiently expel their intracluster medium (ICM), despite significant mass-loss from their stars and their capability to accrete gas from the ambient medium. This raises another question: under which physical conditions can a young cluster get rid of its ICM? To answer this question, we require more comprehensive theoretical studies of the gas dynamics in young clusters, including the key processes, including SN feedback, ram pressure, radiative feedback, cluster potential, and stellar wind heating (Conroy & Spergel 2011; Calura et al. 2015; Gavagnin et al. 2017; Wünsch et al. 2017; Naiman, Ramirez-Ruiz & Lin 2018; Chantereau et al. 2020).

Since a few years, we have focused on the gas accumulation and SG formation in young ($> 30$ Myr) clusters with different masses [$10^5–10^7 M_\odot$] based on the AGB scenario (Calura et al. 2019; Lacchin, Calura & Vesperini 2021; Lacchin et al. 2022; Yaghoobi et al. 2022a, b). In this scenario, the origin of SG stars

is a mixture of the AGB stellar winds and pristine gas accumulated in clusters. D’Ercole, D’Antona & Vesperini (2016) investigate the origin of the pristine gas and its dynamics during the phase of SG formation. Their model follows three requirements: (i) the pristine gas should initially be removed from the cluster by Type II SN explosions within about 30 Myr after the birth of the FG of stars; (ii) the cluster should be able to accrete the pristine gas within about 60 Myr; (iii) the SN II ejecta should not contaminate the pristine gas. They show that these conditions can explain the origin of SG stars only in massive clusters that formed in the discs of galaxies at $z > 2$ (Kravtsov & Gnedin 2005; Kruijssen 2015).

Using 3D hydrodynamical simulations, we investigated this model for clusters of different masses in Calura et al. (2019) and Yaghoobi et al. (2022a; hereafter Paper I). We performed a parameter study of the gas density of the interstellar medium (ISM) with values $10^{-24}$ and $10^{-23} g cm^{-3}$, roughly corresponding to the values observed in star-forming galaxies at low and high redshifts (Marcolini, Brighenti & D’Ercole 2003; D’Ercole, D’Antona & Vesperini 2016; Wardlow et al. 2017) to explore the effect of different environments on the formation of SG stars. We found that the properties of the SG population depend significantly on the surrounding environment of the clusters. In the high-density simulations, massive clusters (with masses $\ge 10^6 M_\odot$) could overcome the ram pressure, accrete pristine gas, and retain their own stellar winds, producing a massive SG. Positive correlations were also found for the SG mass fraction and maximum H enhancement in GCs versus the cluster mass.

In Yaghoobi et al. (2022b, hereafter Paper II), we included ionizing radiation from stars in our simulations and found that the gravitational potential of a very massive cluster ($10^7 M_\odot$, with half-mass radius 30 pc) is strong enough to accumulate both accreted pristine gas and AGB ejecta at its centre. Photoionization heating was found to be efficient in warming up and expanding the gas for only about the first 50 Myrs from the FG formation. After that, the gas cooling dominates the radiative heating and the gas accumulated in the innermost regions results in SF. However, due to ionizing feedback effects, the total SG mass became a factor two lower than in the run without ionizing radiation at the end of simulations (100 Myr after the cluster formation). Nevertheless, the SG mass fraction within the cluster reached 40 per cent in the high-density model, comparable with observed fractions for Galactic GCs (Milone et al. 2017). However, in the low-density case, the gas content and SG fraction were very low (~4 per cent). To further explore the viability of this scenario, our study including the photoionization feedback needs to be expanded to clusters with different structural parameters, since it is shown that the cooling and heating processes depend on the properties of the cluster (Wünsch et al. 2017; Yaghoobi et al. 2022a). Particularly, we are interested in studying the SF in the clusters expected to be equivalent to proto-GCs and YMCs to explore the role of stellar feedback in SG formation. Moreover, we investigate how results are sensitive to the different environments, which are expected to be different for GCs and YMCs.

This study is a follow-up to Paper II to investigate the SF and gas content in young clusters with different properties (mass and half-mass radius). We aim at checking whether these clusters can retain their stellar winds, accrete the pristine gas from the ambient medium, and cool their ICM in the presence of ionizing radiation. The paper is organized as follows. In Section 2, we describe the simulation set-up and main assumptions of our model. In Section 3, we present the results of our simulations, and in Section 4 we discuss them. Finally, in Section 5, we list the main conclusions of our study.
<table><caption>**Table 1.** Main parameters of our simulations for both cases of RT and NoRT. Column description: $M_{FG}$ is the mass of the cluster; $\rho_{pg}$ the gas density of the pristine gas; $t_l$ the time when the cluster reaches the ISM (Section 2.1) at the time reference of this paper (the FG formation). Note that the AGB ejecta in all simulations start at $t_{AGB} = 39$ Myr.</caption><thead><tr><th>Simulation</th><th><i>M</i><sub>FG</sub> (M<sub>☉</sub>)</th><th><i>ρ</i><sub>pg</sub> (g cm<sup>-3</sup>)</th><th><i>t</i><sub>l</sub> (Myr)</th></tr></thead><tbody><tr><td>M5INFALL24</td><td>10<sup>5</sup></td><td>10<sup>-24</sup></td><td>33.9</td></tr><tr><td>M5INFALL23</td><td>10<sup>5</sup></td><td>10<sup>-23</sup></td><td>30.0</td></tr><tr><td>M6INFALL24</td><td>10<sup>6</sup></td><td>10<sup>-24</sup></td><td>42.5</td></tr><tr><td>M6INFALL23</td><td>10<sup>6</sup></td><td>10<sup>-23</sup></td><td>33.9</td></tr></tbody></table>

## 2 SIMULATION SET-UP

The simulation set-up in this study is the same as Paper II for a cluster with the mass of $10^7$ M$_{☉}$ and half-mass radius of 23 pc, except that we study the SG formation within FG clusters with masses of $10^5$ and $10^6$ M$_{☉}$, and half-mass radius of 4 pc, moving through a homogeneous ISM with uniform density. These masses are comparable with the initial masses derived for GCs by calculating the cluster orbits backward in time (Baumgardt et al. 2019). To study the role of environments in SF formation, we consider different values for ISM densities; $10^{-24}$ g cm$^{-3}$ for a typical dwarf galaxy (Marcolini, Brighenti & D’Ercole 2003) and $10^{-23}$ g cm$^{-3}$ as a representative value for star-forming regions in galaxies at high redshifts (Wardlow et al. 2017). As in Paper II, we consider the effects of the ionizing FG and SG radiation using radiation-transfer simulations (RT simulations). We also run their analogues without radiation (NoRT) to better understand the role of stellar radiative feedback in SG formation. The main simulations of this study and their relevant characteristics are listed in Table 1.

We use RAMSES-RT (Rosdahl et al. 2013; Rosdahl & Teyssier 2015), the radiative hydrodynamics version of the RAMSES (Teyssier 2002) code. The code uses a second-order Godunov scheme to solve the Euler equations and a particle-mesh solver to compute the dynamical evolution of particles. In Paper I, the appropriate simulation volume for the $10^5$ and $10^6$ M$_{☉}$ clusters was found to be about 50 pc. However, we find that we need to consider a larger box to achieve convergence in results due to gas expansion resulting from photoionization heating. We find that a cube with a width of 128 pc and inflow boundaries can fulfill this goal. Moreover, we use the Adaptive Mesh Refinement (AMR) strategy with a maximum size of $\Delta x_{max} = 1$ pc ($l_{min} = 7$ in the RAMSES code) and a minimum size $\Delta x_{min} = 0.062$ pc ($l_{max} = 11$). Cells are refined to smaller sizes if their mass exceeds $10$ M$_{☉}$. As in Paper I, we performed a series of tests assuming different resolutions, ranging from $\Delta x_{min} = 1$ pc to $\Delta x_{min} = 0.062$, indicating that the chosen set-up is sufficient for our results to converge.

The simulation set-up of our model has been discussed extensively in Paper II. Here, we briefly recap it and encourage the reader to consult Paper I, Paper II, and Calura et al. (2019) for more details. In our notation, time $t = 0$ corresponds to the birth time of the FG of stars.

### 2.1 Initial conditions

We set the initial conditions of our simulations as proposed in D’Ercole, D’Antona & Vesperini (2016). Type II SN explosions of FG stars have formed a bubble of hot and diffuse gas around the cluster at the end of their activity, assumed to be 30 Myr after the birth of the FG. The bubble radius depends on the cluster mass, ISM density, and velocity with respect to the gas (see Calura et al.

2019). It is assumed that the bubble is filled with extremely diffuse gas in hydrostatic equilibrium with the cluster. Thus, FG stars are not initially in direct contact with the pristine ISM gas until they traverse this radius and reach the ISM. After that, the pristine gas encounters the cluster and starts to dilute the enriched AGB ejecta. We also assume that the injection of AGB ejecta starts at $t_{AGB} = 39$ Myr. Following this model, we consider the velocity of clusters to be $23$ km s$^{-1}$, of the order of the gas velocity dispersion in an ionized medium, and assume it to be constant during the simulations. Accordingly, the time at which the cluster reaches the ISM ($t_l$) can be estimated, as described in Paper I. $t_l$s for our simulations are computed by means of equation (11) of D’Ercole, D’Antona & Vesperini (2016) and presented in the fourth column of Table 1. With the exception of the M6INFALL24 model, the infall starts earlier than stellar winds in all our simulations. This is visible from the comparison of $t_l$ and $t_{AGB}$ in Table 1. All simulations finish at $t = 100$ Myr at the start of FG Type I SNe (D’Ercole et al. 2008).

The FG stars are modelled by a static Plummer (1911) density profile:

$$ \rho_{*,FG}(r) = \frac{3 M_{FG}}{4\pi r_{\text{p}}^3} \left(1 + \frac{r^2}{r_{\text{p}}^2}\right)^{-\frac{5}{2}}, $$

where $r$ is the distance from the centre of the cluster, $r_p = 3$ pc and $M_{FG}$ are the Plummer radius and the cluster mass, respectively. We neglect any changes in the FG stellar distribution due to cluster evolution during the simulation. We perform the simulations in the reference frame of the cluster. Therefore, we assume a static FG cluster at the centre of the box and, after a time of $t_l$, the ISM gas enters into the box from the left side with the velocity of the cluster ($23$ km s$^{-1}$). Hence, the cluster is exposed to the ram pressure of the incoming gas and can accumulate it to form new stars.

### 2.2 Stellar winds of the FG stars

Stellar winds from the AGB stars are assumed to start at $t_{AGB} = 39$ Myr after the formation of the FG stars. Following the models described in Calura et al. (2019), the rate of injected mass by AGB stars can be modelled as a function of time $t$ (expressed in yr) as:

$$ \dot{\rho}_{\text{AGB}}(r) = \alpha \rho_{*,\text{FG}}(r), $$

where $\alpha = 0.065$ t$^{-1.01}$yr$^{-1}$ is the specific injection rate. We add this rate as a source term into the continuity equation for gas density. We specify these ejecta with a different He abundance from the pristine gas. To do that, we define a passive scalar advected with the gas density and follow the evolution of helium (He) mass fraction for all the gas in the box. The He mass fraction of the ejecta is assumed to slightly change from $Y(t = 39$ Myr) = 0.36 to $Y(t = 100$ Myr) = 0.32 (Ventura & D’Antona 2011), whereas for the pristine gas we assume $Y = 0.25$.

### 2.3 Star formation

We use the SF model described in Rasera & Teyssier (2006) to turn gas mass into ‘star’ particles in eligible cells for SF based on the standard Schmidt (1959) law. Thus, a source term is added to the continuity equation as:

$$ \dot{\rho}_{\text{sf}} = -\frac{\rho}{t_{*}}, \qquad (1) $$

where the $t_*$ is the SF time-scale which is assumed to be 0.1 Gyr, as in Calura et al. (2019). In this paper, we adopt four physical criteria
for SF: (i) a gas temperature less than $2 \times 10^4$ K, (ii) a converging gas velocity ($\nabla \cdot \mathbf{v} < 0$), (iii) a local Jeans' length smaller than 4 (finest) cell width, i.e. less than 0.25 pc, and (iv) a density higher than $6 \times 10^{-23}$ g cm$^{-3}$. When these criteria are fulfilled, star particles are formed stochastically. The stellar particle mass is an integer multiple of $m_* = 0.4 M_\odot$, sampled from a Poisson probability distribution as described in Rasera & Teyssier (2006). The star particles are placed at the centre of their parent cell, with a velocity, metallicity, and He abundance equal to their natal gas.

## 2.4 Radiative transfer

The RAMSES-RT code is a momentum-based radiative transfer code that couples the hydrodynamics equations to the radiative transfer equations using the M1 closure for the Eddington tensor. It models the propagation of ionizing radiation and its interplay with gas via non-equilibrium thermochemistry for hydrogen and helium. The frequency range of photons is split into discrete photon groups, and then the radiative transfer equations can be solved separately for each group. In this study, we assume three photon groups (H1, He1, and He II ionizing photons) with mean energies of 18.2, 33.0, and 61.3 eV, as described in Paper II. Thus, in every cell, the evolution of ionization fractions for hydrogen and helium is followed (Rosdahl et al. 2013; Rosdahl & Teyssier 2015). Using this code, we study the effects of photoionization heating, cooling, and radiation pressure of ionizing radiation from the FG and SG stars as our luminous sources. Note though that additional simulations, not included in this paper, show that radiation pressure has no role in gas evolution and SF, and that all effects of the ionizing radiation presented in the next sections are due to photoionization heating. The number of injected photons by the FG and SG is estimated by a spectral energy distribution model, depending on the mass, age, and metallicity of stars. In our simulations, we use the binary population and spectral synthesis (BPASS) code (Eldridge et al. 2017) and assume a metallicity of $Z = 0.001$ and a standard Kroupa (2001) IMF for the FG stars. The same metallicity is considered for the SG stars but assuming a truncated IMF with maximum mass of $m = 8 M_\odot$ to neglect the SN feedback effects, as assumed in Calura et al. (2019) and Yaghoobi et al. (2022b). Moreover, we assume a reduced speed of light $0.002c$ to reduce the computational cost of simulations.

## 2.5 Cooling and heating

In RAMSES-RT, radiative heating and cooling processes contribute with positive and negative rates in the energy equation, respectively. They are functions of the gas temperature, density, photon properties, and ionization fractions, as described in Rosdahl et al. (2013). The cooling processes considered for hydrogen and helium are collisional ionization, excitations, recombinations, dielectronic recombinations, bremsstrahlung, and Compton cooling, and metal cooling.

Following Calura et al. (2019), we also include the heating effects of stellar winds from the AGB stars in the energy equation as:

$$S = 0.5\alpha\rho_{*,\text{FG}}(3\sigma^2 + v^2 + v_{\text{wind}}^2), \quad (2)$$

where $\sigma$ is the one-dimensional velocity dispersion of the cluster, $v_{\text{wind}}$ is the wind velocity of the AGB stars, and $v$ is the local gas velocity (D’Ercole et al. 2008). We assume a wind velocity of $v_{\text{wind}} = 20$ km s$^{-1}$ and an adiabatic index of $\gamma = 5/3$. The stellar winds from other sources, e. g. main-sequence stars, can be neglected due to their low mass-loss rates (Paper II).

## 3 RESULTS

Our recent studies (Calura et al. 2019, Paper I) showed that the gravitational potential of massive clusters ($10^6$ and $10^7 M_\odot$) can overcome ram pressure and stellar wind heating. As a result, clusters can accumulate a significant amount of gas within their central regions, where the accreted pristine gas dilutes the ejecta. In Paper II, we showed that stellar radiative feedback can also not resist against the potential of the cluster. However, it delays the formation of SG stars in massive clusters of $10^7 M_\odot$ and decreases the total SG mass by 50 per cent. In this section, we describe the results of our simulations listed in Table 1 with lower mass clusters and explore how the radiative feedback can affect the gas accumulation, of both retained AGB ejecta and accreted pristine gas. We set our simulations to start at $t = t_{\text{AGB}}$ for M6INFALL24 and at $t = t_1$ for three other simulations. Note that the 64 pc distance between the left side and the centre of the cluster causes a time delay about 1.8 Myr for incoming gas to reach the cluster centre. To take this into account, we allow the incoming gas to come into the simulation box at $t = t_1 - 1.8$ Myr so that it passes the centre at $t_1$.

### 3.1 Low-mass clusters

In Paper I, we found that without the effects of photoionization, the ram pressure exerted by the ambient medium can limit the ability of the gravitational potential of a $10^5 M_\odot$ cluster to retain its stellar winds and accrete the pristine gas. However, it did not suppress the gas accumulation at the centre entirely. As a result, SG stars with high-(≈0.06) and intermediate-He (≈0.12) enhancements were formed within the cluster centre in low- and high-density simulations, respectively. These results were in contrast with the He enhancements observed (≤0.01) in low-mass clusters (Milone et al. 2020). Here, we include the ionizing feedback in our simulations with $10^5 M_\odot$ to explore how this process affects the formation of SG stars.

#### 3.1.1 M5INFALL24

According to the model proposed by D’Ercole, D’Antona & Vesperini (2016), we assume that the cluster in this simulation has already traversed the bubble radius created by SN feedback and reached the ISM. Therefore, the incoming gas reaches the cluster’s centre at $t_1 = 33.9$ Myr (fourth column of Table 1), while the injection of mass and energy from the AGB stars begins at $t_{\text{AGB}} = 39$ Myr. Fig. 1 shows the density (upper panels) and temperature (lower panels) maps in the x-y plane for our RT and NoRT simulations with $10^5 M_\odot$ clusters at the final time (100 Myr). The arrows on the temperature maps indicate the velocity field of the gas. In the NoRT simulation (right column of Fig. 1, also discussed in Paper I), the cluster could retain some of the ejecta in the central regions, accrete a marginal amount of pristine gas, and results in SG stars with a total mass of about $3 \times 10^3 M_\odot$ formed from almost-pure AGB ejecta. The left panels show that the photoionizing feedback of FG stars has increased the temperature and pressure of the gas within the cluster, leading to its expansion in the central regions. In the upper-left panel of the figure, a shell is seen at a distance of about 60 pc from the centre, due to the expanding bubble hitting the incoming pristine gas. The bottom-left panel shows that the gas inside the shell is warm, with $T > 10^4$ K, and therefore ionized. Hence, this shell specifies the region affected by the FG radiation around the cluster. Moreover, the maps display a weak shock in the central regions due to the ICM hitting the ISM. Therefore, gas heating exceeds the gas cooling within the cluster, resulting in no SG SF in this simulation.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Gas density (top panels) and temperature (bottom panels) maps for the M5INFALL24 model in the cases with (RT, left column) and without (NoRT, right column) photoionization feedback at the end of simulation, i.e. at *t* = 100 Myr (after the FG formation). The black arrows represent the gas velocity field. Photoionization heating has caused the gas inside the cluster to expand. In contrast, radiative cooling is insufficient to cool the gas at the centre of the cluster. As a result, no stars can form within ~100 Myr in this low-mass clusters.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Figure 2.** The same as Fig. 1 but for the M5INFALL23 simulation.</figcaption></figure>

### 3.1.2 M5INFALL23

In the M5INFALL23 simulation, the ambient gas density is ten times more than the previous simulation. The model starts with the infall of pristine gas at $t_1 = 30$ Myr, and about 9 Myr later is followed by the injection of the AGB ejecta. $t_1$ is shorter than in the previous model because the FG SN feedback is expected to be less effective in pushing the ICM out into the higher density medium. As can be seen in the right panels of Fig. 2, a denser medium makes ram pressure stronger in the corresponding NoRT case and suppresses SF but does not completely stop the accumulation of the ejecta into the cluster. However, we reported in Paper I that a low-mass SG can be formed

in the case of high-density simulation. This different result comes from the fact that in this study, we have assumed the Jeans' length threshold (0.25 pc) that prevents forming stars in cells in which the Jeans' length is larger than four cell width (see Section 2.3).

The left panels of Fig. 2 correspond to the gas density and temperature maps of the RT simulation of M5INFALL23 at the final time. Again, two shocks are visible in the central regions but not close to the left boundary. The region affected by the ionizing heating is less extended than in the previous simulation due to a denser infall. Within this region, the gas density is much lower than in the other parts of the box due to higher pressure and temperature of the gas inside the shock. Therefore, the gas conditions inside the cluster are again not appropriate for forming SG stars within the first 100 Myr from the formation of the FG stars.

Thus, in low-mass clusters, photoionization heating and ram pressure overcome the FG potential, prevent a significant accumulation of cold gas inside the cluster, and result in no SG SF within 100 Myr of the FG formation.

## 3.2 Massive clusters

From theoretical consideration, a $10^6$ M$_\odot$ cluster can overcome the ram pressure and accumulate gas in the NoRT case (Lin & Murray 2007). We numerically confirmed this conclusion in Paper I and showed that not only could the cluster retain its AGB ejecta and accrete the pristine gas, but also that efficient cooling occurred in the central regions, leading to the formation of a massive SG. The more massive the cluster, the more accreted pristine gas, and the more massive the SG cluster and with lower He abundances. We now study the effects of ionizing radiation on SG formation for this cluster mass.

### 3.2.1 M6INFALL24

A massive cluster in a low-density medium is assumed to create the largest SN bubble among our models. Thus, it takes more time for the incoming gas to reach the cluster's centre. Therefore, this simulation starts with the injection of the AGB ejecta at 39 Myrs after the birth of the FG stars, and after about 3.5 Myr the incoming pristine gas passes through the cluster centre.

Fig. 3 shows the density and temperature maps at different times for our simulations with massive clusters. The white contours on the density maps display the SG stars formed in the simulation box. The left-most panel in the top row shows that the AGB ejecta are already accumulating at the centre of the cluster at $t = 42$ Myr, even before the infall reaches the centre. All the gas in the simulation box is immediately exposed to the FG radiation and becomes ionized. As a result, at distances larger than ~10 pc from the centre of the cluster, the photoionization heating by FG stars pushes the ICM off, as can be seen from the arrows in the temperature map. After the cluster enters the dense ISM, all the gas within the box becomes ionized and warm. As the simulation proceeds, a slight gas accumulation is seen in the cluster centre until the end, but it does not result in sufficient cooling for SF. This situation is more or less unchanged until the end of the simulation. The minimum temperature in this run is $3 \times 10^4$ K. Therefore, the suitable conditions for SF are not met in this case. The two upper-right panels show the final results of the NoRT simulation in which an SG with the mass of $\approx 6 \times 10^4$ M$_\odot$ has been formed. Comparing to the final results of the RT with NoRT simulations, the ionizing radiation suppresses the SF at least for 100 Myr, as predicted by Conroy & Spergel (2011).
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Gas density and temperature maps for the M6INFALL24 (first and second rows from top) and M6INFALL23 (third and fourth rows) simulations at different times for RT runs. The right column shows the results of the corresponding NoRT runs at the final time. In the gas density maps, the outer white and inner orange contours show regions in which the SG stellar densities are $10^{-6}$ (i.e. showing the extent of almost all SG stars) and $10^{-3}$ times the maximum stellar density, respectively.</figcaption></figure>

### 3.2.2 M6INFALL23

In a denser medium, the infall reaches the centre at $t = 33.9$, and the cluster accretes more pristine gas. Eventually, efficient cooling occurs only at the centre of the box and stars start forming at about 40 Myr after the formation of the FG. The two bottom-left panels of Fig. 3 display that two shocks formed at $t = 42$ Myr, as seen in the previous simulations. The white contour on the density map shows the first SG stars formed at this time. Moreover, the dense shell due to the ionizing radiation is formed ahead of the cluster at a distance of about 40 pc from the centre, showing that the radiation of FG stars is not sufficient to ionize all the gas in the simulation box, in contrast to M6INFALL24. As the simulation

proceeds, a cold gaseous tail is formed behind the cluster due to the higher gas accretion rate. However, this tail has a lower density than in the simulation without radiation, resulting in no SF on the tail, unlike the NoRT case (see the last column). As a consequence, SG stars are more concentrated in the cluster centre than in the NoRT case.

The top and bottom panels of Fig. 4 show the cumulative SG mass and the star formation rate (SFR) versus time for the RT and NoRT simulations, in which the mass contributions of the pristine gas (dotted line) and AGB ejecta (dashed line) are shown. Moreover, the mass injection rate from the AGB stars is shown with dash-dotted lines. As seen in the top panel, the amount of recycled AGB
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Top: cumulative stellar mass of SG stars formed in the M6INFALL23 simulation and corresponding NoRT run as a function of time, along with the contribution from processed AGB matter and pristine gas. Bottom: SFR of the SG versus time for these runs. The lines are colour-coded as follows: the solid, dashed, and dotted lines show the final SG stellar mass, the mass formed from AGB ejecta, and the mass formed from the pristine gas, respectively. The dash-dotted lines display the rate of stellar mass return injected into the simulation box.</figcaption></figure>

ejecta is the same for both RT and NoRT simulations at the end of the simulations, whereas the contribution of the accreted pristine gas in SF has decreased by almost 50 per cent in the case of RT. Therefore, the inclusion of ionizing radiation leads to a 30 per cent drop in the final SG mass ($10^5 M_\odot$) due to the lower amount of accreted pristine gas. This result is in agreement with previous analytical studies, showing that the accreted mass decreases with higher temperature (Naiman, Ramirez-Ruiz & Lin 2011). The bottom panel shows that the ionizing radiation produces a 2 Myr delay in SF. In the NoRT simulation, the SFR fluctuates around an almost constant value. This rate is larger than the rate of mass injection from AGB stars and can be attributed to the accretion rate of the pristine gas. In the case of the RT run, the SFR experiences a visible drop at ≈50 Myr related to the ionizing feedback of the first SG stars and then follows approximately the NoRT rate until the end of the simulation.

Fig. 5 shows the He distribution of SG stars formed in the simulations at different times. The SG stars with the highest He abundances are formed in the first Myrs of the simulation when the dilution degree of the AGB ejecta is still low. Later the cluster accretes further pristine gas and the distribution eventually tends towards lower helium abundances. Comparing the final NoRT (dashed line) and RT distributions (grey region), SG stars generally have a higher

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Mass distribution of SG stars versus He abundance $Y$ for the M6INFALL23 simulation at different times. The He mass fraction of the SG stars varies between $Y = 0.26$, which is a little higher than the helium mass fraction of the pristine ISM gas (and FG stars, $Y = 0.245$) at the end of the simulation, and $Y = 0.34$, corresponding to the SG stars formed from less-diluted ejecta at the beginning of the simulations.</figcaption></figure>

He mass fraction in the RT simulation because the ionizing radiation limits the pristine gas accretion. While no extremely He-rich ($Y = 0.36$) stars are formed in this run, the maximum He abundance of SG stars is 0.34.

# 4 DISCUSSION

In a previous study (Paper II), we investigated the effects of ionizing feedback on the formation of SG stars in a very massive cluster ($10^7 M_\odot$) in the context of the AGB scenario of D'Ercole, D'Antona & Vesperini (2016). To explore the role of different environments in SG SF, we assumed that the cluster moves through a medium of different densities, i.e. $10^{-24}$ and $10^{-23}$ g cm$^{-3}$. In the low- and high-density simulations, we found that the total SG masses formed within 100 Myr after the FG formation were a fraction of 0.06 and 0.6 of the FG mass, respectively. In this work, we expand the study to FG clusters with the different masses of $10^5$ and $10^6 M_\odot$ and assume a half-mass radius of 4 pc, estimated to be typical for proto-GCs and YMCs (Baumgardt et al. 2019; Krumholz, McKee & Bland-Hawthorn 2019). In both our low- and high-density simulations with the cluster mass $10^5 M_\odot$, the radiative heating dominates the ram pressure and gravitational potential of the cluster. Hence, it does not allow the gas inside the cluster to accumulate at its centre, and as a result, SF is suppressed throughout the entire time of the simulation (100 Myr, Figs 1 and 2). For $10^6 M_\odot$ clusters, our results are density-dependent. In the low-density simulation (M6INFALL24), the FG ionizing radiation heats and expands the gas inside the cluster, resulting in no SG stars formed until the end of the simulation (the first two rows in Fig. 3). On the other hand, in the high-density case (M6INFALL23), a higher accretion rate on the cluster occurs, causing the ICM to cool efficiently at the cluster's centre. Eventually, a more compact SG cluster than in the NoRT case is formed, with a mass of 0.1 of FG mass (Fig. 3).

To explore the role of ionizing feedback in the SG formation, we assume the maximum possible luminosities for our stellar populations, including binary stars (BPASS model; Eldridge et al. 2017). As in Paper II, we assume a bottom-heavy IMF for the SG stars so that the mass of SG stars is truncated at $8 M_\odot$. Hence, Type
<table><caption>**Table 2.** Main results obtained at the end of M6INFALL23 simulation in this study and in Paper II with a $10^7 M_\odot$ cluster. From left to right, the columns show the name of the model, the mass of the FG cluster, the final SG mass, the final SG mass fraction, the fraction of AGB ejecta, pristine gas incorporated in SG stars, and the minimum and maximum He abundance of SG stars.</caption><thead><tr><th>Model</th><th>FG mass (M<sub>☉</sub>)</th><th>SG mass (M<sub>☉</sub>)</th><th>f<sub>SG/FG</sub></th><th>f<sub>AGB</sub></th><th>f<sub>P</sub></th><th>Y<sub>min</sub></th><th>Y<sub>max</sub></th></tr></thead><tbody><tr><td>M6INFALL23</td><td>10<sup>6</sup></td><td>1.0 × 10<sup>5</sup></td><td>0.1</td><td>0.9</td><td>0.1</td><td>0.26</td><td>0.34</td></tr><tr><td>M7INFALL23</td><td>10<sup>7</sup></td><td>64.0 × 10<sup>6</sup></td><td>0.64</td><td>0.1</td><td>0.9</td><td>0.247</td><td>0.36</td></tr><tr><td>M7INFALL24</td><td>10<sup>7</sup></td><td>6.0 × 10<sup>5</sup></td><td>0.06</td><td>0.6</td><td>0.4</td><td>0.345</td><td>0.36</td></tr></tbody></table>

II SNe for SG stars are not included in our simulations, and SG luminosities are not as high as the FG ones.

In Table 2, we report the final properties of the SG population obtained at the end of our RT simulations with cluster masses $10^6$ (this work) and $10^7 M_\odot$ (Paper II). The SG masses in our high-density simulations exhibit an increasing trend with FG cluster mass. The reason is that more massive clusters can accrete more pristine gas, leading to the formation of more massive SG clusters. Assuming that the FG and SG stars have the same IMF at 100 Myr, the SG mass fraction and SG number ratio are equivalent. Therefore, the AGB model can qualitatively explain the observed positive correlation between the SG fraction and cluster mass for the GCs hosting MSPs (see Paper I). However, the observed SG fractions (0.1–0.9) are higher than the ratios computed in our simulations. This difference might come from the fact that our simulations can account only for the first 100 Myr of the evolution, whereas the observed systems have undergone several Gyrs of dynamical evolution. The SG-to-FG ratio is expected to increase during the dynamical evolution (D’Ercole et al. 2008; Vesperini et al. 2021) due to the depletion of FG stars. It is found that most SG stars are retained within clusters because they are located in the central regions. In Paper I, we showed that assuming factors between 5 and 20¹ for the ratio between the initial mass of the cluster and the present-day mass, allows us to account for the observed correlation between SG fraction and cluster mass. A similar conclusion can be obtained for our RT simulations of different masses. The dilution of the AGB ejecta in the RT-M6INFALL23 simulation results in the formation of SG stars with He mass fractions of $0.26 < Y < 0.34$. The contribution of pristine gas in forming the SG is about 10 per cent of the FG mass, which is consistent with the value proposed by D’Ercole, D’Antona & Vesperini (2011). In fact, based on the observed Na–O anticorrelation, they show that the amount of pristine gas cannot be larger than ~$0.1 M_{FG}$. On the other hand, the $10^7 M_\odot$ clusters can form SG stars with extreme-He mass fractions (≈0.36, see last two columns of Table 2), showing that some of the SG stars originate out of pure AGB ejecta. It is worth noting that the dilution process of the enriched ejecta has a longer delay in massive clusters, due to the formation of a larger hot bubble by FG Type II SNe around the cluster. This causes some SG stars to form from enriched materials only before the dilution of stellar winds. This result is consistent with what was proposed by

D’Ercole, D’Antona & Vesperini (2016) to explain the existence of extreme-He fraction in SG stars in massive GCs. From our high-density simulations with cluster masses of $10^6$ and $10^7 M_\odot$, the maximum He enhancement shows an increasing trend with cluster mass. In principle, this increasing trend can explain the observed positive correlation between maximum He enhancement and cluster mass (Milone et al. 2020). The capability of reproducing this trend is a promising result for the AGB scenario. Moreover, we found in Paper I that in NoRT simulations corresponding to M5INFALL24, M5INFALL23, and M6INFALL24, SG stars with extreme-He fractions could be formed, which is not in agreement with observations. This study shows that the radiative feedback suppresses the SF in these simulations and prevents the formation of these extreme-He stars. Moreover, it suppresses SF in low-mass clusters, regardless of the density of the external gas. Therefore, our results indicate that the ionizing stellar feedback is not a severe problem for SG formation. Rather, it can help the AGB scenario to account for some observables, such as the lack of extreme He-rich stars in lowest mass GCs. However, other issues, such as the mass budget problem and the assumption of a bottom-heavy IMF for SG stars, still need further investigation.

Our study confirms that only massive GCs can form a new generation from their AGB ejecta and accreted pristine gas, as already proposed for this model (D’Ercole et al. 2008; D’Ercole, D’Antona & Vesperini 2016). Moreover, we introduce another requirement for this scenario; these massive clusters should be located in dense environments to accrete sufficient pristine gas, in order to fulfill two requirements. First, to provide a significant gas accumulation within clusters, for the efficient cooling required for SF. Secondly, to supply the pristine gas required for appropriate dilution of the AGB ejecta, so that the chemical compositions of SG stars can be comparable to observations. Therefore, our simulations with photoionization heating confirm the key role of pristine gas accretion in SG formation in diluting the AGB ejecta and solve the abundance problem. In addition, we confirm its undeniable role in efficient cooling within young clusters.

The density of $10^{-23} g cm^{-3}$, used in our high-density simulations, is found in star-forming disc galaxies observed at redshift $z > 2$, in agreement with the D’Ercole, D’Antona & Vesperini (2016) assumption that the formation of SG stars is limited to these redshifts. On the other hand, this point can explain why no strong evidence has been found so far for the presence of MSPs in YMCs (Mucciarelli et al. 2008, 2014; Martocchia et al. 2018). Probably, they could not be able to accumulate a considerable amount of pristine gas at their centre, because they are located in lower density environments than the ones assumed here for the formation of MSPs, and their ionizing feedback becomes effective in heating the ICM. As a result, the necessary conditions for SG SF may not be satisfied in these clusters. This can also explain why YMCs do not show evidence of ongoing SF and why they are gas free (Bastian et al. 2013b). We can check the gas content within young clusters in our simulations with

<footnote n="1">¹As for the factor 5, Larsen, Strader & Brodie (2012) showed that the mass of the field stars of the Fornax dwarf galaxy can at most be 4–5 times larger than that of his GCs. In case the field stars have all been lost from the GCs, this result imposes an upper limit on their initial mass. However, an important quantity currently unknown is the mass of Fornax, which in the past could have been greater than today. The factor 20 comes from the argument that it is possible to show that approximately 5 per cent (i.e. 1/20) of the FG mass comes out with the right composition, in terms of He abundance and p-process elements, to produce SG stars (Renzini 2013), assuming a 100 per cent SF efficiency.</footnote>
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Gas mass within different radii from 1 to 5 pc in the clusters studied in this paper as a function of time.</figcaption></figure>

the caveat that the assumed metallicity in our model is lower than in local galaxies, with possible effects on radiative cooling.

Fig. 6 shows the evolution of the ICM mass in the radii of 1 (dash-dotted line), 2 (thin solid line), and 5 pc (thick solid line) from the cluster's centre for the simulations of this study. To better follow how gas mass inside the cluster depends on the radius, the coloured shaded regions indicate the areas covered by lines corresponding to each simulation. The gas mass quickly increases after the infall passes through the cluster centre and then settles on to an approximately constant value. This translates into an equilibrium state between the FG gravity, ram pressure, and FG ionizing feedback. Bastian et al. (2013b) estimate the mass of ionized gas within a sphere of radius 1 pc in YMCs according to their H β flux. They find that the gas mass within YMCs should be lower than $100 M_{\odot}$. As can be seen in Fig. 6, the gas mass within 1 pc for those simulations in which no SG stars form (M5INFALL24, M5INFALL23, and M6INFALL24), including both low-density simulations, is much lower than this value (grey dashed line), while it is larger in the M6INFALL23 simulation. Therefore, our study suggests that radiative feedback can prevent the accumulation of gas and SF inside clusters of masses $\le 10^6 M_{\odot}$ in local galaxies. However, further simulations specific to YMCs should be performed to confirm this result. This topic represents an interesting subject for future work.

## 5 CONCLUSIONS

Using hydrodynamical simulations, we have investigated the effects of photoionization feedback on SF and gas content in young clusters with masses of $10^5$ and $10^6 M_{\odot}$ within 30–100 Myr from the birth of clusters. We assume that the clusters move in a uniform ISM with gas densities of $10^{-24}$ and $10^{-23}$ g cm$^{-3}$ to explore the effect of different environments on the amount of gas accretion on to clusters and SG formation. These simulations are designed to check the predictions of the AGB scenario in the formation of a new generation from the stellar winds of AGB stars and accreted pristine gas in the presence of photoionization heating. Our main results can be summarized as follows:

(i) In the simulations with a $10^5 M_{\odot}$ cluster, the radiative heating dominates the ram pressure and gravitational potential of the cluster. Then, it does not allow the gas inside the cluster to accumulate at its centre and create any low-mass stellar aggregate. As a result, SF is suppressed until the end of simulations.

(ii) For $10^6 M_{\odot}$ clusters, we find that our results are density-dependent. In the low-density simulation (M6INFALL24), the FG ionizing radiation heats and expands the gas inside the cluster, resulting in no SG stars formed until $t = 100$ Myr. On the other hand, in the high-density case (M6INFALL23), a higher accretion rate on the cluster occurs, causing the ICM to cool efficiently at the cluster's centre. Eventually, a more compact SG cluster than in the NoRT case is formed, with a mass of 0.1 of FG mass. Due to ionizing effects, the cluster accretes less pristine gas, and therefore the new generation has a lower mass (30 per cent) and higher He abundances in comparison to the simulation without radiation.

(iii) We confirm the key role of pristine gas accretion in SG formation and introduce another requirement for the AGB scenario: massive clusters should be located in dense environments to accrete sufficient pristine gas to fulfill two goals. First, to provide a significant gas accumulation within clusters for the efficient cooling required for SF. Secondly, to supply the pristine gas required to dilute the AGB ejecta and produce an increasing trend between maximum Y variation and mass.

(iv) From the results of our study, it is apparent that the ionizing stellar feedback is not a severe problem for SG formation. Moreover, it can also help the AGB scenario to account for some observables, such as the lack of extreme He-rich stars in lowest mass GCs. Moreover, the SG masses and maximum He enhancement in high-density simulations exhibit increasing trends with FG cluster mass.

(v) In our low-density simulations, the gas content is too low so that it can be compared with estimated values for YMCs. In addition, these simulations do not show SF, as seen in YMCs. Hence, our high-density simulations can better be matched to observations corresponding to GCs hosting MSPs, while our low-density ones can quantitatively explain the properties of young clusters in the local Universe. Therefore, the AGB model can provide a reason why YMCs do not show evidence of SG formation and any significant ICM.

## ACKNOWLEDGEMENTS

We acknowledge support and computational resources from the PSMN (Pôle Scientifique de Modélisation Numérique) of the ENS de Lyon. FC acknowledges support from PRIN INAF 1.05.01.85.01, the National Recovery and Resilience Plan (Piano Nazionale di Ripresa e Resilienza, PNRR) Project ID CN 00000013 ‘Italian Research Center on High-Performance Computing, Big Data and Quantum Computing’ funded by MIUR Missione 4 Componente 2 Investimento 1.4: Potenziamento strutture di ricerca e creazione di ‘campioni nazionali di R&S (M4C2-19)’ – Next Generation EU (NGEU). The research activities described in this paper have been co-funded by the European Union – NextGenerationEU within PRIN 2022 project no. 20229YBSAN – Globular clusters in cosmological simulations and in lensed fields: from their birth to the present epoch. AY acknowledges the support of Ferdowsi University of Mashhad.

## DATA AVAILABILITY

The data supporting this study's findings are available from the corresponding author upon reasonable request.

## REFERENCES

Bastian N., Lardo C., 2018, ARA&A, 56, 83

Bastian N., Strader J., 2014, MNRAS, 443, 3594
Bastian N., Lamers H. J. G. L. M., de Mink S. E., Longmore S. N., Goodwin S. P., Gieles M., 2013a, MNRAS, 436, 2398

Bastian N., Cabrera-Ziri I., Davies B., Larsen S. S., 2013b, MNRAS, 436, 2852

Bastian N., Hollyhead K., Cabrera-Ziri I., 2014, MNRAS, 445, 378

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

Bekki K., 2017, MNRAS, 467, 1857

Cabrera-Ziri I. et al., 2015, MNRAS, 448, 2224

Cadelano M., Dalessandro E., Webb J. J., Vesperini E., Lattanzio D., Beccari G., Gomez M., Monaco L., 2020, MNRAS, 499, 2390

Calura F., Few C. G., Romano D., D'Ercole A., 2015, ApJ, 814, L14

Calura F., D'Ercole A., Vesperini E., Vanzella E., Sollima A., 2019, MNRAS, 489, 3269

Carretta E. et al., 2009, A&A, 505, 117

Chantereau W., Biernacki P., Martig M., Bastian N., Salaris M., Teyssier R., 2020, MNRAS, 493, 1306

Conroy C., Spergel D. N., 2011, ApJ, 726, 36

D'Ercole A., Vesperini E., D'Antona F., McMillan S. L. W., Recchi S., 2008, MNRAS, 391, 825

D'Ercole A., D'Antona F., Vesperini E., 2011, MNRAS, 415, 1304

D'Ercole A., D'Antona F., Vesperini E., 2016, MNRAS, 461, 4088

Decressin T., Meynet G., Charbonnel C., Prantzos N., Ekström S., 2007a, A&A, 464, 1029

Decressin T., Charbonnel C., Meynet G., 2007b, A&A, 475, 859

Denissenkov P. A., Hartwick F. D. A., 2014, MNRAS, 437, L21

Denissenkov P. A., VandenBerg D. A., Hartwick F. D. A., Herwig F., Weiss A., Paxton B., 2015, MNRAS, 448, 3314

Dickens R. J., Croke B. F. W., Cannon R. D., Bell R. A., 1991, Nature, 351, 212

Eldridge J. J., Stanway E. R., Xiao L., McClelland L. A. S., Taylor G., Ng M., Greis S. M. L., Bray J. C., 2017, Publ. Astron. Soc. Aust., 34, e058

Gavagnin E., Bleuler A., Rosdahl J., Teyssier R., 2017, MNRAS, 472, 4155

Gratton R., Sneden C., Carretta E., 2004, ARA&A, 42, 385

Gratton R. G., Carretta E., Bragaglia A., 2012, A&AR, 20, 50

Hollyhead K., Bastian N., Adamo A., Silva-Villa E., Dale J., Ryon J. E., Gazak Z., 2015, MNRAS, 449, 1106

Johnson C. I., Pilachowski C. A., 2010, ApJ, 722, 1373

Khalaj P., Baumgardt H., 2015, MNRAS, 452, 924

Krause M., Charbonnel C., Decressin T., Meynet G., Prantzos N., 2013, A&A, 552, A121

Krause M. G. H., Charbonnel C., Bastian N., Diehl R., 2016, A&A, 587, A53

Kravtsov A. V., Gnedin O. Y., 2005, ApJ, 623, 650

Kroupa P., 2001, MNRAS, 322, 231

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

Krumholz M. R., McKee C. F., Bland-Hawthorn J., 2019, ARA&A, 57, 227

Lacchin E., Calura F., Vesperini E., 2021, MNRAS, 506, 5951

Lacchin E., Calura F., Vesperini E., Mastrobuono-Battisti A., 2022, MNRAS, 517, 1171

Larsen S. S., Strader J., Brodie J. P., 2012, A&A, 544, L14

Larson R. B., 1996, in Morrison H. L., Sarajedini A., eds, ASP Conf. Ser. Vol. 92, Formation of the Galactic Halo...Inside and Out. Astron. Soc. Pac., San Francisco, p. 241

Lin D. N. C., Murray S. D., 2007, ApJ, 661, 779

Marcolini A., Brighenti F., D'Ercole A., 2003, MNRAS, 345, 1329

Martocchia S. et al., 2018, MNRAS, 473, 2688

Milone A. P., Marino A. F., 2022, Universe, 8, 359

Milone A. P. et al., 2017, MNRAS, 464, 3636

Milone A. P. et al., 2020, MNRAS, 491, 515

Mucciarelli A., Carretta E., Origlia L., Ferraro F. R., 2008, AJ, 136, 375

Mucciarelli A., Dalessandro E., Ferraro F. R., Origlia L., Lanzoni B., 2014, ApJ, 793, L6

Naiman J. P., Ramirez-Ruiz E., Lin D. N. C., 2011, ApJ, 735, 25

Naiman J. P., Ramirez-Ruiz E., Lin D. N. C., 2018, MNRAS, 478, 2794

Piotto G. et al., 2015, AJ, 149, 91

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

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

Prantzos N., Charbonnel C., Iliadis C., 2007, A&A, 470, 179

Rasera Y., Teyssier R., 2006, A&A, 445, 1

Renzini A., 2013, Mem. Soc. Astron. Ital., 84, 162

Renzini A. et al., 2015, MNRAS, 454, 4197

Rosdahl J., Teyssier R., 2015, MNRAS, 449, 4380

Rosdahl J., Blaizot J., Aubert D., Stranex T., Teyssier R., 2013, MNRAS, 436, 2188

Schmidt M., 1959, ApJ, 129, 243

Teyssier R., 2002, A&A, 385, 337

Ventura P., D'Antona F., 2011, MNRAS, 410, 2760

Vesperini E., Hong J., Giersz M., Hypki A., 2021, MNRAS, 502, 4290

Wardlow J. L. et al., 2017, ApJ, 837, 12

Wünsch R., Palouš J., Tenorio-Tagle G., Ehlerová S., 2017, ApJ, 835, 60

Yaghoobi A., Calura F., Rosdahl J., Haghi H., 2022a, MNRAS, 510, 4330

Yaghoobi A., Rosdahl J., Calura F., Khalaj P., Haghi H., 2022b, MNRAS, 517, 4175

Yong D., Grundahl F., Norris J. E., 2015, MNRAS, 446, 3319

This paper has been typeset from a TeX/LATEX file prepared by the author.