Long-term orbital evolution of Galactic satellites and the effects on their star formation histories

We investigate the orbital motions of dwarf spheroidal galaxies (dSphs) in the halo of the Milky Way (MW) to understand their possible effects on the diversity of the star formation histories seen in these MW satellites. In this work, we explicitly consider a time-varying gravitational potential due to the growth of the MW's dark halo mass to calculate the long-term orbital evolutions of the dSphs, guided with {\it Gaia} DR2 proper motions, over the past 13.5 billion years. We find that the infall time of a satellite, defined at which the galaxy first crosses within the growing virial radius of the MW's halo, coincides well with the time when the star formation rate (SFR) is peaked for the sample of classical dSphs. On the other hand, ultra-faint dSphs already finished their SF activity prior to their infall times as already suggested in previous works, but there is a signature that their earlier SF histories are affected by interaction with the growing MW's halo to some extent. We also find, for classical dSphs, that the relative fraction of stars formed after the peak of the SFR to the current stellar mass is smaller for the smaller pericentric radius of the galaxy at its first infall. These results suggest that the infalling properties of the dSphs into the MW and the resultant environmental effects such as ram-pressure stripping and/or tidal disturbance in the MW's dark halo containing hot gas play important roles in their star formation histories.


INTRODUCTION
Our understanding of how dwarf spheroidal galaxies (dSphs) in the Milky Way (MW) have formed and evolved to what we see today remains still far from being complete. The distribution of the member stars in a color-magnitude diagram indicates that each of the MW's dSphs shows a variety of time evolution of star formation rate (SFR) (e.g., Mateo 1998;Grebel 1999). The diversity of these SF histories may be summarized into the following different groups: (1) SF occurred only in the early phase of the dSph, leading to the dominance of old stellar populations as old as 12 Gyrs, (2) SF occurred mostly in the near past, as inferred from the dominance of relatively young stellar populations as young as ∼ 5 Gyrs, (3) SF has been occurred over several Gyrs in the middle of the galaxy history, and (4) the galaxy had experienced episodic SF events (e.g., Buonanno et al. 1999;Grebel & Gallagher 2004;Tolstoy et al. 2009;de Boer et al. 2012;McConnachie 2012;de Boer et al. 2014).
In particular, deep photometric observations of Galactic satellites with Hubble Space Telescope (HST) by Weisz et al. (2014) clearly indicate that the so-called classical dSphs in the MW, having a V -band absolute magnitude, M V , brighter than −8 mag, show diverse SF histories at given M V . On the other hand, ultra-faint dwarf galaxies (UFDs) with M V > −8 mag contain only old member stars, whereby the SFR was peaked and quenched at early epochs (Okamoto et al. 2008;Brown et al. 2012Brown et al. , 2014Weisz et al. 2014;Simon 2019). Stars in dSphs also show their characteristic iron abundance as well as abundance-ratio distributions (e.g., Tolstoy et al. 2004;Koch et al. 2006;Battaglia et al. 2006;Tolstoy et al. 2009;Kirby et al. 2011Kirby et al. , 2013Ishigaki et al. 2014;Tsujimoto et al. 2015;Simon 2019). These chemical properties of dSphs must be intimately related to the diversity of their SF histories, but what causes this diversity is unsolved yet, including their intrinsic properties or external effects (e.g., Mayer et al. 2006;Revaz et al. 2009;McConnachie 2012;Okayasu & Chiba 2016;Bermejo-Climent et al. 2018;Revaz & Jablonka 2018;Escala et al. 2018).
One of the key ingredients that controls the SF histories in Galactic dSphs is their past and current environment within the halo of the MW. This includes the tidal effects from the gravitational field of the MW, the rampressure stripping of cold gas in their progenitor galaxies in the presence of the MW's hot halo gas, the photoionizing effect of UV light from the MW as well as from the Universe, and so on. Indeed, recent hydrodynamical simulations of galaxy formation in the framework of Λdominated cold dark matter theory suggest that all of these physical processes are actually at work on each of dark-matter subhalos falling into a host, MW-sized halo, which may eventually become currently observed luminous satellites (e.g., Wetzel et al. 2015;Macciò et al. 2017;Frings et al. 2017;Simpson et al. 2018;Buck et al. 2019;Genina et al. 2019;Garrison-Kimmel et al. 2019;Rodriguez Wimberly et al. 2019;Fillingham et al. 2019).
To assess these environmental effects on the SF histories of Galactic dSphs, it is important to investigate when and how these dSphs are falling into and eventually orbiting in the gravitational field of the MW's dark halo and how these orbital evolutions are associated with the past SF events in each of the dSphs. This approach is possible only when the reliable kinematical information are available for many of the different dSphs showing different stellar populations. Several previous works have already suggested, based on the calculations of the orbits of Galactic dSphs, that their close passages to the MW can be linked with the SF histories, including the peak and the subsequent time evolution of the SFR (e.g., Sohn et al. 2007;Pasetto et al. 2011;Rocha et al. 2012;Fillingham et al. 2019;Rusakov et al. 2020).
In this respect, Gaia DR2 has revolutionized both precisions and amounts in the astrometric data of Galactic stars, so the estimation of the precise spatial motions of many different dSphs, including both classical dSphs and UFDs, is now possible. Based on the calibration of the proper motions from Gaia DR2, Helmi et al. (2018) calculated the orbits of the nine classical dSphs and one UFD as well as 75 globular clusters for three different models of static Galactic gravitational potentials and showed the distribution of their orbital properties. Fritz et al. (2018) further derived the Gaia DR2 proper motions of more than 39 dwarf galaxies located out to 420 kpc from the Galactic Center and integrated their orbits in static canonical MW potentials. Similar studies for the orbits of 17 UFDs in a static Galactic potential are made by Simon (2018) based on the calibration of the Gaia DR2 proper motions for these satellites. Kallivayalil et al. (2018) showed that some of these satellites are associated with the Large Magellanic Cloud (LMC) (see also Yozin & Bekki 2015;Pardy et al. 2019;Erkal & Belokurov 2019;Patel et al. 2020).
However, we note that these orbital calculations of Galactic dSphs are made under the assumption that the Galactic gravitational field is fixed with the currently observed form. Although this assumption is appropriate for the relatively short period of a few Gyrs, it cannot be applied to the whole history of the MW, especially at the early stage of the first passage of satellites to the host, MW's dark halo, whose mass is growing through hierarchical accretion processes. Thus, it is not fully understood yet whether the orbital motions of dSphs play a major role in their SF histories through the relevant environmental effects.
This work intends to relax this assumption of a static Galactic potential and explores the long-term orbital evolutions of Galactic dSphs under the situation of a growing mass of the MW's dark halo. Similar calculations for the past orbits of the LMC or star clusters in a time-varying Galactic potential were made by Zhang et al. (2012) and Haghi et al. (2015) and we adopt here the method in these previous works.
Recently  performed the Phat ELVIS suite of 12 high-resolution, dissipationless simulations of forming MW-sized halos  to follow the orbital evolutions of subhalos relative to their host halo. They derived the infall time at which a subhalo first crosses within the virial radius of a growing MW-sized host halo. These subhalos are then matched with 37 Galactic satellites by comparing both in the diagram of binding energy vs. distance from host (Rocha et al. 2012), whereby the infall time inferred for each satellite is compared with its SF history.
Here, instead of using such extensive cosmological simulations, we adopt an analytically tractable, timevarying mass of a MW-sized host halo as in Zhang et al. (2012) and Haghi et al. (2015) to directly integrate the long-term orbital evolution of each of Galactic dSphs and investigate their first infall and subsequent orbital motions in comparison with a growing virial radius of the host halo. We then infer the possible relation of these orbital motions with their SF histories. As detailed below, we find an intriguing coincidence of the first infalling time with the peak of the SFR for many of Galactic dSphs.
The paper is organized as follows. Section 2 shows the method for the calculations of the orbital motions of the sample dSphs over many dynamical times of the Galaxy. Section 3 is devoted to the results of these calculations. In Section 4, we discuss the physical mechanism behind our calculated results by comparing with recent high-resolution simulations of galaxy formation and the conclusions are made.
For all the relevant calculations in what follows, we adopt the set of the cosmological parameters based on WMAP7 (Larson et al. 2011): Ω m = 0.266, Ω Λ = 0.734 and H 0 = 71 km s −1 Mpc −1 .

METHOD
We simply assume that the form of the Galactic potential at each epoch is spherically symmetric and is given by the so-called NFW profile (Navarro, Frenk & White 1996) where M vir , r vir and c denote the virial mass, virial radius and concentration parameter of the MW's dark halo, respectively. In this work, we explicitly consider the time evolution of these quantities over the past 13.5 Gyrs to follow the corresponding long-term orbital evolution of Galactic satellites. Regarding the growth of the MW's virial mass, M vir (t), Wechsler et al. (2002) obtained the approximate formula as a function of redshift, z, based on their cosmological N-body simulations, where a c controls the growth time scale of the MW's mass, which is given as a c = 0.34. In this case, the growth is slow, where the epoch when the MW's mass was half of the current value was about 8 Gyrs ago. For the current virial mass of the MW, we set M vir (z = 0) = 1.54 +0.75 −0.44 ×10 12 M ⊙ taken from the recent measurement of the spatial motions of globular clusters based on Gaia DR2's proper motions (Watkins et al. 2019) and also consider the effect of these 1 σ uncertainties in the M vir (z = 0) value on the calculation of the satellites' orbits. We note that Krumholz & Dekel (2012) also showed the mass accretion history of the MW's halo and the corresponding time evolution for the MW's dark halo mass is found to be nearly the same as that Equation (2) provides.
Given M vir at each epoch, the virial radius is estimated as the radius within which the mean density of the corresponding halo is 200 times as large as the critical density of the Universe, ρ c (z), namely, M vir = (4π/3)r 3 vir 200ρ c . This yields r vir (z), For the time evolution of c, we adopt the relation, c(M vir , z), given by Prada et al. (2012).
The proper motions of the sample of the dSphs that we use here are based on Gaia DR2 catalog. Here, we adopt the list of 6D data in the work by Riley et al. (2019), who assembled the distances, line-of-sight velocities and proper motions for 38 Galactic satellites. Among these, we select 8 classical dSphs (Carina, Draco, Fornax, Leo I, Leo II, Sculptor, Sextans and Ursa Minor) and 8 UFDs (Bootes I, Coma Berenices, Canes Venatici I (CVn I), Canes Venatici II (CVn II), Reticulum II, Segue I, Ursa Major I and Ursa Major II), which covers a wide range of orbital properties as shown below and for most of which SF histories are available from the HST observations by Weisz et al. (2014) and Brown et al. (2014). In this sample selection for the orbit calculation, we avoid the satellites having large uncertainties in the measured proper motions and thus 3D velocities, such as Leo IV and Hercules, and those being thought as LMC satellites in recent studies (Kallivayalil et al. 2018;Patel et al. 2020).
With the above set-up, we calculate the past orbit of each satellite using galpy (Bovy 2015) with a time step of 10 Myr. The Galactocentric distance of the Sun, its circular velocity and the local solar motion for this calculation are adopted as 8 kpc, 220 km s −1 and (11.1, 12.24, 7.25) km s −1 (Schönrich et al. 2010). We note that the virial radius is increasing with time as the virial mass of the MW grows. Thus, at the specific epoch, an infalling satellite first crosses within this increasing virial radius of the MW with time, and we define this time when it occurs as "the infall time" hereafter denoted as t first infall . We note that in these orbit calculations of Galactic satellites, the choice for the form of the Galactic potential is significantly affecting the resulting time-evolution of the orbit compared with the effects stemmed from the measurement errors of the observed quantities. The principal importance of the host halo mass in the orbits of the associated subhalos is also suggested from the analysis of their orbital properties from cosmological N-body simulation (Wetzel et al. 2011). Taking this regard into account, we here estimate the range of t first infall resulting from the choice of M vir (z = 0) within its 1 σ uncertainties, namely, 1.10 × 10 12 M ⊙ to 2.29 × 10 12 M ⊙ .  at early epochs and thus have executed many orbital oscillations around the MW, (3) Leo I is now located beyond r vir after having crossed within it about 2 Gyr ago and arrived at the pericentric radius about 1 Gyr ago (Sohn et al. 2007).
For these sample satellites, we derive the infall time, t first infall , in the currently adopted, time-varying Galactic potential. The comparison is then made with the infall time obtained in the simulation work of Fillingham et al. (2019), who analyze the orbits of subhalos in 12 MW-sized halos from the Phat ELVIS simulation ) and matching with the observed dSphs is made. As shown in Figure 2, the currently derived infall times are roughly in agreement with those in the work of Fillingham et al. (2019) within the uncertainties, except the case of CVn II perhaps due to the difference in the adopted mass models. This suggests that the current simplified treatment of the satellites' obits given in the Galactic potential of Equation (2) provides generally consistent results with those based on the high-resolution simulations of evolving MW-sized dark halos.
The relation between the infall time and the binding energy of each satellite is shown in Figure 3. The general trend of this relation is again consistent with that obtained from the cosmological simulations of MWsized halos (Rocha et al. 2012;Fillingham et al. 2019), namely a satellite that has fallen earlier is more strongly bound to the MW's halo located at the smaller Galactocentric radius.
To assess the effect of considering the growing mass of the MW on these orbits, we also calculate the case when the Galactic potential is fixed in the form of Equation (1) for all of our sample dSphs. The plots for all of these orbits are presented in Appendix. In short, we find that while these orbits show only a simple oscillation in r for a static Galactic potential, those in the currently timevarying Galactic potential start to notably deviate from these simple oscillations at the epochs about 4 Gyrs ago. Thus, to derive the past orbits of the dSphs more than 4 Gyrs ago, we need to explicitly take into account the time variation of the Galactic potential as adopted here.   (lower panel) and their SF histories (upper panel), where the latter is expressed in terms of the differential SFR as a function of look-back time using the cumulative SF history and total stellar mass available from the Weisz et al. (2014) work. Figure 6 is the same as these figures but for two UFDs, CVn I and CVn II, for which SF histories are again available from the same reference (Weisz et al. 2014), where Leo IV and Hercules therein are excluded in our analysis because of the large uncertainties in their measured 3D velocities. We also note that we confine ourselves to use the work of Weisz et al. (2014) for the source of the SF histories to avoid any systematics in their derivations associated with difference methods.

Comparison with star formation histories
It is evident from Figure 4 and 5 that the timing of the first crossing through r vir for each classical dSph, defined as the infall time here, occurs nearly at the same timing as when the SFR is peaked. For Sculptor, although its first infall did not cross r vir so the infall time is delayed somewhat, the epoch when it first reached its pericentric radius occurs nearly at the peak time of the SFR.
In contrast to these classical dSphs, two UFDs shown in Figure 6 had the peak of the SFR well before their infall times. This result is in agreement with that in Fillingham et al. (2019), which showed that the UFDs shut down their SF prior to the infall to the MW. However, it is worth remarking that in CVn I, the second peak of the SFR occurs nearly at the infall time, whereas in CVn II, the peak of the SFR is realized at around its first passage to the MW, although the pericentric radius then did not cross the virial radius of the MW. This suggests that even in UFDs, their SF histories may be affected by the early stages of the host, MW halo to some extent.
In Figure 7, we show the difference between the time when the SFR is peaked and the time of the first infall, t SF peak − t first infall , as a function of V-band absolute magnitude, M V , of each galaxy. The error bar for the time difference stems from the effect of the range of the adopted M vir (z = 0) on the estimate of t first infall . For UFDs other than CVn I and CVn II, we use the results of HST observation by Brown et al. (2014), which estimate the ages of the dominant old stars in their sample of UFDs. It is clear that in the classical dSphs, this time difference is confined only within a few Gyrs. Sculptor shows somewhat a large time difference because its first infall does not cross within the virial radius of the MW. Instead of using the infall time for this galaxy, we also plot, with the filled diamond, the difference between the time when the SFR is peaked and the time when it first reached the pericentric radius, which is now found to be small. In fact, the time of the first arrival at the pericentric radius occurs just after the infall time as long as this first infall crosses within the virial radius of the MW. In contrast, for UFDs, the peak of the SFR occurred much earlier than their first infall to the MW's virial radius, in agreement with Fillingham et al. (2019). This comparison between the orbital evolution and SF history of each classical dSph shown in Figure 4 and 5 suggests that not only the similarity between the infall time and the peak time of the SFR, but also the notable properties of the subsequent SF activity after the first pericentric passage are inferred, such that when the pericentric radius of the first infall is small, the SFR after its peak appears to be considerably reduced.
This property is presented in Figure 8, which shows the relation between the fraction of stars formed after  Figure 6. The same as Figure 4 and 5 but for the two UFDs, (a) CVn I and (b) CVn II for which star formation histories are available from Weisz et al. (2014). In contrast to classical dSphs, these UFDs finished star formation before the infall time. We note that in CVn I, the second peak of the SFR occurs nearly at the infall time, whereas in CVn II, the peak of the SFR occurred at around its first passage to the MW. Classical UFD Figure 7. The difference between the time when SFR is peaked and the time of the first infall, t SF peak − t first infall as a function of V-band absolute magnitude, MV , for the classical dSphs (filled blue squares) and UFDs (filled black squares). For Sculptor, the filled diamond denotes the difference between the time when SFR is peaked and the fime when it first reached the pericentric radius.
the peak of the SFR relative to each dSph's current stellar mass, M * , and the first pericentric radius. Again, the error bar for the latter quantity stems from the effect of the range of the adopted M vir (z = 0) on the estimate of t first infall . The clear tight correlation between these quantities is found, such that the relative amount of stars formed after its peak is reduced for the smaller pericentric radius at the first infall. This result is basically consistent with that in Fillingham et al. (2019), which showed that the quenching timescale of SF, i.e., the difference between the quenching time, defined as the look-back time when 90 % of the current stellar mass is formed, and infall time is short for the dSphs having high orbital eccentricities.  Figure 8. The relation between the fraction of the stars formed after the peak of the SFR relative to each dSph's current stellar mass, M * , and the galaxy's pericentric radius at its first infall to the MW. There is a remarkable correlation such that the smaller pericentric radius at the first infall leads to the reduction of star formation after its peak of the SFR.
We have investigated the long-term orbital motions of Galactic dSphs in the course of the growing mass of the MW's dark halo over the past 13.5 Gyrs. Their current motions are taken from the recent compilation (Riley et al. 2019), which adopts the high-precision measurements of proper motions from Gaia DR2. We have compared these orbital motions with the SF histories of the dSphs. It is found that the infall time of each classical dSph first crossing within the time-varying MW's virial radius coincides remarkably well with the time when its SFR is peaked. Also, in the classical dSph whose pericentric radius after the first infall is small, the SFR afterwards is reduced. Finally, we have confirmed that the formation of stars in UFDs is finished well before they first enter into the virial radius of the MW's dark halo as already suggested in previous works. For these UFDs however, we have found a signature that their earlier SF histories are subject to environmental effects provided by the early stages of the MW halo to some extent. The adopted prescription for the form of the Galactic potential and its time dependence given in Equation (1) and (2) is admittedly simplistic. For instance, the accretion process of dark matter onto a host MW-like halo is no longer spherically symmetric and continuous, but rather anisotropic and sporadic through merging/accretion of many subhalos. Also, the time evolution of the MW halo and the resultant orbits of the satellites can be affected by the environment associated with the formation of the Local Group, including the falling and binding Magellanic Clouds to the Galactic potential (Bekki & Chiba 2005;Yozin & Bekki 2015;Kallivayalil et al. 2018;Pardy et al. 2019;Erkal & Belokurov 2019;Patel et al. 2020). Nonetheless, it is interesting to note that the infall times of Galactic satellites derived here are generally in agreement of those based on dissipationless cosmological simulations by Fillingham et al. (2019) (Figure 2), suggesting that the adopted simplified form of the Galactic potential for the calculation of satellites' orbits is not far from reality. Said that, for a more comprehensive comparison with the SF histories, it will be important to consider the above other dissipationless processes that we ignore here as well as the effect of the later baryonic infall that deepens the Galactic potential, and also to refine the determination of the current MW mass and its distribution (e.g., Eilers et al. 2019;Hammer et al. 2020).
The correlation between the first infall time of a MW satellite and the time of its maximum SFR has been suggested from recent hydrodynamical simulations of galaxy formation, APOSTLE (Genina et al. 2019) and Auriga (Simpson et al. 2018). Indeed, these simulations show that the very efficient star formation is achieved when a subhalo containing cold interstellar gas is first approaching to its pericentric radius. The snapshots of gas densities shown in the Genina et al. (2019) suggest that gas inside a subhalo is compressed and removed due to ram pressure from hot gas in a MW-sized halo and that this effect is strongest around the pericentric radius, whereby star formation may be induced from gas compression. Simpson et al. (2018) from their highresolution, zoom-in cosmological simulations also show that the orbital motion of a subhalo having cold interstellar gas within a MW-sized dark halo affects both the star formation history and the time evolution of its gas fraction in each subhalo: the epoch when 90% of the mass of the currently observed stars in a subhalo is formed, τ 90 in their notation, as well as the timing of the reduction of its gas fraction appears to be well correlated with the first infall time into a host halo. It is also suggested from their simulations that a smaller pericentric radius at the first infall of a subhalo seems to yield a large reduction of its interstellar gas, thereby suggesting the suppression of subsequent star formation (see also the relevant simulation works by Macciò et al. 2017;Frings et al. 2017;Buck et al. 2019). These properties are naturally understood if interstellar gas in a subhalo is efficiently compressed and then removed through rampressure stripping and/or tidal disturbance. Indeed, the presence of a hot gas in the MW and its effect on satellite galaxies has been suggested and studied from both the observations of the spectra of background QSOs and theoretical models (e.g., Miller & Bregman 2013, 2015Emerick et al. 2016).
Finally, star formation activities in UFDs may be entirely determined by the first collapse of a subhalo in the early Universe and the subsequent quenching or suppression of star formation by reionization of the Universe (Bovill et al. 2009;Brown et al. 2014). Thus, the formation of stars in UFDs is already finished when they enter into the MW's dark halo. These properties for UFDs are actually suggested from the ELVIS simulations (Rodriguez Wimberly et al. 2019;Fillingham et al. 2019).
These dynamical effects of the orbits of Galactic satellites, especially the timing of the first infall in comparison with their star formation histories, can be imprinted in the chemical abundances of stars in their member stars (e.g., Koch et al. 2006;Kirby et al. 2011Kirby et al. , 2013 and also the density distribution of dark matter halos associated with these dSphs through external tides (e.g. Walker et al. 2009;Peñarrubia et al. 2010;Walker & Peñarrubia 2011;Bullock & Boylan-Kolchin 2017;Hayashi et al. 2017;Kaplinghat et al. 2019). For instance, in Carina, the multiple maxima of the SFR after the first infall seem to occur at the subsequent pericentric radii, and these phenomena may trigger gaseous infall and outflow, whereby governing the metallicity distribution of the member stars (Koch et al. 2006). Also, these star formation activities and associated feedback effects from supernovae can modify the density profile of a dark halo in dSphs (Read & Gilmore 2005;Ponntzen & Governato 2012).
These chemical and dynamical information of member stars in Galactic dSphs have been biased toward an inner part of each dSph compared to its nominal tidal radius, so that the global chemo-dynamical state, especially in an outer part or up to an outer boundary of each dSph, which is more sensitive to environmental effects, is yet largely unknown. Extensive spectroscopic data of stars in Galactic satellites out to their outskirts will be available from Prime Focus Spectrograph (PFS) to be attached to Subaru Telescope (Takada et al. 2014;Tamura et al. 2016). Using this fiber-fed, wide field of view spectrograph, Galactic Archaeology survey will allow us to get new insights into these subjects and thus understand the formation histories of the dSphs in the MW. Figure A2. The same as Figure A1 but for UFDs.