Open Access
Issue
A&A
Volume 692, December 2024
Article Number L1
Number of page(s) 10
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/202451381
Published online 28 November 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

In the hierarchical scenario of galaxy formation, structures grow inside out (Gott 1975). This means that newcomers have lower binding energies than satellites that entered a main galaxy host at early epochs. Because of the host mass growth, first-comers are naturally most strongly bound (Rocha et al. 2012; Boylan-Kolchin et al. 2013). Cosmological simulations recovered a tight linear correlation between the lookback infall time and the logarithm of the binding energy, showing an evolution of more than 1 dex, and with a scatter of only 0.13 dex (Rocha et al. 2012). The slope of this correlation is consistent with the slope estimated for the Milky Way (MW) accretion history. Assuming the MW mass model from Eilers et al. (2019, total mass of 8.3 × 1011 M), Hammer et al. (2023, see their Fig. 6) determined this relation for the MW. The authors accounted for globular clusters (GCs) associated with the bulge, Kraken, Gaia-Sausage-Enceladus (GSE), and Sgr infall and adopted the associations made by Malhan et al. (2022) and Kruijssen et al. (2020). According to the latter study, these events occurred 12.5 ± 0.5, 11.5 ± 0.5, 9 ± 1, and 5 ± 1 Gyr ago, respectively. The relation between their lookback infall time and the logarithm of their binding energy was found to be linear, which agrees with cosmological simulations.

The relation was used by Hammer et al. (2023) to infer that the binding energies of MW dwarf galaxies are far lower than the binding energy of GSE GCs (by a factor of 6 on average). This prevented them from having been accreted ∼9 Gyr ago. Extrapolating the relation from bulge to Sgr GCs toward the dwarf energy regime, Hammer et al. (2023) concluded that most dwarf infall lookback times should be shorter than 3 Gyr. The relation could also constrain the MW mass, because the more massive a galaxy, the deeper its potential well. This allows it to capture satellites with lower binding energies (e.g., Leo I could be bound if the MW were very massive; Boylan-Kolchin et al. 2013). Recent Gaia measurements of the MW rotation curve (RC) provided 8.3 × 1011 M (Eilers et al. 2019) and 1.4 × 1012 M (Jiao et al. 2023) for an adopted Navarro et al. (1997, Navarro, Frenk & White, hereafter NFW) and Einasto (1965, see also Retana-Montenegro et al. 2012) mass profile, respectively. However, the NFW mass profile is too shallow in the disk outskirts, and it is rejected at 3σ level by RC measurements from the 3rd Gaia data release (hereafter Gaia DR3), which also revealed a decline beyond 19 kpc that is consistent with Keplerian expectations (Wang et al. 2023; Jiao et al. 2023; Ou et al. 2024). The consideration of the MW gravitational potential as being almost equivalent to a point mass1 beyond 19 kpc has profoundly impacted the cosmological community, for whom the MW halo-limiting radius was generally assumed to range from 150 to 300 kpc.

However, the MW is not entirely axisymmetric or in dynamical equilibrium, as is assumed when resolving the Jeans equation to derive the rotation curve. In particular, the disk shows many substructures, including ridges, warps, and flares, in particular, in its outer range, with different upward or inward velocities, and a significant difference in the stellar rotational velocity above and below the disk (Gaia Collaboration 2021, and references therein). Recently, Koop et al. (2024) suggested that the MW disk is so perturbed by passages of satellites that its rotation curve cannot be used to predict its mass. However, the regularity of the MW rotation curve until ∼19 kpc suggests that nonequilibrium mechanisms are not sufficiently strong to perturb its dynamics significantly, although they could be a serious concern beyond this radius.

The goal of this paper is to verify which MW mass can be consistent with the relation between lookback infall time and the binding energy of satellites after we compare the mass to the masses obtained in cosmological simulations. Furthermore, we assess the accretion epoch of dwarf galaxies, and we verify whether satellites may impact the disk stability when they pass near the disk. Section 2 tests the accretion history of the MW, including that of dwarf galaxies, and compares it to the history determined from high-resolution cosmological simulations. Section 3 discusses the pertinence of cosmological simulations for a retrieval of the MW properties and its accretion history, and for identifying whether satellites may have affected its disk dynamics. We then summarize our main conclusions.

In this paper, R200c is the virial radius for which the enclosed DM mass, dubbed M200c, corresponds to an overdensity of 200 times the critical density. The total mass includes M200c, to which we added the baryonic component, in a similar way as Eilers et al. (2019). The orbital energies correspond to the addition of the kinetic energy provided by Gaia proper motions and radial velocities with the potential energy of the considered mass distribution. To calculate the potential, we followed Rocha et al. (2012), who assumed that the DM-halo potential rises to zero at 1 Mpc. We therefore assumed that the DM potential associated with a given virial mass M200c must reach zero at the rescaled radius of 1 (M200c/1.4 × 1012)1/3 Mpc, where 6.2 × 1010 M is the DM-halo mass M200c of Rocha et al.2. The total mass (or potential) is the sum over the DM and the bulge and disk components. As an example, detailed equations are given in Appendix A for an NFW DM-halo component.

2. Cosmological simulations compared to the MW RC and accretion history

The top panel of Figure 1 compares the binding energies of the simulated subhalos according to Rocha et al. (2012)3 (black dots) with those of GCs that are associated with the bulge, GSE, and Sgr (red dots). The Rocha et al. (2012) model has a high total mass (M200c = 9.12 × 1011 M). However, this model faces two fundamental difficulties: (1) The binding energies of the simulated subhalos are lower by two to three times than those of Sgr and the GSE GCs, and to reconcile them, extremely high lookback time values are required for these MW substructures. (2) The dark matter (DM) halo is far too massive to fit the MW RC (see the top panel of Fig. 2 and compare the dotted black line to the long-dash red line). The addition of a baryonic component to the simulations of Rocha et al. (2012) would not help us to reproduce the locations of bulge, Kraken, GSE, and Sgr points (see the bottom panel of Fig. 1) because a concentrated mass component would worsen the offset by increasing their binding energy.

thumbnail Fig. 1.

Relation between the infall lookback time and the logarithm of the binding energy, for which black dots represent the subhalos that are associated with the Rocha et al. (2012) simulation, the solid line represents the best fit, and the dashed lines show its one σ scatter. The binding energies of dwarf galaxies (blue triangles) and of GCs associated with MW events (red dots with corresponding labels in magenta) are calculated using the observed kinetic energy together with the potential of the main halo. The blue triangles indicate the location of dwarf galaxies with a very good accuracy (< 0.1 dex) in binding energy, and for which we adopted the predicted dwarf infall time of Barmentloo & Cautun (2023). Top panel: Relation as derived from the Rocha et al. (2012) dark-matter-only simulations. Bottom panel: Same as the top panel, but we added a baryonic mass (1.4 × 1012 M) after a slight rescaling of the dark matter halo mass from 1.34 × 1012 M to M200c = 1.4 × 1012 M (see Appendix A). In the bottom right corner of the panel, the numbers indicate the predicted time for GCs associated with the bulge, GSE, and Sgr when they are inserted into the halo potential, as well as the probability that they lie in the distribution of simulated subhalos.

In Appendix A we tried to rescale the simulations by Rocha et al. (2012) to different values of the total mass within the virial radius. We failed to find a value for which they fit the locations of GCs associated with Sgr, GSE, Kraken, and the bulge, however. Thus, while Rocha et al. (2012) reproduced the slope of the relation of infall time versus binding energy well, the kinetic energy of the simulated subhalos is far higher than that of typical GC systems associated with Sgr, GSE, Kraken, and the bulge even after the total mass is rescaled.

It is then useful to compare the MW accretion history to high-resolution simulated halos with a full treatment of the baryonic physics. We examined and extracted infall times and binding energies for satellites around Milky Way-mass galaxies in different high-resolution zoom-in hydrodynamical cosmological simulations: Seven from FIRE (Wetzel et al. 2023), and 15 from Auriga (Grand et al. 2024) simulated MW-mass galaxies (see details in Appendix B). Most FIRE-simulated MW-mass galaxies do not reproduce the RC slope at large radii (see Figs. 2 and B.3), which is expected because their halo radii are large (R200c from 168 to 251 kpc), which implies significant amounts of mass outside of ∼20 kpc. They do not reproduce the slope of the relation between infall epoch and binding energy either, which ranges from 27 to 37 Gyr km−2 s2, that is, two to three times larger than observed (12.2 Gyr km−2 s2). This motivated us to also investigate Auriga-simulated MW-mass galaxies, including their low MW-mass subsample. Their low MW mass halos 3 and 8 show a rotation velocity that is consistent with observed values around 15–20 kpc (see Fig. B.5). However, their RCs are too flat and cannot provide a reasonable fit of the observed MW RC. Furthermore, even low MW-mass halos (e.g., numberes 1 and 9) underestimate the MW RC velocities, which would require a more concentrated DM component. Most of the Auriga-simulated MW-mass galaxies show very to extremely steep slopes for the relation between infall epoch and binding energy and predict binding energies that are ∼0.4–0.8 dex too low at the epochs of the bulge formation and of the Kraken, GSE, and Sgr merging events. This is far more problematic.

thumbnail Fig. 2.

Rotation curve and MW accretion history for the FIRE model m12c. Top panel: Rotation curve from an isolated MW-like galaxy from FIRE (halo m12c) compared to Gaia data. The blue points show DR2 (Eilers et al. 2019) without systematics, and the red points show DR3 (Jiao et al. 2023), which includes measurements from Wang et al. (2023) plus systematic errors. The contributions from the DM halo, bulge, and disk are plotted as long-dash red, short-dash green, and short-dash blue lines, respectively. To illustrate the impact of the very massive halo from Rocha et al. (2012), we also added its contribution (dotted black line), which far exceeds the FIRE m12c DM halo. Bottom panel: Same as the bottom panel of Figure 1, but for the FIRE halo m12c. The binding energies of dwarf galaxies and of GCs associated with MW events were calculated using the observed kinetic energy together with the potential of each main halo. The top left corner of the panel displays the combined probability that the nine dwarfs have infall lookback times longer then 8 Gyr, as predicted by Barmentloo & Cautun (2023). In the bottom right corner of the panel, the numbers indicate the predicted time for GCs associated with the bulge, GSE, and Sgr when they are inserted into the m12c halo potential, as well as the probability that they lie in the distribution of subhalos.

Among all the 22 simulated halos, FIRE m12c (Garrison-Kimmel et al. 2019) was found to reproduce both the MW RC and its accretion history best (see Fig. 2). Its total mass is × 1012 M, which is 10% higher than in the models adopted by Eilers et al. (2019) or by Bovy (2015), with which it shares many similarities. The top panel of Figure 2 suggests that this model reproduces the Gaia DR3 MW RC relatively well (associated χ2 probability P = 0.33). When this is compared to the RC from Eilers et al. (2019), the associated χ2 probability becomes extremely low because the latter authors did not account for systematic errors. A better estimate requires a significant improvement of the systematic error estimates, which is the subject of a future paper (Jiao et al. in preparation). Furthermore, the m12c predicted slope of the relation between the infall epoch and binding energy (31 Gyr km−2 s2) is more than twice what is observed in the MW. When we account for the fitting of the MW accretion history, the m12c halo might be discarded with a probability of P = 2 × 10−4, which is much better than any other tested simulated halo, however. We also note that FIRE halos m12b and m12m might be rescaled in total mass and would then reproduce the relation of infall time versus energy with the same significance as halo m12c, but their slopes would again be too high in this relation.

Barmentloo & Cautun (2023) pioneered the use of the accretion history relation to infer the dwarf galaxy infall times. They used 1628 simulated MW halos from EAGLE (Grand et al. 2017) with total masses from 0.5 to 2 M200c =1.3 × 1012 M. We adopted their predicted infall lookback times for 14 dSphs, which include Carina, Draco, Fornax, Sculptor, Sextans, and Ursa Minor (in addition to Bootes I and II, Coma Berenices, Grus I, Hercules, Segue II, Triangulum II, and Ursa Major II; see the blue triangles in Fig. 1), for which the binding energy is accurately determined with an error smaller than 0.1 dex. Figure 1 confirms the prediction by Rocha et al. (2012) that most dSph infall lookback times are within 8–9 Gyr. This requires a very high MW mass that is similar to the mass derived when dwarf orbits are assumed to be at equilibrium with the MW potential (e.g., 3 × 1012 M from Li et al. 2017). However, a fit of the MW accretion history relation requires a lower mass, such as that of FIRE m12c. The bottom panel of Figure 2 shows that an infall lookback time of 5–9 Gyr would lead for most dwarfs to an extremely low probability (see the value in the top left corner), while values below 3 Gyr are very likely, except for Segue II. In other words, when a cosmological simulated halo is able to fit the MW relation between the binding energy and the infall time, it automatically implies very low probabilities for most dwarfs to have entered the halo at a similar epoch as the GSE. The infall times of Barmentloo & Cautun (2023) have very large uncertainties, and their prediction of dwarf infall times is affected by the use of many halos that do not fit the MW accretion history and RC and by their use of several additional parameters that are insensitive to the infall time (see the discussion in the Appendix of Hammer et al. 2024).

3. Discussion and summary

We have presented the relation between infall epoch and binding energy, which is representative of the MW accretion history. We discuss whether this relation can be sufficiently robust to be predictive. For example, Pagnini et al. (2023) have questioned the reliability of associating GCs with past merger events in the MW (Kruijssen et al. 2019, 2020; Malhan et al. 2022) after they simulated mergers of two galaxies which had their own GC systems. They found GCs that progressively lost their energy after they were stripped in a 1:10 minor encounter, which prevents their use for a dating of the epoch of their simulated mergers. However, most GCs are expected to be formed during strong star formation events that occurred 12–9 Gyr ago (Haywood et al. 2016), that is, during the Kraken and GSE merger events (De Lucia et al. 2024; Valenzuela et al. 2024). It is doubtful whether the 1:10 mass ratio adopted by Pagnini et al. (2023) is representative of the two latter events. For example, Naidu et al. (2021) considered higher mass ratios4 of 1:2 to 1:4. Only a modest fraction of MW GCs can then be accreted through the mechanism proposed by Pagnini et al. (2023), which might provide a complementary channel for GCs that are not identified inside a structure in the plane of total energy – angular momentum (Hammer et al. 2023, see their Fig. 5). In addition, D’Souza & Bell (2022) argued that the relation of the infall epoch and binding energy might become very uncertain for halos with very active merger histories. The authors analyzed the ELVIS suite of 48 simulated host halos with a total mass ranging from 1 to 2.05 × 1011 M (Garrison-Kimmel et al. 2014). Among them, D’Souza & Bell (2022, see their Fig. 9) identified 3 host halos with a rich merger history such that the scatter of the relation almost equaled its amplitude. However, the MW is known for its relatively quiescent merger history when compared to spiral galaxies of similar masses (Hammer et al. 2007). For example, while an average spiral galaxy experimented its last major merger 6 Gyr ago, this occurred 9 Gyr ago (GSE) for the MW. The latter is thus unlikely to be one of the few galaxies with the richest merger history, which agrees with the well-identified linear relation (Hammer et al. 2023).

A comparison of the MW accretion history relation to that of cosmological simulations may provide us an additional constraint on the MW mass. However, the mass constraint only applies on the farthest considered GCs, that is, those attached to Sgr, with an average distance of ∼21 kpc. Within the latter radius, we find a mass of 2.05 × 1011 M for the m12c simulated halo (see Fig. 2). This value agrees well with estimates for GCs by Watkins et al. (2019, 2.15 × 1011 M) and with the values extracted from Gaia DR2 Eilers et al. (2019, 1.95 × 1011 M) and DR3 Jiao et al. (2023, 1012 M) RCs.

Cosmological simulations that can reproduce the binding energy of GSE GCs (see Fig. 2 and the three first panels of Fig. B.3) inferred that dwarf galaxies cannot have been bound to the MW 8–10 Gyr ago. It confirms the conjecture by Hammer et al. (2023) that they are newcomers (infall lookback time < 3 Gyr). This also suggests that they are likely out-of-equilibrium systems (see the detailed explanations and modeling in Hammer et al. 2024 and in Wang et al. 2024). This suggests that a derivation of the mass from dwarf galaxy orbits after considering them as long-term MW satellites systematically overestimates the total MW mass. This is confirmed by the fact that an MW mass in excess of 6 × 1010 M can neither reproduce the binding energy of the bulge, GSE, and Sgr GCs (see Figs. 1 and B.3) nor its RC.

We compared our results to cosmological hydrodynamical simulations from FIRE and Auriga and found the results at odds with both the MW RC and its accretion history relation. Many simulated galaxies predicted an RC and accretion history that are very different from the observed events. In particular, Rocha et al. (2012) and Auriga-simulated galaxies show dynamically hot subhalo systems with binding energies that are lower by three to nine times than those of the bulge, GSE, and Sgr GCs. This suggests that cosmological simulations would benefit from investigating different initial conditions. For example, it could be useful to investigate different initial mass fluctuations especially toward the low-mass end to allow the possibility of more recent major mergers, as are expected for understanding the formation of spiral galaxies and the acquisition of their angular momentum (Hammer et al. 2009; Hopkins et al. 2010; Puech et al. 2012).

Our comparisons of the relation of infall epoch and binding energy are limited to cosmological simulations that might provide a too shallow halo-mass distribution, that is, those following an NFW mass profile. If the Gaia DR3 RC with a Keplerian decline (Jiao et al. 2023; Ou et al. 2024) is confirmed, simulations of more concentrated DM halos with lower total masses would be required, and we would need to verify whether they can reproduce the MW RC and the accretion history. However, it would also be necessary to generate simulated halos that are sharply truncated, as suggested by the MW RC.

Finally, the relation of the infall time to energy may help us to verify whether the MW RC can be significantly altered by the passages of heavy satellites through the disk, as suggested by Koop et al. (2024)5. Their study compared such an impact after assuming the model of Laporte et al. (2018) of Sgr with an initial mass of 3.8 × 109 M. However, a much lower initial mass for Sgr (M200m =1.9 × 1012 M) is required6 to reproduce the Sgr stream (Vasiliev et al. 2021), which according to the latter authors would need a reanalysis of the role of Sgr in seeding the spiral pattern in phase space that was discovered in the MW disk (Gaia Collaboration 2021). Concerning the impact of Sgr on the MW RC, we may consider that decreasing a perturbing mass by a factor ∼16 would reduce the predicted changes in the rotation velocity by Koop et al. (2024) by about four times when the perturbations are in the linear regime. This implies that the changes in rotation velocity of 10–15 km s−1 reported by Koop et al. (2024, see their Fig. 11) would become consistent with the analyses of the systematic uncertainties (∼2% of the rotation velocity below R = 22 kpc) made by Jiao et al. (2023). Koop et al. (2024) also considered the six halos from Auriga (see Fig. B.4) to also infer large systematics caused by the passage of satellites similar to Sgr. However, Auriga subhalos are four to six times dynamically hotter than Sgr for a common infall lookback time of 5 Gyr, or in other words, the velocity dispersion of the subhalo system in Auriga is greater than the typical velocity of Sgr relative to the MW. The impact of the satellites on the MW disk would thus become quite modest and consistent with the expected systematic uncertainties (Sylos Labini et al. 2023; Jiao et al. 2023; Ou et al. 2024). A significant advantage of the above analysis is that it naturally explains why all Gaia RCs show a smooth decline from R = 8 to 18 kpc (Eilers et al. 2019; Mróz et al. 2019; Jiao et al. 2023; Ou et al. 2024), that is, why they do not show strong local variations in velocity that would be caused by massive satellite impacts on the MW disk.


1

The best-fit Einasto index found by Jiao et al. (2023) is consistent with a Gaussian cutoff of the dark matter density profile.

2

The M200c value was converted from the Rocha et al. (2012) × 1011 M, which corresponds to an overdensity of 200 times the mass density.

3

The infall time versus binding energy relation of Rocha et al. (2012) is based on a high-resolution (4100 M for dark matter particles) dark-matter-only simulation based on Via Lactea II (Diemand et al. 2007).

4

Notice that these larger mass ratios are necessary to explain the origin of both thin and thick disk of a spiral galaxy like the MW (Hammer et al. 2009, 2018; Hopkins et al. 2010; Athanassoula et al. 2016; Sauvaget et al. 2018).

5

Koop et al. (2024) also suggested a truncation of the MW disk profile from their Jeans equation analysis of Gaia data with radial velocities. Their analysis was based on Bayesian techniques for estimating distances, but their statistics beyond R = 17 kpc were lower than in Ou et al. (2024). We defer this to a discussion in a paper (Jiao et al., in preparation) that further analyzes the systematics associated with the Gaia DR3 data.

6

Both studies considered a similar MW total mass.

Acknowledgments

We are very grateful for the discussion we have had with Robert Grand. We warmly thank the referee whose comments have helped us to improve the paper. We are grateful for the support of the International Research Program Tianguan, which is an agreement between the CNRS in France, NAOC, IHEP, and the Yunnan Univ. in China. Y.-J.J. acknowledges financial support from the China Scholarship Council (CSC), No.202108070090. I. A. would like to thank the Graduate Program in Astrophysics of the Paris Sciences et Lettres (PSL) University for funding.

References

  1. Athanassoula, E., Rodionov, S. A., Peschken, N., & Lambert, J. C. 2016, ApJ, 821, 90 [NASA ADS] [CrossRef] [Google Scholar]
  2. Barmentloo, S., & Cautun, M. 2023, MNRAS, 520, 1704 [Google Scholar]
  3. Bovy, J. 2015, ApJS, 216, 29 [NASA ADS] [CrossRef] [Google Scholar]
  4. Boylan-Kolchin, M., Bullock, J. S., Sohn, S. T., Besla, G., & van der Marel, R. P. 2013, ApJ, 768, 140 [NASA ADS] [CrossRef] [Google Scholar]
  5. De Lucia, G., Kruijssen, J. M. D., Trujillo-Gomez, S., Hirschmann, M., & Xie, L. 2024, MNRAS, 530, 2760 [Google Scholar]
  6. Diemand, J., Kuhlen, M., & Madau, P. 2007, ApJ, 667, 859 [NASA ADS] [CrossRef] [Google Scholar]
  7. D’Souza, R., & Bell, E. F. 2022, MNRAS, 512, 739 [CrossRef] [Google Scholar]
  8. Eilers, A.-C., Hogg, D. W., Rix, H.-W., & Ness, M. K. 2019, ApJ, 871, 120 [Google Scholar]
  9. Einasto, J. 1965, Trudy Astrofizicheskogo Instituta Alma-Ata, 5, 87 [NASA ADS] [Google Scholar]
  10. Gaia Collaboration (Antoja, T., et al.) 2021, A&A, 649, A8 [EDP Sciences] [Google Scholar]
  11. Garrison-Kimmel, S., Boylan-Kolchin, M., Bullock, J. S., & Lee, K. 2014, MNRAS, 438, 2578 [NASA ADS] [CrossRef] [Google Scholar]
  12. Garrison-Kimmel, S., Hopkins, P. F., Wetzel, A., et al. 2019, MNRAS, 487, 1380 [NASA ADS] [CrossRef] [Google Scholar]
  13. Gott, J. R. 1975, ApJ, 201, 296 [Google Scholar]
  14. Grand, R. J. J., Gómez, F. A., Marinacci, F., et al. 2017, MNRAS, 467, 179 [NASA ADS] [Google Scholar]
  15. Grand, R. J. J., Fragkoudi, F., Gómez, F. A., et al. 2024, MNRAS, 532, 1814 [Google Scholar]
  16. Hammer, F., Puech, M., Chemin, L., Flores, H., & Lehnert, M. D. 2007, ApJ, 662, 322 [NASA ADS] [CrossRef] [Google Scholar]
  17. Hammer, F., Flores, H., Puech, M., et al. 2009, A&A, 507, 1313 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Hammer, F., Yang, Y. B., Wang, J. L., et al. 2018, MNRAS, 475, 2754 [NASA ADS] [Google Scholar]
  19. Hammer, F., Li, H., Mamon, G. A., et al. 2023, MNRAS, 519, 5059 [CrossRef] [Google Scholar]
  20. Hammer, F., Wang, J., Mamon, G. A., et al. 2024, MNRAS, 527, 2718 [Google Scholar]
  21. Haywood, M., Lehnert, M. D., Di Matteo, P., et al. 2016, A&A, 589, A66 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Hopkins, P. F., Bundy, K., Croton, D., et al. 2010, ApJ, 715, 202 [Google Scholar]
  23. Jiao, Y., Hammer, F., Wang, H., et al. 2023, A&A, 678, A208 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Koop, O., Antoja, T., Helmi, A., Callingham, T. M., & Laporte, C. F. P. 2024, A&A, in press, https://doi.org/10.1051/0004-6361/202450911 [Google Scholar]
  25. Kruijssen, J. M. D., Pfeffer, J. L., Reina-Campos, M., Crain, R. A., & Bastian, N. 2019, MNRAS, 486, 3180 [Google Scholar]
  26. Kruijssen, J. M. D., Pfeffer, J. L., Chevance, M., et al. 2020, MNRAS, 498, 2472 [NASA ADS] [CrossRef] [Google Scholar]
  27. Laporte, C. F. P., Johnston, K. V., Gómez, F. A., Garavito-Camargo, N., & Besla, G. 2018, MNRAS, 481, 286 [Google Scholar]
  28. Li, Z.-Z., Jing, Y. P., Qian, Y.-Z., Yuan, Z., & Zhao, D.-H. 2017, ApJ, 850, 116 [NASA ADS] [CrossRef] [Google Scholar]
  29. Malhan, K., Ibata, R. A., Sharma, S., et al. 2022, ApJ, 926, 107 [NASA ADS] [CrossRef] [Google Scholar]
  30. Mróz, P., Udalski, A., Skowron, D. M., et al. 2019, ApJ, 870, L10 [Google Scholar]
  31. Naidu, R. P., Conroy, C., Bonaca, A., et al. 2021, ApJ, 923, 92 [NASA ADS] [CrossRef] [Google Scholar]
  32. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 [Google Scholar]
  33. Ou, X., Eilers, A.-C., Necib, L., & Frebel, A. 2024, MNRAS, 528, 693 [NASA ADS] [CrossRef] [Google Scholar]
  34. Pagnini, G., Di Matteo, P., Khoperskov, S., et al. 2023, A&A, 673, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Puech, M., Hammer, F., Hopkins, P. F., et al. 2012, ApJ, 753, 128 [Google Scholar]
  37. Retana-Montenegro, E., van Hese, E., Gentile, G., Baes, M., & Frutos-Alfaro, F. 2012, A&A, 540, A70 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Rocha, M., Peter, A. H. G., & Bullock, J. 2012, MNRAS, 425, 231 [NASA ADS] [CrossRef] [Google Scholar]
  39. Sauvaget, T., Hammer, F., Puech, M., et al. 2018, MNRAS, 473, 2521 [Google Scholar]
  40. Sylos Labini, F., Chrobáková, Ž., Capuzzo-Dolcetta, R., & López-Corredoira, M. 2023, ApJ, 945, 3 [NASA ADS] [CrossRef] [Google Scholar]
  41. Valenzuela, L. M., Remus, R. S., McKenzie, M., & Forbes, D. A. 2024, A&A, 687, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Vasiliev, E., Belokurov, V., & Erkal, D. 2021, MNRAS, 501, 2279 [NASA ADS] [CrossRef] [Google Scholar]
  43. Wang, H.-F., Chrobáková, Ž., López-Corredoira, M., & Sylos Labini, F. 2023, ApJ, 942, 12 [NASA ADS] [CrossRef] [Google Scholar]
  44. Wang, J., Hammer, F., Yang, Y., et al. 2024, MNRAS, 527, 7144 [Google Scholar]
  45. Watkins, L. L., van der Marel, R. P., Sohn, S. T., & Evans, N. W. 2019, ApJ, 873, 118 [NASA ADS] [CrossRef] [Google Scholar]
  46. Wetzel, A., Hayward, C. C., Sanderson, R. E., et al. 2023, ApJS, 265, 44 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Re-scaling binding energy and orbital times for halos of different mass

Figure A.1 illustrates the major problem for the simulations made by Rocha et al. (2012) if it is used to reproduce the MW accretion history after comparing the location of bulge, GSE and Sgr GCs to that of the simulated subhalos. We have added a baryonic component (bulge and disk) very similar to that used by Jiao et al. (2023) and Ou et al. (2024), as done in the bottom panel of Figure 1. The two panels show that for any re-scaling up or down of the Rocha et al. (2012) main halo, the subhalos are less bound to their host galaxy (rescaled to a possible Milky Way mass) than the bulge, Kraken, GSE or Sgr GCs are to the Milky Way. Since we kept the kinetic terms of the binding energies of these four systems from the bulk line-of-sight velocities and Gaia bulk proper motions, while adopting the re-scaled potential of the simulation for their potential terms, any rescaling of the simulation will push the positions of these four systems (red points) in the same direction as are pushed the subhalos (black points), although by a smaller amount.

thumbnail Fig. A.1.

Same as Fig. 1, but for different mass re-scalings of the Rocha et al. (2012) simulations. In the bottom right of each panel the numbers indicate the predicted time for GCs associated with the bulge, GSE, and Sgr when they are inserted into the m12c halo potential, as well as the probability they can lie in the distribution of subhalos. (Top:) We have scaled down the Rocha et al. (2012) main halo mass to 7.4× 1010 M and add a baryonic component of 6.2× 1011 M, similar to that used in Jiao et al. (2023) and Ou et al. (2024). (Bottom panel:) We have scaled up the Rocha et al. (2012) main halo mass to 19.4× 1010 M, also adding a baryonic mass of 6.2× 1011 M.

Here we show that the binding energies and infall times of subhalos in the simulations of Rocha et al. (2012) can be converted for halos having similar concentrations but different masses. Indeed, the binding energy at the ‘virial’ radius R200c follows Ebinding(R200c)∝V200c2, where V200c2 = GM200c/R200c is the square of the circular velocity at the virial radius. Since, at any epoch, the mean density within the virial radius is the same for all halos, it leads to R 200 c M 200 c 1 / 3 $ R_{\mathrm{200c}} \propto M_{\mathrm{200c}}^{1/3} $, and V 200 c R 200 c M 200 c 1 / 3 $ V_{\mathrm{200c}} \propto R_{\mathrm{200c}} \propto M_{\mathrm{200c}}^{1/3} $. One then finds that E binding ( R 200 c ) M 200 c 2 / 3 $ E_{\mathrm{binding}}(R_{\mathrm{200c}}) \propto M_{\mathrm{200c}}^{2/3} $. Since the concentration of halos varies weakly with halo mass (c ∝ M−0.1, Navarro et al. 1997), decreasing the mass by a multiplicative factor μ = 1.9 leads to only a 6 per cent increase of the concentration. Thus, to first order, re-scaling the mass is self-similar and should not affect the orbital parameters of satellites, if considered in virial units. Therefore, decreasing the mass by a factor of μ causes their binding energies to be divided by μ2/3. On the other hand, the circular orbital times of dwarfs around the MW at its virial radius, torb, circ ∝ R200c/V200c, will be independent of μ. Therefore, the orbital times for general orbits within the MW halo will be independent of μ, and so will be the infall times, as well as the lookback times from entry to a given fraction of the current virial radius.

In summary, re-scaling the virial mass of the MW by a factor μ should shift the location of the infall lookback time vs. binding energy of dSphs by a linear factor of (2/3) log(μ) along the logarithmic energy axis, without affecting the infall lookback time. This will conserve the slope of the correlation. Since Rocha et al. (2012) assumed a potential that rises to zero at 1 Mpc, we have scaled all the considered halo potentials in this paper to reach zero at a radius Rlim equal to 1 Mpc multiplied by (M200c/1.4 × 1012)1/3. In Appendix B we verify whether models having a full treatment of the baryonic physics may help to identify simulated galaxies having properties similar to that of the MW.

Here, we give an example of how the virial mass (M200c) and the potential are calculated for the Rocha et al. (2012) main halo mass of M200c= 14× 1011 M that is calculated within the virial radius R200c= 231 kpc.

The DM (spherical) mass profile for our assumed NFW model is:

M ( R ) = M 200 c ln ( R / a + 1 ) R / a / ( R / a + 1 ) ln ( c + 1 ) c / ( c + 1 ) , $$ \begin{aligned} M(R)= M_{200c}\, \frac{\ln (R/a+1)-R/a\,/\,(R/a+1)}{\ln (c+1)-c/(c+1)}\ , \end{aligned} $$(A.1)

where c = 11.48 = R200c/a is the concentration, while a = 20.12 kpc is the NFW scale radius. The gravitational potential corresponding to the NFW model is given by:

Φ ( R ) = G M 200 c / R 200 c ln ( c + 1 ) c / ( c + 1 ) × [ ln ( R / a + 1 ) R / R 200 c ( R 200 c R lim ) ln ( R lim a + 1 ) ] = G M 200 c / R ln ( c + 1 ) c / ( c + 1 ) × [ ln ( R a + 1 ) ( R R lim ) ln ( R lim a + 1 ) ] . $$ \begin{aligned} \Phi (R)&= - \frac{G\,M_{200c}/R_{200c}}{\ln (c+1)-c/(c+1)} \nonumber \\&\qquad \times \left[ \frac{\ln (R/a+1)}{R/R_{200c}} - \left(\frac{R_{200c}}{R_{\rm lim}}\, \right) \ln \left(\frac{R_{\rm lim}}{a}+1\right) \right] \nonumber \\&= -\frac{G\,M_{200c}/R}{\ln (c+1)-c/(c+1)} \nonumber \\&\qquad \times \left[ \ln \left(\frac{R}{a}+1\right) - \left(\frac{R}{R_{\rm lim}}\, \right) \ln \left(\frac{R_{\rm lim}}{a}+1\right) \right] \ . \end{aligned} $$(A.2)

The last term in the square brackets accounts for the potential going to zero at Rlim= 1 Mpc. The total mass and total potential can then be calculated by adding baryonic mass and potential such as described in Jiao et al. (2023) and Ou et al. (2024).

Appendix B: Cosmological models compared to the observed MW RC and accretion history

Cosmological hydrodynamical simulations bring widely-used resources providing important predictions for galaxy formation and evolution. In particular "zoom-in" simulations allow to investigate numerical resolutions associated with 104 − 105 (103 − 104) M per DM (baryonic) particle (Wetzel et al. 2023; Grand et al. 2024), respectively. These simulations include many ingredients of baryonic physics allowing comparison to very detailed observations of the MW substructures and its halo satellites. FIRE (Wetzel et al. 2023) and Auriga (Grand et al. 2024) have considered MW-like (or MW-mass) galaxies with total mass ranging from M200c= 7.3 to 14× 1011 M, and from 5.3 to 19 5 × 104 M, respectively. The above mass ranges, while consistent with the value found by Gaia DR2 RC (Eilers et al. 2019), are larger than that derived from the Gaia DR3 RC (Jiao et al. 2023; Sylos Labini et al. 2023; Ou et al. 2024). We have tested their predictions on the observed RC and accretion history of the MW, to verify whether this could bring additional constraints on the MW halo mass. We have selected halos from data release of both the FIRE (Wetzel et al. 2023) and Auriga (Grand et al. 2024). We have preselected the 7 isolated halos from FIRE, and 15 Auriga halos, including 6 high-mass MW with high-resolution, plus the 9 low-mass MW halos. Our goal was to account for isolated simulated galaxies as well as to account for the simulations with the highest resolution as well as considering the largest possible range for the simulated halo masses.

Each subhalo is represented by its rotation curve (solid black lines in the top panels) that is compared with Gaia DR3 RCs, and its accretion history relation (bottom panels). In the latter, infall lookback times have been retrieved from the data base for FIRE simulations, while we have derived them for Auriga simulations.

To do this, we have generated the whole orbital history of subhalos and compare their orbits to the evolution of R200. Rocha et al. (2012) used the virial radius R200m for estimating the infall times as well as it has been done by FIRE, while Auriga provides the virial radius R200c at which the mean density is 200 times the critical density. For consistency reasons, we have calculated the virial radius R200m of 15 MW-like simulations from Auriga at z = 0 (the last snapshot) with Ωm = 0.30966 (Planck Collaboration VI 2020). We perform a linear fit on these radii and then roughly convert all R200c values of all simulations at each redshift to R200m using this linear relation as presented in Figure B.1. Finally we assume the infall epoch to be that of the last entry of the subhalo within R200m as shown in the top panel of Figure B.2. Simulation predictions for both RC and MW accretion history are given in Figure B.3 for FIRE, and in Figures B.4 (6 high-mass MW with high-resolution) and B.5 (9 low-mass MW with low-resolution) for Auriga.

thumbnail Fig. B.1.

Linear relation of R200m and R200c for different MW-like galaxies realized from Auriga simulations.

During this exercise, we have realized that some subhalos were lacking of information about their orbital positions (e.g. bottom panel in Figure B.2), which render unreliable the infall time estimate. Going one step further, we found this to occur essentially for high-resolution (resolution 3, 3.5 × 104 M per DM particle) simulated subhalos with less than 100 particles. Therefore, we have adopted a secure mass limit for subhalos to be 107 M, a value we have also adopted for FIRE simulated subhalos (4 × 105 M per DM particle). We have also use the same mass limit for low-resolution (resolution 4, 1012 M per DM particle) low-mass, Auriga subhalos because otherwise their numbers become too small to allow an efficient fit of the MW accretion history relation. It means that the latter relation is less precise for low-mass Auriga halos L1 to L10 (the last 10 pair of panels in Figure B.5).

thumbnail Fig. B.2.

Example of calculation of infall time for Auriga simulations. Blue points and dash lines indicate the distance of subhalos 1 (top panel) and 1749 (bottom) from Auriga MW-like Level 3 resolution Halo 27 to the main host galaxy center. Red points and lines indicate the evolution of the virial radius (R200m) of the main host galaxy at different epochs using the linear relation in Fig B.1. Vertical red dashed lines indicate the last infall time of the subhalo.

thumbnail Fig. B.3.

Rotation curves (top panels) and accretion history (bottom panels) from isolated MW-mass galaxies from FIRE. Same symbols for points and lines than in Figure 2. For simplicity, here we have only show the six classical dwarfs (blue triangles) and bulge, GSE and Sgr GCs (red crosses). The top (bottom) panels include simulated galaxies with total mass M200c smaller (higher) than .

thumbnail Fig. B.4.

Rotation curves (top) and accretion history (bottom) from high MW-mass galaxies from Auriga. Same symbols for points and lines than in Figure B.3. Some Auriga simulations have been considered to be representative of the infall of GSE or of the LMC, which is indicated in the bottom panel of the corresponding simulations.

thumbnail Fig. B.5.

Rotation curves (top) and accretion history (bottom) from low MW-mass galaxies from Auriga. Same symbols for points and lines than in Figure B.4.

All Figures

thumbnail Fig. 1.

Relation between the infall lookback time and the logarithm of the binding energy, for which black dots represent the subhalos that are associated with the Rocha et al. (2012) simulation, the solid line represents the best fit, and the dashed lines show its one σ scatter. The binding energies of dwarf galaxies (blue triangles) and of GCs associated with MW events (red dots with corresponding labels in magenta) are calculated using the observed kinetic energy together with the potential of the main halo. The blue triangles indicate the location of dwarf galaxies with a very good accuracy (< 0.1 dex) in binding energy, and for which we adopted the predicted dwarf infall time of Barmentloo & Cautun (2023). Top panel: Relation as derived from the Rocha et al. (2012) dark-matter-only simulations. Bottom panel: Same as the top panel, but we added a baryonic mass (1.4 × 1012 M) after a slight rescaling of the dark matter halo mass from 1.34 × 1012 M to M200c = 1.4 × 1012 M (see Appendix A). In the bottom right corner of the panel, the numbers indicate the predicted time for GCs associated with the bulge, GSE, and Sgr when they are inserted into the halo potential, as well as the probability that they lie in the distribution of simulated subhalos.

In the text
thumbnail Fig. 2.

Rotation curve and MW accretion history for the FIRE model m12c. Top panel: Rotation curve from an isolated MW-like galaxy from FIRE (halo m12c) compared to Gaia data. The blue points show DR2 (Eilers et al. 2019) without systematics, and the red points show DR3 (Jiao et al. 2023), which includes measurements from Wang et al. (2023) plus systematic errors. The contributions from the DM halo, bulge, and disk are plotted as long-dash red, short-dash green, and short-dash blue lines, respectively. To illustrate the impact of the very massive halo from Rocha et al. (2012), we also added its contribution (dotted black line), which far exceeds the FIRE m12c DM halo. Bottom panel: Same as the bottom panel of Figure 1, but for the FIRE halo m12c. The binding energies of dwarf galaxies and of GCs associated with MW events were calculated using the observed kinetic energy together with the potential of each main halo. The top left corner of the panel displays the combined probability that the nine dwarfs have infall lookback times longer then 8 Gyr, as predicted by Barmentloo & Cautun (2023). In the bottom right corner of the panel, the numbers indicate the predicted time for GCs associated with the bulge, GSE, and Sgr when they are inserted into the m12c halo potential, as well as the probability that they lie in the distribution of subhalos.

In the text
thumbnail Fig. A.1.

Same as Fig. 1, but for different mass re-scalings of the Rocha et al. (2012) simulations. In the bottom right of each panel the numbers indicate the predicted time for GCs associated with the bulge, GSE, and Sgr when they are inserted into the m12c halo potential, as well as the probability they can lie in the distribution of subhalos. (Top:) We have scaled down the Rocha et al. (2012) main halo mass to 7.4× 1010 M and add a baryonic component of 6.2× 1011 M, similar to that used in Jiao et al. (2023) and Ou et al. (2024). (Bottom panel:) We have scaled up the Rocha et al. (2012) main halo mass to 19.4× 1010 M, also adding a baryonic mass of 6.2× 1011 M.

In the text
thumbnail Fig. B.1.

Linear relation of R200m and R200c for different MW-like galaxies realized from Auriga simulations.

In the text
thumbnail Fig. B.2.

Example of calculation of infall time for Auriga simulations. Blue points and dash lines indicate the distance of subhalos 1 (top panel) and 1749 (bottom) from Auriga MW-like Level 3 resolution Halo 27 to the main host galaxy center. Red points and lines indicate the evolution of the virial radius (R200m) of the main host galaxy at different epochs using the linear relation in Fig B.1. Vertical red dashed lines indicate the last infall time of the subhalo.

In the text
thumbnail Fig. B.3.

Rotation curves (top panels) and accretion history (bottom panels) from isolated MW-mass galaxies from FIRE. Same symbols for points and lines than in Figure 2. For simplicity, here we have only show the six classical dwarfs (blue triangles) and bulge, GSE and Sgr GCs (red crosses). The top (bottom) panels include simulated galaxies with total mass M200c smaller (higher) than .

In the text
thumbnail Fig. B.4.

Rotation curves (top) and accretion history (bottom) from high MW-mass galaxies from Auriga. Same symbols for points and lines than in Figure B.3. Some Auriga simulations have been considered to be representative of the infall of GSE or of the LMC, which is indicated in the bottom panel of the corresponding simulations.

In the text
thumbnail Fig. B.5.

Rotation curves (top) and accretion history (bottom) from low MW-mass galaxies from Auriga. Same symbols for points and lines than in Figure B.4.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.