We present the first comprehensive analysis of the internal dynamics of multiple stellar populations MPs in 28 Galactic Globular Clusters GCs across a wide field of view extending from the innermost regions to the clusters outskirts. Using astrophotometric catalogs from groundbased observations Gaia and the Hubble Space Telescope HST we identify first 1P and secondpopulation 2P stars and study the internal dynamics of MPs using highprecision Gaia DR3 and HST proper motions. Our results reveal that while the 1P transitions from isotropy to slight tangential anisotropy toward the outer regions 2P stars become increasingly radially anisotropic beyond the halflight radius. We also explore the connection between the dynamics of MPs and the clusters structural and dynamical properties finding statistically significant differences in the anisotropy profiles of dynamically young and nonrelaxed clusters particularly beyond the 12 halflight radii. In these regions 1P stars transition from isotropic to slightly tangentially anisotropic motion while 2P stars become more radially anisotropic. In contrast dynamically older clusters with mixed MPs exhibit weaker relative differences. Furthermore clusters with orbits closer to the Galactic center exhibit larger dynamical differences between 1P and 2P stars than those with larger periGalactic radii. These findings are consistent with a scenario where 2P stars form in a more centrally concentrated environment where the interaction with the Milky Way tidal field plays a crucial role in the dynamical evolution of MPs especially of 1P.
['Astrophysics - Astrophysics of Galaxies', 'Astrophysics - Solar and Stellar Astrophysics']
Internal Dynamics of Multiple Populations in 28 Galactic Globular Clusters A WideField study with Gaia and the Hubble Space Telescope
{'Internal Dynamics of Multiple Populations in 28 Galactic Globular Clusters: A Wide-Field study with Gaia and the Hubble Space Telescope': "- G. Cordoni 1 , 2 , ★ L. Casagrande 1 , 2 , A. P. Milone 3 , 4 , E. Dondoglio 3 , A. Mastrobuono-Battisti 5 , S. Jang 6 ,\n- A. F. Marino 4 , 7 , E. P. Lagioia 8 , M. V. Legnardi 3 , T. Ziliotto 3 , F. Muratore 3 , V. Mehta 1 ,\n- E. Lacchin 3 , 9 , 10 , 11  \n, M. Tailo 3  \n- 1 Research School of Astronomy and Astrophysics, The Australian National University, Canberra, ACT 2611, Australia\n- 2 ARC Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO 3D), Australia\n- 3 Dipartimento di Fisica e Astronomia 'Galileo Galilei' - Univ. di Padova, Vicolo dell'Osservatorio 3, Padova, IT-35122\n- 4 Istituto Nazionale di Astrofisica - Osservatorio Astronomico di Padova, Vicolo dell'Osservatorio 5, Padova, IT-35122\n- 5 GEPI, Observatoire de Paris, PSL Research University, CNRS, Place Jules Janssen, 92195 Meudon, France\n- 6 Center for Galaxy Evolution Research and Department of Astronomy, Yonsei University, Seoul 03722, Korea\n- 7 Istituto Nazionale di Astrofisica - Osservatorio Astrofisico di Arcetri, Largo Enrico Fermi, 5, Firenze, IT-50125\n- 8 South-Western Institute for Astronomy Research Yunnan University, Kunming, 650500, P.R. China\n- 9 Osservatorio di Astrofisica e Scienza dello Spazio di Bologna, Via Gobetti 93/3, 40129 Bologna, Italy\n- 10 Institut für Theoretische Astrophysik, ZAH, Universität Heidelberg, Albert-Ueberle-Straße 2, D-69120, Heidelberg, Germany\n- 11 INFN - Padova, Via Marzolo 8, I-35131 Padova, Italy  \nAccepted XXX. Received YYY; in original form ZZZ", 'ABSTRACT': "We present the first comprehensive analysis of the internal dynamics of multiple stellar populations (MPs) in 28 Galactic Globular Clusters (GCs) across a wide field of view, extending from the innermost regions to the clusters' outskirts. Using astrophotometric catalogs from ground-based observations, Gaia, and the Hubble Space Telescope (HST), we identify first- (1P) and second-population (2P) stars, and study the internal dynamics of MPs using high-precision Gaia DR3 and HST proper motions. Our results reveal that while the 1P transitions from isotropy to slight tangential anisotropy toward the outer regions, 2P stars become increasingly radially anisotropic beyond the half-light radius. We also explore the connection between the dynamics of MPs and the clusters' structural and dynamical properties, finding statistically significant differences in the anisotropy profiles of dynamically young and non-relaxed clusters, particularly beyond the 1-2 half-light radii. In these regions, 1P stars transition from isotropic to slightly tangentially anisotropic motion, while 2P stars become more radially anisotropic. In contrast, dynamically older clusters, with mixed MPs, exhibit weaker relative differences. Furthermore, clusters with orbits closer to the Galactic center exhibit larger dynamical differences between 1P and 2P stars than those with larger peri-Galactic radii. These findings are consistent with a scenario where 2P stars form in a more centrally concentrated environment, where the interaction with the Milky Way tidal field plays a crucial role in the dynamical evolution of MPs, especially of 1P.  \nKey words: stars: Hertzsprung-Russell and colour-magnitude diagrams; stars: kinematics and dynamics; Galaxy: globular clusters: general; Galaxy: kinematics and dynamics", '1 INTRODUCTION': "More than two decades of Hubble Space Telescope (HST) observations have revealed that nearly all Galactic Globular clusters (GCs) host Multiple Populations (MPs) with at least two main groups of stars: first-population (1P) and second-population (2P Gratton et al. 2012; Piotto et al. 2015; Renzini et al. 2015; Milone et al. 2017; Bastian & Lardo 2018). These groups exhibit specific chemical compositions, with 2P stars being for example enriched in Sodium, Nitrogen and Helium and depleted in Carbon and Oxygen (see e.g. Bastian & Lardo 2018; Milone & Marino 2022, for recent reviews). Despite several models have been proposed to explain the formation of  \nMPs, none can simultaneous match the numerous observational constraints (see e.g. Renzini et al. 2015; Bastian & Lardo 2018; Milone & Marino 2022). Some scenarios predict that 2P stars formed from the ejecta of short-lived massive 1P stars (Ventura et al. 2001; Decressin et al. 2007; D'Ercole et al. 2010; Denissenkov & Hartwick 2014). As the expelled gas accumulated toward the cluster center, 2P stars would dominate the composition of the central cluster regions. Alternatively, the chemical composition of 2G stars could be the result of the accretion of chemically enriched material emitted from 1P stars in binary configuration (de Mink et al. 2009; Bastian et al. 2013; Gieles et al. 2018).  \nWhile some clusters indeed exhibit a more centrally concentrated 2P, there are numerous examples where 1P and 2P stars are spatially mixed, and even clusters with a more centrally concentrated 1P are  \nobserved (see e.g. Dalessandro et al. 2019; Leitinger et al. 2023). Recent studies suggest that the extent of this mixing between 1P and 2P stars may correlate with the dynamical age of the cluster, indicating that dynamically older clusters tend to show a more complete mixing of MPs. Nonetheless, despite the numerous observational constraints, a comprehensive picture of how MPs interact form and evolve over time remains elusive.  \nIn recent years, a promising new path has emerged: the study of the internal kinematics of cluster stars. This approach is crucial for uncovering the physical processes behind the formation of multiple populations. Specifically, N -body simulations indicate that the dynamical evolution of 2P stars formed in different environment should differ significantly from that of more spatially extended 1P stars. Such differences may still be detectable in present-day GC kinematics (Mastrobuono-Battisti & Perets 2013, 2016; Vesperini et al. 2013; Vesperini et al. 2021; Hénault-Brunet et al. 2015; Tiongco et al. 2019; Lacchin et al. 2022), provided the populations are not completely relaxed.  \nOver the past decade, a number of studies have analyzed the spatial distribution and the internal kinematics of MPs in various GCs, using HST photometry and proper motions (see e.g. Libralato et al. 2023; Leitinger et al. 2023), MUSE/APOGEE line-of-sight (LoS) velocities (Szigeti et al. 2021; Martens et al. 2023), ground-based photometry coupled with Gaia proper motions, or a combination of these methods (Milone et al. 2018; Cordoni et al. 2020a; Cordoni et al. 2020b; Cordoni et al. 2023; Cadelano et al. 2024). Such works presented the first evidence of significant dynamical differences between 1P and 2P stars in some GCs, qualitatively aligning with 𝑁 -body simulations' results. However, the limited number of clusters and the small field of view prevented a complete and detailed characterization of the phenomenon.  \nIn this work we aim to overcome the limitation of previous studies by investigating the internal dynamics of MPs in a sample of 28 GCs, combining photometric information from HST (Nardiello et al. 2018), ground-based (see e.g. Stetson et al. 2019; Jang et al. 2022) and Gaia (Mehta et al. 2024) and accurate HST (Libralato et al. 2022) and Gaia (Gaia Collaboration et al. 2023a) proper motions. This work presents the first comprehensive dynamical investigation of MPs in a large sample of clusters, extending from the central arcminutes to nearly the tidal radius. The paper is structured as follows: in Sec. 2 we introduce the dataset and MPs selection process. Section 3 presents the study of the internal dynamics in the analyzed clusters and the global dynamical profiles of MPs. Finally, the discussion and summary are presented in Sec. 4 and 5, respectively.", '2 DATA': 'To investigate the internal dynamics of MPs in Galactic GCs over a wide field of view, we combined HST photometry (Milone et al. 2017), ground-based UBVI photometry (Stetson et al. 2019; Jang et al. 2022), Gaia XP spectro-photometry (Mehta et al. 2024) 1 , with HST(Libralato et al. 2022) and Gaia DR3 proper motions (Gaia Collaboration et al. 2023a). We refer to the different works for a detailed description of the astro-photometric datasets. The total sample used in this work consists of 28 GCs for which we could reliably select MPs stars over a wide field of view.  \nTo identify 1P and 2P among RGB stars, we have made use of  \nthe Chromosome Maps (ChMs) presented in Milone et al. (2017) for HST data, photometric diagrams made with the 𝐶 UBI index 2 , (Cordoni et al. 2020a; Jang et al. 2022) for ground-based observations, and the ChMs introduced in Mehta et al. (2024) from Gaia spectro-photometry. ChMs effectively and accurately separate MPs in GCs, with 1P stars defining a clump of stars centered around (0, 0), and 2P stars spread over larger absolute values of Δ BI and Δ CUBI ( Δ F275W , F814W and Δ CF275W , F336W , F438W in HST filters). In the next section, we discuss in detail the procedure used to identify 1P and 2P stars in the ground-based ChMs of 28 GCs, whereas we refer to Milone et al. (2017) and Mehta et al. (2024) for the identification of MPs in HST and Gaia XP photometry.  \nThe final selection of 1P and 2P stars for each cluster and the analyzed field of view are shown in App. A as online supplementary material. 3', '2.1 Multiple Populations from ground-based photometry': 'In a nutshell, the 𝐼 vs. 𝐵 -𝐼 and 𝐼 vs. 𝐶 UBI = 𝑈 -2 𝐵 + 𝐼 CMDsofRed Giant Branch stars (RGB) have been used to derive the verticalized Δ BI and Δ CUBI colors, that are respectively shown on the 𝑥 and 𝑦 axis of the ChMs. We refer to Jang et al. (2022) for a detailed discussion on how these quantities have been derived.  \nTo select 1P and 2P stars from the ChMs we adapted the procedure introduced in Milone et al. (2017) on HST ChMs to our ground-based ChMs, as done in Jang et al., in preparation. Figure 1 illustrates the procedure for the cluster NGC 2808. Briefly, we first selected by eye the bulk of 1P stars among the stars with ( Δ BI , Δ CUBI ) ∼ ( 0 , 0 ) , and fitted a straight line to 1P stars (black solid line in Fig. 1, leftmost panel). We then rotated the ChMs by the angle 𝜃 defined by slope of the 1P best-fit line. We refer to the new rotated reference frame as Δ x and Δ y (Fig. 1, middle panel). By construction, 1P stars define a horizontal distribution in Δ y vs Δ x, and we can use the Δ y distribution to disentangle 1P and 2P stars. We used the scikit-learn python package (Pedregosa et al. 2011) and employed Gaussian Mixture Models (GMM) to fit 2 Gaussians to the Δ y distribution. We remind here that GMM is applied directly to the data points, so that the results do not depend on the bin choice. Finally, 1P stars (red points in the rightmost panel) have been selected as stars with Δ y < Δ sep y , with Δ sep y being the intersection between the two best-fit Gaussians, indicated by the horizontal line in the middle panel of Fig. 1. 2P stars are marked with blue points. Additionally, stars with unusual values of Δ BI and Δ CUBI have been excluded from the selection.  \nThe ground-based ChMs with selected 1P and 2P stars of the 28 analyzed Galactic GCs are shown in App. A as online supplementary material.  \nAdditionally, we employed the 𝑉 vs. 𝐶 UBI pseudo-CMD to extend the regions analyzed in Jang et al. (2022), including stars at larger distance from the cluster centers. We refer to Cordoni et al. (2020a) for a detailed description of the selection of MPs.  \nFigure 1. 1P and 2P star selection process for the GC NGC 2808 using ground-based ChM. Left panel. Ground-based ChM from Jang et al. (2022). Middle panel. ChM rotated by an angle 𝜃 (as shown in the left panel) with a GMM fit applied to Δ y. Right panel. The same ChM as in the left panel, with 1P and 2P stars highlighted in red and blue, respectively.  \n<!-- image -->', '2.2 Multiple populations from HST photometry': 'To study MPs among RGB stars in the innermost cluster regions we exploited the ChMs and selection carried out in Milone et al. (2017). In a nutshell, HST ChMs have been derived by means of the verticalized 𝑚 F814W vs. 𝑚 F275W -𝑚 F814W (for the abscissa) and 𝑚 F814W vs. 𝐶 F275W , F336W , F438W 4 CMDs (on the ordinate). 1P and 2P stars have been selected with the procedures detailed in Milone et al. (2017, see e.g. Sec. 3.3 for a detailed description). Each HST ChM has then been cross matched with the proper motion catalogs of Libralato et al. (2022).', '2.3 Multiple populations from Gaia XP synthetic photometry': "The recent Gaia DR3, (Gaia Collaboration et al. 2023a) publicly released low-resolution spectra, namely XP spectra, for ∼ 200 million sources (De Angeli et al. 2023), allowing to derive synthethic photometry in virtually any photometric systems (see e.g. Gaia Collaboration et al. 2023b). In a recent work, Mehta et al. (2024) exploited Gaia XP spectra to compute synthetic photometry special filters designed to maximize the separation between 1P and 2P stars in GCs' RGB. Their analyzed field of view extends well beyond that covered by the ground-based photometric dataset of Stetson et al. (2019), reaching their tidal radius. We refer to Mehta et al. (2024) for a detailed description of the procedure used to derive and separate the populations, and for the validation of such identification. In a nutshell, exploiting synthethic spectra with the typical composition of 1P and 2P stars, Mehta and collaborators defined a series of photometric filters to maximize the separation between 1P and 2P stars in 5 Galactic GCs. To validate their MPs identification, Mehta et al. (2024) compared the classification with available spectroscopic information, finding consistent results. In this work, we exploit their identification to extend our ground-based dataset and investigate the dynamics of 1P and 2P stars up to the clusters' tidal radii.", '3 INTERNAL DYNAMICS OF MULTIPLE POPULATIONS': "After carefully selecting MPs in the 28 analyzed clusters, we exploited HST (Libralato et al. 2022) and Gaia DR3 (Gaia Collaboration et al. 2023a) proper motions to investigate the internal dynamics of stars belonging to different populations. Specifically, we  \nfollowed the procedure described in Libralato et al. (2022); Libralato et al. (2023, for HST) and Bianchini et al. (2018) 5 . We first transformed the celestial coordinates and their proper motion components ( 𝛼, 𝛿, 𝜇 ∗ 𝛼 , 𝜇 𝛿 ) 6 into a Cartesian reference frame ( 𝑥, 𝑦, 𝜇 x , 𝜇 y ) , using the orthographic projection (eq. 2 in Gaia Collaboration et al. 2018), and then the proper motions into their sky-projected radial and tangential components, defined as 𝜇 R = 𝑥𝜇 x + 𝑦𝜇 y and 𝜇 T = 𝑦𝜇 x -𝑥𝜇 y, with 𝜇 R pointing outward (positive) or inward (negative), and 𝜇 T positive when counterclockwise. 𝜇 x and 𝜇 x here indicate the proper motion of each star relative to the motion of the cluster, determined by Vasiliev & Baumgardt (2021). Furthermore, considering that the sample of analyzed stars is different from that analyzed in Vasiliev & Baumgardt (2021), we repeated the analysis determining the mean cluster motion from our sample of stars, and accounting for systematic errors as in Vasiliev & Baumgardt (2021) 7 , finding consistent results.  \nThe uncertainties on the radial and tangential components of the proper motions have been determined from the uncertainties on 𝜇 ∗ 𝛼 and 𝜇 𝛿 , accounting for the correlation between the proper motions. Finally, the radial proper motions have been corrected for the perspective expansion/contraction due the bulk motion of the cluster along the line of sight, as in Bianchini et al. (2018).  \nThe mean motion and velocity dispersion of 1P and 2P stars have then been computed by dividing the field of view into different concentric annuli containing approximately the same number of stars. For each annulus we determined the mean motion ( 𝜇 R / T ) and dispersion ( 𝜎 R / T ) from Gaia DR3 radial and tangential components by minimizing the negative log-likelihood defined in Bianchini et al. (2018) including the covariance term as in Sollima et al. (2019). Uncertainties on the determination of the mean motion and dispersion have been determined using Markov Chain Monte Carlo algorithm ( emcee Foreman-Mackey et al. 2013).  \nThe dynamical profiles have been re-computed using different number of bins and bins with the same radial width to determine the effect of arbitrary bin choice, with consistent results. Concerning the HST proper motions, we used the same likelihood without the covariance term, (see e.g. Equation 1 in Libralato et al. 2022) and  \n6  \n7  \n∗  \n𝛼  \n=  \n𝛼  \nroutines  \navailable  \nat  \nhttps://github.com/  \nGalacticDynamics-Oxford/GaiaTools .  \n𝜇  \n𝜇  \ncos  \n𝛿  \nmaximizing only for the radial/tangential dispersion. This is because HST proper motions are relative to the cluster's bulk motion and do not provide individual stellar absolute mean motion.  \nIn the following sections, we present and discuss the results for selected individual clusters as well as the global dynamical profiles derived by combining data from all clusters. This approach enables us to explore the general dynamical properties of 1P and 2P stars. Our analysis focuses on the dispersion profiles, while we defer a detailed investigation of the rotation of 1P and 2P stars to a forthcoming paper that is currently in preparation.", '3.1 Individual dynamical profiles': 'In this section, we present the study of the dynamical profiles of 1P and 2P stars in a few individual clusters with sufficiently large samples, namely NGC 0104, NGC 2808, NGC 5904, and NGC 6205. The dynamical profile of the remaining clusters are shown in App. B as online supplementary material.  \nFigure 2 shows, from top to bottom, the mean radial and tangential profiles, the radial and tangential dispersion, and the anisotropy profiles, defined as 𝛽 = 𝜎 T / 𝜎 R -1. The red and blue solid lines indicate the dynamical quantities, determined as described in Sec. 3, while the shaded rectangles indicate the relative uncertainties and the extension of each radial bin. The distance from the cluster center are in units of 𝑅 h indicated in the middle panel (Harris 1996, revision of 2010). Black diamonds indicate the values of the 1-dimensional dispersion from Libralato et al. (2022), while the regions with available HST data is indicated by the vertical green lines. The top and right axis (shown in gray) of each panel display distances in units of parsec and km/s, converted considering cluster distances determined in Baumgardt & Hilker (2018).  \nA visual inspection of the individual profiles reveals that, as expected, the mean radial motions of cluster stars are consistent with 0, as we already subtracted the perspective contraction/expansion. Consistently with Vasiliev & Baumgardt (2021), we find that all four clusters exhibit rotation in the plane of the sky. The mean tangential profiles of NGC 0104 and NGC 5904 exhibit clear negative and positive values, respectively indicating counterclockwise and clockwise rotation. Possible differences in the rotation profile are consistent with the inferred uncertainties. NGC 2808 and NGC 6205, show possible hints of stellar rotation, even though the signal too weak to carry out a comparison.  \nConcerning the dispersion profiles, 1P stars exhibit a larger tangential and similar radial dispersion compared with 2P stars, especially outside ∼ 1 -2 𝑅 h , while 1P and 2P have similar radial dispersion profiles. As a result, 1P stars exhibit an isotropic motion in the central regions, possibly shifting toward tangentially anisotropic motion in the outer regions. On the other hand, the motion of 2P stars is consistent with being isotropic in the cluster center, becoming radially anisotropic beyond 1 𝑅 h . A similar pattern is also observed in NGC3201, consistently with Cadelano et al. (2024). The mean, dispersion and anisotropy profiles of the remaining clusters are displayed in App. B. Overall, we find a good agreement between HST and Gaia dynamical profiles in overlapping or neighbouring regions.  \nThe observed dynamical profiles of 1P and 2P stars are qualitatively consistent with the predictions from 𝑁 -body simulations where 2P stars form in a more centrally concentrated environment. We refer to Tiongco et al. (2019); Vesperini et al. (2021); Pavlík & Vesperini (2022); Pavlík et al. (2024) for a detailed description of the assumptions adopted in the simulations and the results. Additional discussion is carried out in Sec. 4.', '3.2 Global dynamical profiles': "Given the low number of 1P and 2P stars in many clusters, a detailed analysis of the individual dynamical profiles as a function of cluster radius is often unfeasible. Therefore, similarly to Libralato et al. (2023), we opted to investigate the internal dynamics in all clusters simultaneously, combining together the dispersion profiles of the individual clusters 8 . In order to properly compare different clusters, we normalized the radial coordinates to the cluster 𝑅 h from Harris (1996, revision of 2010) and the dispersion profiles to the central dispersion determined in Libralato et al. (2022) from HST data. Four clusters do not have determination of the central velocity dispersion from HST data, namely NGC 1904, NGC 4147, NGC 6712 and NGC7492. For these clusters, we adopted the central dispersion determined in Vasiliev & Baumgardt (2021). We remind here that such approach can be used only to compare dispersion profiles, but not mean motions, as the latter are not on a relative scale.  \nThe normalized global dispersion and anisotropy profiles as a function of distance from cluster center for the 28 analyzed GCs are shown in Fig. 3. Global profiles for 1P and 2P stars are shown in the first and second columns, while the comparison between the average trend is shown in the third column. The average profiles were determined using the Locally Weighted Regression (LOESS) as implemented in Cappellari et al. (2013). To estimate the uncertainties in the average trend, we repeated the LOESS fitting on 1000 bootstrapped samples, drawing the values of 𝜎 R / T in each realization from a Gaussian distribution centered on the observed value with a dispersion equal to the observed uncertainties. The uncertainties were determined as the 16 th and 84 th percentiles of the 1000 LOESS fits. This is a more robust approach than simply bootstrapping the sample, as it accounts for observed uncertainties as well. Finally, we estimated the significance of the observed differences between 1P and 2P stars accounting for uncertainties in both data and fitting. The details of the procedure are outlined in App. C. In a nutshell, for each of the 1000 realizations of the LOESS fit, we computed the difference between 1P and 2P as a function of the radial coordinate. The 1, 2 and 3 𝜎 confidence regions are shown as gray shaded regions in the rightmost panels, where we display the difference between 1P and 2P a function of normalized radius. Additionally, we quantified the statistical significance of each point as the fraction of simulations which returned a difference smaller than the observed one, at the same radial location. The color of the line in the rightmost panels is indicative of of the statistical significance, in units of 𝜎 , as indicated in the right colorbar.  \nWhile there are no differences between the global radial dispersion profiles of 1P and 2P stars, 2P stars exhibit a lower tangential dispersion with respect to 1P stars, especially in the outer regions, e.g. 𝑅 > 1 -2 𝑅 h . As a consequence, the 1P displays an isotropic/tangentially anisotropic motion in the outer regions, while the 2P is radially anisotropic. Due to the large uncertainties of the mean trend, the observed differences are significant beyond the 3 𝜎 level only between 2 and 4 𝑅 h . Nonetheless, the magnitude of the differences is consistent with the predictions from the theoretical simulations (see e.g. Tiongco et al. 2019; Vesperini et al. 2021). Additionally, the global profiles within 𝑅 h are consistent with the results of Libralato et al. (2023).  \nFigure 2. Dynamical profiles of NGC 0104, NGC 2808, NGC 5904 and NGC 6205. From top to bottom: mean radial and tangential profile, radial and tangential dispersion profiles and anisotropy profile. Red/blue lines refer to 1P and 2P stars respectively, while the corresponding uncertainties are represented by the shaded boxes. The region where HST data, if available, have been used is indicated by the vertical green line, while the dot-dashed and dashed and lines mark the core and half-light radii respectively. Parsec and km/s units are shown respectively in the top and right axis with gray colors.  \n<!-- image -->  \nFigure 3. Global dynamical dispersion and anisotropy profiles for the 28 analyzed clusters. Velocity dispersion profiles are normalized to the central dispersion value from Libralato et al. (2022), and radial coordinates are normalized to the clusters' half-light radius. LOESS trends and uncertainties are represented by solid lines and shaded regions. The first and second columns show the LOESS fit for 1P and 2P separately, while the third column displays their comparison. The fourth column presents the difference between the 1P and 2P global profiles along with its statistical significance. The 1, 2, and 3 𝜎 confidence intervals are indicated by gray shaded regions, and the significance of the difference is represented by the color of the lines. Points with significance lower than 1 𝜎 are shown in the same color for clarity.  \n<!-- image -->", '4 DISCUSSION': "In the last decade, several theoretical works based on 𝑁 -body simulations have investigated the dynamical evolution of MPs in GCs and their potential for distinguishing between different formation scenarios. The velocity dispersion and anisotropy profiles of MPs can provide fundamental insights into the origins of these stellar populations (see e.g. Hénault-Brunet et al. 2015; Mastrobuono-Battisti & Perets 2013, 2016; Mastrobuono-Battisti & Perets 2021; Tiongco et al. 2019; Vesperini et al. 2021; Lacchin et al. 2022; Hypki et al. 2022; Pavlík & Vesperini 2022; Pavlík et al. 2024). However, due to dynamical evolution and spatial mixing, possible differences in the internal dynamics of MPs disappear as the clusters reach relaxation. As shown in Fig.2 for individual clusters and Fig.3 for all clusters, the tangential dispersion and anisotropy profiles of 1P and 2P stars exhibit some differences, with 1P being isotropic across the entire cluster field and 2P stars shifting from isotropy to radially anisotropy beyond 𝑅 h .  \nTo investigate how these dynamical differences depend on the clusters' properties evolutionary state, we divided our sample of 28 clusters into different groups based on various internal and external factors: dynamical ages (determined as the ratio between cluster age and half-mass relaxation time, see e.g. Libralato et al. 2022), tidal filling factor defined as R = 𝑅 h / 𝑅 J (Baumgardt et al. 2010; Shin et al. 2013) with 𝑅 J being the Jacobi radius derived in Balbinot & Gieles (2018) 9 , spatial mixing of MPs (as determined in Leitinger et al. 2023, and Jang et al. in preparation), tidal interaction with the Milky Way (determined on the basis of the peri-Galactic radius, see e.g. Milone et al. 2020; Libralato et al. 2023) 10 , clusters' origin (e.g. in situ or accreted) as identified in Massari et al. (2019) and escape velocity ( 𝑣 esc, see e.g. Baumgardt & Hilker 2018; MastrobuonoBattisti & Perets 2021) 11 . The groups are shown in the table provided in App. C. We find worth mentioning that different groups have overlap, as some clusters' properties are not independent one from the other (e.g. the spatial mixing of 1P and 2P is connected with the cluster dynamical age Dalessandro et al. 2019; Leitinger et al. 2023). Additionally, we also divided clusters based on their internal bulk rotation (as determined in Vasiliev & Baumgardt 2021), but, to keep the figure simpler, we not show the results in Fig. 4. Nonetheless, the difference between 1P and 2P anisotropy profiles for the groups and the statistical significance are included in App. C.  \nThe dynamical profiles for different groups of clusters are displayed in Fig. 4. Each group is indicated in the top inset, together with the number of clusters in each group. The global profiles (solid lines) have been computed with the LOESS algorithm, while the uncertainties (shaded regions, 1 𝜎 ) have been determined by bootstrapping with replacements a 1000 times and accounting for the observational uncertainties (see Sec. 3.2 for a detailed description). The significance of the differences between the anisotropy profiles of 1P and 2P for each group is presented in App. C, with the description of the procedure adopted to estimate it.  \nWe find significant differences, i.e. above the 3 𝜎 level, in the av-  \nge anisotropy profiles of many of the analyzed groups. Overall, we find that the 2P is more radially anisotropic beyond 𝑅 h , while 1P is isotropic and become slightly tangentially anisotropic beyond 23 𝑅 h . There are important differences in the 1P-2P relative dynamics amongsomeofthe analyzed groups. For example, while dynamically young clusters and those with centrally concentrated 2P stars show clearly different dispersion and anisotropy profiles, these differences are less pronounced in intermediate-age and mixed-population clusters. These findings are qualitatively consistent with the simulations by Tiongco et al. (2019, see, e.g., their Fig. 9) and Vesperini et al. (2021, see, e.g., their Figs. 13-17), which suggest that due to different initial spatial configurations (such as centrally concentrated 2P stars) or initial isotropy, 1P stars tend to develop isotropic motion, while 2P stars display more radially anisotropic motion, especially in the outer regions. Over time, these dynamical differences gradually diminish as spatial and dynamical mixing of MPs occurs.  \nThe comparison between clusters with small and large periGalactic radii (shown in the first two columns of the bottom row of Fig. 4) reveals some interesting patterns. In clusters with small peri-Galactic radii, 1P stars exhibit tangential anisotropy beyond 𝑅 h , while in clusters with larger peri-Galactic radii, 1P stars are either isotropic or slightly radially anisotropic. In both cases, 2P stars become increasingly radially anisotropic as they move outward. These findings suggest that the interaction with the Milky Way's tidal field plays a significant role in the dynamical evolution of MPs. Moreover, in clusters with small peri-Galactic radii, 1P stars appear to behave like tidally filling systems, where the outer regions develop slight tangential anisotropy, whereas 2P stars resemble tidally underfilling systems with radial anisotropy. In contrast, weaker interactions with the Milky Way's tidal field result in a more isotropic 1P population. This scenario aligns qualitatively with the theoretical models presented in Tiongco et al. (2019); Hypki et al. (2022); Pavlík et al. (2024).  \nTo further explore this hypothesis, we also examined the internal dynamics of MPs in clusters classified as tidally underfilling ( R < 0 . 05) and tidally filling ( R > 0 . 05). While tidally underfilling clusters exhibit relative differences consistent with other groups, such as isotropic 1P and radially anisotropic 2P, MPs in tidally filling clusters display a different pattern. Specifically, 1P stars are tangentially anisotropic, while 2P stars exhibit isotropic motion across the entire field. Interpreting these results is challenging, as the tidal filling factor R is closely linked to other properties, such as 𝑣 esc, 𝑅 peri , and cluster mass. Specifically, tidally underfilling clusters tend to be more massive, with larger 𝑣 esc and smaller 𝑅 peri .  \nAqualitatively similar pattern is observed in clusters with low and high escape velocities, with the latter often exhibiting tangentially anisotropic 1P stars. This difference may be due to the preferential loss of 1P stars on radial orbits in clusters with higher 𝑣 esc (see, e.g., Mastrobuono-Battisti & Perets 2021; Lacchin et al. 2024). However, the relationship between 𝑣 esc and the fraction of MPs remains complex and is highly dependent on the specific formation scenario adopted for MPs.  \n1P and 2P stars exhibit comparable relative differences in clusters, regardless of whether they are in situ or accreted, and whether they are rotating or non-rotating (see, e.g. App. C). In these cases, 2P stars remain radially anisotropic, while 1P stars display isotropic motion across the entire cluster field.  \nIn summary, the results presented in Fig. 4 indicate that clusters in a less dynamically evolved state show significant dynamical differences among MPs. Additionally, the interaction with the Galaxy appears to play a crucial role in shaping the evolution of different populations. These findings align with the conclusions of Libralato  \net al. (2023) regarding the internal dynamics of MPs in the innermost regions.", '5 SUMMARY AND CONCLUSION': "In this study, we present the first comprehensive investigation of the internal dynamics of multiple populations across a wide field of view, from the innermost arcmin to the clusters' outskirts, in a large sample of 28 Galactic GCs. Using HST, ground-based and XP synthethic photometry, we identified first- and second-population stars among RGB stars from the innermost regions to the clusters' outskirts. To achieve this, we exploited the pseudo CMD dubbed Chromosome maps as well as the 𝐶 UBI index, where 1P and 2P stars form well-separated groups. The internal dynamics of MPs was investigated using HST and Gaia DR3 proper motions, allowing us to determine the mean, dispersion, and anisotropy profiles as function of distance from the cluster center. The analyzed cluster regions range from 0.2 𝑅 h to more than 10 𝑅 h . The analysis presented in this work confirms and extends the results of previous studies based on HST and Gaia proper motions (see e.g. Richer et al. 2013; Bellini et al. 2015, 2018; Milone et al. 2018; Cordoni et al. 2020a; Cordoni et al. 2020b; Libralato et al. 2023; Cadelano et al. 2024).  \nOur analysis of individual clusters reveals distinct dynamical differences in the anisotropy profiles of 1P and 2P stars in several cases. On average, 2P stars are more radially anisotropic beyond 𝑅 h , whereas 1P stars generally exhibit isotropic motion, with some showing tangential anisotropy beyond 2-3 𝑅 h . These anisotropy differences between 1P and 2P stars are primarily driven by a lower tangential dispersion of 2P, with no significant differences observed in the radial component. These results agree with the analysis of Libralato et al. (2023) for the innermost cluster regions.  \nStudying the global trends in dispersion and anisotropy profiles, derived by combining the dynamical profiles of various clusters, we observe significant differences between 1P and 2P stars, especially outside 𝑅 h where the differences are significant beyond the 3 𝜎 level. In dynamically young clusters, 1P stars are isotropic in the inner regions but become slightly tangentially anisotropic toward the outskirts. In contrast, 2P stars are isotropic in the cluster centers and become radially anisotropic in the outer regions. These patterns are also evident in clusters where 2P stars are more centrally concentrated and in clusters with escape velocities exceeding 20 km / s. However, dynamical evolved clusters, with spatially mixed MPs and lower escape velocities do not exhibit these dynamical differences. Rotating/non-rotating clusters, and in situ/accreted clusters show similar relative differences between 1P and 2P stars.  \nWe also explore the influence of the Milky Way's tidal field on the dynamical properties of MPs by analyzing clusters with different peri-Galactic radii. Our findings reveal a strong connection between the dynamical behavior of 1P stars and the strength of the Milky Way's tidal field. In clusters with orbits closer to the Galactic center, where the tidal field is stronger, 1P stars tend to exhibit tangential anisotropy beyond 1-2 𝑅 h . Conversely, in clusters with weaker interactions with the Milky Way's tidal field, 1P stars display isotropic motion. This suggests that the Milky Way's tidal field plays a crucial role in the dynamical evolution of MPs. Further supporting this conclusion, our analysis shows that tidally underfilling and filling clusters exhibit distinct relative patterns in the dynamical profiles of 1P and 2P stars.  \nThe observed differences in the internal dynamics between 1P and 2P stars qualitatively align with the predictions of 𝑁 -body and theoretical simulations by Mastrobuono-Battisti & Perets (2016);  \nTiongco et al. (2019); Vesperini et al. (2021); Pavlík et al. (2024). These studies indicate that the distinct dynamical properties of 1P and 2P stars are indicative of 2P stars forming in a more centrally concentrated environment.  \nOverall, the analysis presented in this paper offers key insights into the formation scenarios of multiple stellar populations and their relationship with both internal factors (such as 𝑣 esc and dynamical age) and external influences (such as interaction with the Milky Way's tidal field).", 'ACKNOWLEDGMENTS': "SJ acknowledges support from the NRF of Korea (2022R1A2C3002992, 2022R1A6A1A03053472). EPL acknowledges support from the 'Science & Technology Champion Project' (202005AB160002) and from the 'Top Team Project' (202305AT350002), all funded by the 'Yunnan Revitalization Talent Support Program'. TZ acknowledges funding from the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie Grant Agreement No. 101034319 and from the European Union - Next Generation EU. This work has received funding from 'PRIN 2022 2022MMEB9W Understanding the formation of globular clusters with their multiple stellar generations ' (PI Anna F. Marino), and from INAF Research GTO-Grant Normal RSN2-1.05.12.05.10 - (ref. Anna F. Marino) of the 'Bando INAF per il Finanziamento della Ricerca Fondamentale 2022'. EL acknowledges financial support from the European Research Council for the ERC Consolidator grant DEMOBLACK, under contract no. 770017", 'DATA AVAILABILITY': 'Relevant data underlying this work is available in the article. All other data will be shared upon reasonable request to the corresponding author.', 'REFERENCES': 'Figure 4. Comparison between dispersion and anisotropy profiles of clusters with different properties, as indicated in the insets. Details about each group and the median trend is discussed in the text. Red and blue colors indicate 1P and 2P stars. The statistical significance of the difference between 1P and 2P trend is shown in App. C.  \n<!-- image -->  \nDalessandro E., et al., 2019, The Astrophysical Journal Letters, 884, L24 De Angeli F., et al., 2023, A&A, 674, A2  \nDecressin T., Meynet G., Charbonnel C., Prantzos N., Ekström S., 2007, åp, 464, 1029  \nDenissenkov P. A., Hartwick F. D. A., 2014, MNRAS, 437, L21  \nForeman-Mackey D., Hogg D. W., Lang D., Goodman J., 2013, PASP, 125,  \n306  \nGaia Collaboration et al., 2018, A&A, 616, A12  \nGaia Collaboration et al., 2023a, A&A, 674, A1  \nGaia Collaboration et al., 2023b, A&A, 674, A33  \nGieles M., et al., 2018, MNRAS, 478, 2461  \nGratton R. G., Carretta E., Bragaglia A., 2012, A&ARv, 20, 50  \nHarris W. E., 1996, AJ, 112, 1487  \nHypki A., Giersz M., Hong J., Leveque A., Askar A., Belloni D., Otulakowska-Hypka M., 2022, Monthly Notices of the Royal Astronomical Society, 517, 4768  \nHénault-Brunet V., Gieles M., Agertz O., Read J. I., 2015, MNRAS, 450,', 'APPENDIX A: SELECTION OF MULTIPLE POPULATIONS': 'In this Appendix we show the selection of MPs for all the 28 analyzed clusters. As discussed in Sec. 2, we used a combination of HST and ground-based photometry (Milone et al. 2017; Nardiello et al. 2018; Stetson et al. 2019; Jang et al. 2022), and Gaia XP synthetic photometry (Mehta et al. 2024). In Fig. A1- A3 we display the analyzed field of view in the central panel, together with the photometric diagrams used to separate 1P and 2P, always shown with red and blue colors respectively. Specifically, ground-based photometry is shown in the leftmost panels, with the ChMs from Jang et al. (2022) shown in the top one. HST and Gaia XP ChM are instead displayed in the right  \npanels, respectively in the top and bottom one. The coloured circles in the central panels indicate the analyzed regions.  \nConcerning ground-based observations, we first used the selection in the ChM whenever available, and we adopted the 𝐶 UBI selection for stars without ChM information. MPs outside the orange circles are selected by means of Gaia XP synthethic photometry. We remind here that Gaia XP synthethic photometry is only available for four clusters, namely NGC 0104, NGC 3201, NGC 6121 and NGC 6752 (see Mehta et al. 2024, for a complete discussion).', 'APPENDIX B: INDIVIDUAL PROFILES OF ALL CLUSTERS': "We present in this section the individual dynamical profiles of the analyzed clusters. The profiles are displayed in a similar fashion as Fig. 2. Mean trends are shown in Fig. B1-B2, while dispersion and anisotropy profiles are presented in Fig. B3-B5. Values computed from HST/Gaia data are displayed as crosses and circles respectively, while their uncertainties are indicated by the shaded rectangles. Mean motions and dispersion profiles have been converted to km / s adopting the clusters' distances derived in Baumgardt & Hilker (2018). Red and blue colors indicate 1P and 2P stars, while black diamonds represent the values of the 1D-dispersion computed in Libralato et al. (2022). Finally, the locations of the 𝑅 c , 𝑅 h , 𝑅 t radii are indicated by the dash-dotted, solid and dashed lines respectively. 𝑅 t is indicated only if contained within the analyzed field of view.", 'APPENDIX C: STATISTICAL SIGNIFICANCE OF THE OBSERVED GLOBAL PROFILES.': 'To assess the statistical significance of the differences in the global dynamical profiles of 1P and 2P shown in Fig. 3 and 4, we adopted the following procedure. We refer for simplicity to the case with all clusters. First we created a 1000 realizations of bootstrapped samples of 𝑦 = 𝜎 R , 𝜎 T , 𝛽 , where each value has bee drawn from a Gaussian with dispersion equal to the observed uncertainty. For each realization we repeated the LOESS fit for 1P and 2P, and computed the difference the between the two fits. Finally, we estimate the 1 , 2 , 3 𝜎 confidence intervals from the distribution of simulated difference at each 𝑅 / 𝑅 h . The three confidence intervals are shown as gray shaded regions in the top panels of Fig. C1. In addition to this analysis, we also directly quantified the fluctuations introduced solely by the uncertainties. To do this, we repeated the same analysis on two identical distributions, testing the null hypothesis that 1P and 2P stars share the same dynamical profiles, and determined the fraction ( 𝑓 ) of simulations with differences larger than the observed ones. Such fraction indicates the probability that uncertainties alone can reproduce a difference between 1P and 2P as large as the observed ones, and thus we computed the significance of the difference as 𝑝 = 1 -𝑓 . The value of 𝑝 is indicated by the color of the lines in the top panels of Fig. C1, as indicated by the colorbar. We stress here that, instead of determining one single value of significance for each 1P/2P dynamical profiles, we compute the statistical significance for each radial coordinate.  \nThis paper has been typeset from a T E X/L A T E X file prepared by the author.  \nFigure A1. Photometric diagrams used to identify 1P and 2P stars across the cluster field of view. The left panels display photometric diagrams built with 𝑈𝐵𝑉𝐼 ground-based photometry from Stetson et al. (2019), while ChMs built with space-based observations are shown on the right side (HST ChMs in the top panel and Gaia XP ChMs in the bottom panel). The central panels indicate the analyzed field of view, with different circles indicating the extension of each dataset. We only show panels whenever a dataset is available. Details of the selection are discussed in Sec. 2.  \n<!-- image -->  \nFigure A2. Same as Fig. A1.  \n<!-- image -->  \nFigure A3. Same as Fig. A1.  \n<!-- image -->  \nFigure B1. Mean motion of the radial component, determined as discussed in the Sec. 3. Red and blue colors represent the profiles of 1P and 2P stars. The radial coordinate is normalized to the half-light radius, indicated in the inset in units of arcmin.  \n<!-- image -->  \nFigure B2. Same as Fig. B1 for the tangential component.  \n<!-- image -->  \nFigure B3. Dispersion profile along the radial direction. The values of the 2-dimensional dispersion computed in the cluster center and at the core and half-light radii in Libralato et al. (2022) are indicated by black diamonds, while values and uncertainties computed in this work for 1P and 2P stars are shown with red and blue colors, respectively.  \n<!-- image -->  \nFigure B4. Same as Fig. B4 for the tangential component.  \n<!-- image -->  \nFigure B5. Anisotropy profile, defined as 𝛽 = 𝜎 T / 𝜎 R -1. Red and blue colors indicate results for 1P and 2P stars.  \n<!-- image -->  \nFigure C1. Difference between the 1P and 2P global anisotropy profiles for the groups discussed in Sec. 4 and its statistical significance. In each row, the bottom panels display the global 1P and 2P profiles, along with their uncertainties (red and blue lines and shaded regions). The top panels illustrate the difference between the profiles and its uncertainties. The color of the lines indicates the significance of the difference at each 𝑅 / 𝑅 h , as shown by the colorbar. For clarity, points with statistical significance below 1 , 𝜎 are represented with the same color.  \n<!-- image -->  \nTable C1. Clusters groups introduced and analyzed in Sec. 4.'}