The Magellanic Stream at 20 kpc: A New Orbital History for the Magellanic Clouds

We present new simulations of the formation of the Magellanic Stream based on an updated first-passage interaction history for the Magellanic Clouds, including both the Galactic and Magellanic Coronae and a live dark matter halo for the Milky Way. This new interaction history is needed because previously successful orbits need updating to account for the Magellanic Corona and the loosely bound nature of the Magellanic Group. These orbits involve two tidal interactions over the last 3.5 Gyrs and reproduce the Stream's position and appearance on the sky, mass distribution, and velocity profile. Most importantly, our simulated Stream is only $\sim$20 kpc away from the Sun at its closest point, whereas previous first-infall models predicted a distance of $100-200$ kpc. This dramatic paradigm shift in the Stream's 3D position would have several important implications. First, estimates of the observed neutral and ionized masses would be reduced by a factor of $\sim$5. Second, the stellar component of the Stream is also predicted to be $<$20 kpc away. Third, the enhanced interactions with the MW's hot corona at this small distance would substantially shorten the Stream's lifetime. Finally, the MW's UV radiation field would be much stronger, potentially explaining the H$\alpha$ emission observed along most of the Stream. Our prediction of a 20 kpc Stream could be tested by searching for UV absorption lines towards distant MW halo stars projected onto the Stream.


INTRODUCTION
The Magellanic System is essential to our understanding of the ongoing formation and evolution of the Local Group. It consists of the two closest massive dwarf galaxies to the Milky Way (MW), the Large and Small Magellanic Clouds (LMC and SMC), and the Magellanic Stream, a massive network of gaseous filaments trailing behind the Clouds (see D'Onghia & Fox 2016 for a review). Significant theoretical and computational effort has gone into furthering our understanding of the origin and formation of the Magellanic Stream, beginning with simple analytical models (e.g. Fujimoto & So-fue 1977), eventually including simple hydrodynamics and self-gravity (e.g. Moore & Davis 1994;Gardiner & Noguchi 1996). Modern high-resolution simulations include live N-body models for both the LMC and the SMC as well as self-consistent hydrodynamics with radiative cooling and star formation (Besla et al. 2012, hereafter B12;Hammer et al. 2015;Pardy et al. 2018, hereafter P18;Wang et al. 2019). Most recently, Lucchini et al. (2020, hereafter L20) showed that including a Magellanic Corona of warm gas surrounding the LMC and SMC can explain the ionized gas component of the Magellanic Stream (Fox et al. 2014).
However, as these models have improved, one piece of the puzzle has remained unconstrained: the exact past orbits of the LMC and SMC. Proper motion (PM) measurements for the Magellanic Clouds have become very precise (Kallivayalil et al. 2013;Zivick et al. 2018), but the total mass of the MW and the LMC are still imprecisely known (Bland-Hawthorn & Gerhard 2016;D'Onghia & Fox 2016). The PM measurements favor a first-infall scenario (Besla et al. 2007), which is supported by the LMC's wake in the dark matter (DM) distribution of the MW halo (Conroy et al. 2021). The largest uncertainties in the orbits of the Clouds come from hydrodynamical effects including ram pressure and tidal energy losses, which are difficult to include in analytical integrators. Given the recent indications for the Magellanic Corona (see L20), as well as the need to include the MW's hot circumgalactic medium (CGM), these hydrodynamical effects will play a significant role in the orbital history of the Clouds. Upon the inclusion of the gaseous halos, we find that the Clouds can survive fewer recent interactions than previously thought (B12) if they are to remain separated at the present day. Additionally, evidence for a Magellanic Group (e.g. D'Onghia & Lake 2008; Nichols et al. 2011) suggests that, being loosely bound, its two largest members should expect only a couple of direct interactions within the past ∼5−7 Gyr.
In this Letter we explore the large-scale structure and location of the Stream resulting from an alternate firstpassage interaction history between the Clouds as motivated by the existence of the Magellanic Corona. Our new simulations are consistent with the observed PMs of the Clouds and have dramatic implications on the 3D location of the Magellanic Stream. In Section 2 we outline the methods and initial conditions used in our simulations. In Section 3 we discuss our main results, and in Section 4 we dissect the significant outcomes and implications of the model.

METHODS
For this work, we used the gizmo massively parallel, multiphysics code (Hopkins 2015;Springel 2005). gizmo employs a Lagrangian meshless finite-mass (MFM) hydrodynamics scheme that is ideal for simulations with large bulk velocities and large dynamic ranges in density. The MFM scheme provides the ability to track individual fluid elements while still capturing Kelvin-Helmholtz instabilities and shocks (Hopkins 2015). We used adaptive gravitational softening lengths for gas particles, and softening lengths of 350 and 100 pc for the DM and stellar components, respectively. Additionally, the default cooling (see Appendix B of Hopkins et al. 2018) and star formation (Springel & Hernquist 2003) routines in gizmo were included.
We do not include metal-line cooling, time-variable ionizing radiation from the MW disk, or UV background radiation in our model, because although these complex mechanisms would influence the Stream's thermal state and ionization level, we expect they would not affect its location, which is the focus of this Letter. Following L20, we assume the cold material (the H I Stream) originates from the disks of the LMC and SMC (the Magellanic ISM), whereas the warm ionized material originates from the Magellanic Corona.

Initial Conditions
The simulation presented here contains the same components as those used in L20, except we have added a live DM halo for the MW. In brief, those components are galaxies with stellar and gaseous exponential disks embedded in live Hernquist profile DM halos following the methods outlined in Springel et al. (2005), and gaseous coronae around the LMC and MW. Table 1 outlines the parameters used to generate the initial setup of the simulation.
Based on recent measurements (Peñarrubia et al. 2016;Erkal et al. 2019), and the mounting evidence for a Magellanic Group (e.g. D'Onghia & Lake 2008; Nichols et al. 2011), we include a high-mass LMC in our model. It has a total DM mass of 1.8 × 10 11 M . Our SMC is consistent with previous works (P18; L20) with a DM mass of 1.9 × 10 10 M . We use the MW model of D'Onghia & Aguerri (2020) with a gaseous component added to the disk. The total number of particles used for each galaxy is given in row 10 of Table 1, and leads to masses per particle ranging from 1.8 × 10 5 to 7.0 × 10 5 M for DM, 4.2 × 10 3 to 1.6 × 10 4 M for stars, and 4.3 × 10 3 to 1.8 × 10 4 M for gas.
Following L20, we included a Magellanic Corona around the LMC and SMC by extracting the radial density profile from an LMC analog in the Auriga cosmological simulations (Grand et al. 2017), although we used slightly different selection criteria. Specifically, we selected all gas with T > 2.5 × 10 4 K, the lowest temperature at which the galaxy's disk gas was excluded from the selection. The corona was then added around our LMC with its total mass scaled by the ratio between the Auriga LMC analog's total mass and our LMC's total mass, leading to a mass of 8.3 × 10 9 M in the Magellanic Corona and a density of 10 −4 cm −3 at 50 kpc from the LMC. Finally the corona's temperature was set to 2.4 × 10 5 K, the expected virial temperature of the LMC. These values are also consistent with other recent cosmological simulations (Hafen et al. 2019;Jahn et al. 2021). Run in isolation for 4 Gyr, the LMC and Magellanic Corona remain stable.
Around the MW, we initially included a gaseous corona following the "fiducial" density profile in Salem Notes. The resultant galaxies have rotation curve peaks of ∼240, ∼120, and ∼65 km s −1 at 12, 8.5, and 5 kpc for the MW, LMC, and SMC respectively. a The DM concentration parameter, c, is defined as the ratio of the virial radius, R vir , to the Navarro-Frenk-White (NFW) scale radius, Rs: R vir = cRs. The NFW scale radius is then converted into a Hernquist scale radius with a = Rs 2 (log(1 + c) − c/(1 + c)). b Proper motions given as (µ W , µ N ) following the convention in Kallivayalil et al. (2013). et al. (2015): a β-profile with n 0 = 0.46 cm −3 , r c = 0.35 kpc, and β = 0.559. The maximum density was capped within 13 kpc from the Galactic Center where the corona overlapped with the gaseous disk. Additionally, the profile exponentially declines for r > r vir = 166 kpc. This initial simulation was unable to reproduce the velocity profile of the Stream (see Section 3), so in the final simulation presented here, the total mass was increased by a factor of two to 4 × 10 10 M , which solved the kinematic discrepancy. This increased the MW halo density at 50 kpc from 1.1 × 10 −4 to 2.0 × 10 −4 cm −3 , still consistent with current data (e.g. Anderson & Bregman 2010;Li & Bregman 2017). We acknowledge that the Galactic and Magellanic Coronae included in these simulations are simplified when compared with the complex multiphase intricacies known to exist in circumgalactic media (van de Voort et al. 2019).

Orbits
The orbits used in previous simulations of Magellanic Stream formation were established without the Galactic or Magellanic Coronae included (B12; P18). While the Cloud orbits will clearly require modification due to the changes in total masses of the galaxies due to these added components, the increased friction and ram pressure that the Clouds experience as they move through these media also play a significant role. Therefore, to match the observed positions and velocities of the Magellanic Clouds while including the Galactic and Magellanic Coronae in our simulations, we need to determine an alternate orbital history with fewer recent interactions between the Clouds.
To do this, we analytically integrated the orbits of the LMC, SMC, and MW backward in time starting from their present-day observed positions and velocities. We included radially extended Hernquist DM halos including dynamical friction and represented the stellar and gaseous disks as point particles. By varying the presentday velocities within the observed errors, we obtained a suite of 1458 possible orbits. We filtered these to select first-passage orbits with multiple interactions between the Clouds. These orbits are generally consistent with the results found with previous analytic integrations af-ter accounting for differences in LMC mass and live versus static MW DM halos (Besla et al. 2007;Kallivayalil et al. 2013;Garavito-Camargo et al. 2019). We then chose 10 of these orbits sampling the parameter space with varying morphologies to run in full N -body hydrodynamical simulations. Due to the effects of tidal stripping, ram pressure, and friction, the initial conditions required modification to match the present-day observations of the Clouds. For 9 of the 10 chosen orbits, we were unable to modify them such that they reproduced the present-day observations. However, after a few iterations on one of the orbits, we found a solution that could match the Clouds' observed kinematic properties, which we present in this Letter. The initial positions and velocities used for the orbits in this simulation are given in rows 11 and 12 of Table 1.
Given the backward-integrated orbits, we had to choose when to begin the N -body simulation. We chose the apocenter between the Clouds after their second interaction (3.5 Gyr ago). While the Clouds would have another encounter if we integrate these orbits further back in time, this interaction would occur 7 Gyr ago (5.5 Gyr before their next interaction) at a distance of 0.9 Mpc away from the MW. Any gas stripped in this interaction would be tidally thrown out to great distances from the Clouds and become too diffuse and distant to contribute to the Stream today (see material stripped from the first interaction in Figure 2 of P18).

Comparison with Previous Orbits
In the previous orbital model of B12 and P18 (B12's "Model 2"), the LMC and the SMC experience three interactions in isolation (without any MW influence) over ∼6 Gyr. They are then rotated into the correct orientation, placed at the virial radius of the MW, and allowed to continue evolving until they reach their present-day positions (∼1 Gyr). During this infall they experience an additional direct collision that forms the Magellanic Bridge.
In the model presented here, all the interactions between the Clouds and the MW occur in a single, selfconsistent simulation lasting 3.5 Gyr. The LMC and the SMC have two interactions, the second of which has a very low impact parameter and forms the Magellanic Bridge. The main differences between these two orbital models are (see Figure 1b): • the number of interactions (4 versus 2), • the length of the simulation (7 versus 3.5 Gyr), • the maximum separation between the Clouds (100 versus 150 kpc), and • the sense of the SMC's orbit around the LMC (see Section 4 and Figure 4).
While there have been many other proposed orbital models of the interactions between the Clouds (e.g. Růžička et al. 2010;Diaz & Bekki 2011;Guglielmo et al. 2014), none have looked at this explicit combination, especially when considering a first-passage scenario with the inclusion of an MW CGM and the Magellanic Corona.

RESULTS
The simulation ran for 3.46 Gyr, when the positions and velocities of the Clouds matched current observations just after their first pericentric passage. When the LMC and SMC are at their present-day sky positions in our model, they are at distances of 52.1 ± 1.7 kpc and 78.0 ± 7.8 kpc respectively (with errors calculated as err res discussed in footnote c of Table 1; err are negligible), with proper motions and radial velocities listed in Table 1. The kinematics of the Clouds are fully consistent with the observed values: within 1σ for the proper motions and radial velocity for both the LMC and the SMC (with the exception of µ W for the LMC at 1.18σ). Additionally, the relative velocity between the Clouds (66 km s −1 ) matches within 2σ.
The present-day disk gas masses 1 for the LMC and SMC are 4.4 × 10 8 and 3.3 × 10 7 M respectively (with peak column densities of 10 21.8 and 10 21.4 cm −2 ). A full exploration of the structure of the Clouds themselves will be performed in an upcoming paper with more complete star formation and feedback routines (Lucchini, S. et al. 2021, in preparation).
Throughout their 3.5 Gyr interaction history, a trailing Stream is formed through tidal interactions between the Clouds. The orbital paths of the LMC and SMC along with the resultant Stream at the present day are shown projected onto the y-z plane relative to the MW disk in Figure 1a. The LMC and SMC experience two close encounters shown as minima in Figure 1b. Their first interaction, 1.4 Gyr ago with an impact parameter of 9.9 kpc, provides the tidal forces necessary to strip material from the SMC to create the bulk of the H I Stream. Their second interaction, 295 Myr ago, has a significantly lower impact parameter of 3.0 kpc, and this direct collision forms the Magellanic Bridge. The present-day position of the Magellanic System can also be seen in 3D in Figure 2. To test the dominant stripping mechanism in our model, we ran two additional simulations: one consisting of just the MW, its hot corona, and the SMC (without the LMC), and another with the MW, its hot corona, and the LMC with the Magellanic Corona (without the SMC). When the SMC alone passes through the MW's CGM, we find only 8% of the total gas mass stripped from the full model. When the LMC alone passes through, we find negligible stripping of gas from the LMC, but the Magellanic Corona still sees 95% of the stripping of the full model. Therefore, we conclude that the neutral Stream is stripped through tidal forces (consistent with previous findings; e.g. Salem et al. 2015), whereas the ionized component of the Stream (originating in the Magellanic Corona) is stripped mostly through ram pressure against the MW hot corona, even during the first passage.
As in L20, the bulk of the mass of the Magellanic Stream is composed of ionized gas originating in the Magellanic Corona. However, the distribution of the ionized Magellanic Corona gas on the sky is substantially different in our new model (compare Figure 3(c) with L20 Figure 2(a)), but the neutral Stream's appearance on the sky in our model is generally consistent with previous models and observations (Figure 3(a), (b), (d)). While the simulated Stream is longer than observed and slightly offset spatially, the morphology of the Stream in this model more accurately reproduces the turbulent, filamentary nature of the data. Due to its interactions with the Magellanic Corona and the MW's hot CGM, instabilities fragment and distort the Stream, leading to a significant improvement in its appearance and morphology when compared to models that do not include these gaseous components (e.g. P18).
The velocity profile of our simulated Stream also matches observations (see Figure 3e). Previous models found a velocity gradient too shallow when compared to the data (L20), but this is resolved in our new model. As stated in Section 2.1, the density profile and total mass of the MW hot corona were increased by a factor of 2 over the Salem et al. (2015) values to better match the velocity gradient along the Stream. If the MW's hot halo is not massive enough, the stripped gas from the Clouds is accelerated toward the MW and in some cases ends up with greater velocities than the LMC and SMC themselves. A higher gas density around the MW provides the ram pressure forces to slow down the trailing Stream to match the observed velocity gradient. This is in contrast to previous works that placed upper limits on the MW coronal mass such that a Leading Arm can form (Tepper-García et al. 2019;L20).
Additionally, as seen in Figure 3d, this model does not self-consistently reproduce the Leading Arm gas. This is because the Clouds only have two close encounters so there is not enough time for gas to be tidally thrown ahead of the Clouds in their orbits. This is in contrast to previous studies (Tepper-García et al. 2019;L20) where the lack of a Leading Arm was due to the MW hot coronal density being too high. In this new model, even without an MW hot corona, a Leading Arm is not formed. The true nature of the Leading Arm is one of the biggest outstanding questions in the Magellanic  (Putman et al. 1998) and the metallicities (although they vary with position; Fox et al. 2018), several works have proposed alternative, non-Magellanic sources. Hammer et al. (2015) and Tepper-García et al. (2015) suggested the Leading Arm structures could be remnants from dwarf spheroidal satellites of the MW whose gas has been stripped from the MW hot corona, and a non-Magellanic explanation for the Leading Arm remains a possibility.

System (D'Onghia & Fox 2016). While a Magellanic origin is supported by the kinematics
The most notable implication of this new model is that the Stream is significantly closer to us than previously thought (Figure 3f). While some past models have predicted a close Stream (via multiple passages around the MW; e.g. Moore & Davis 1994;Diaz & Bekki 2012), all previous first-passage orbits have resulted in the Stream flowing behind the Clouds out to distances of 100 − 200 kpc or greater (see Figure 4e; B12; P18; L20). The first-passage model presented here forms a tidal Stream that reaches as close as ∼20 kpc away from the Sun with a column density-weighted average distance of 24.7 kpc (between Magellanic longitudes of −25 • and −150 • ). See Section 4 for an in-depth discussion of this finding.
We emphasize that in all our simulations that formed a Stream (7 of the 10 selected from the backward integration, see Section 2.2), that Stream was <50 kpc away from us. While the distances and kinematics of the Clouds, and the Stream morphology, vary greatly between these individual simulation runs, the finding that the Stream remains nearby is a robust prediction of our first-infall models. Additionally, we have run several convergence tests at various numerical resolutions, and despite minor differences on small scales, all runs produced accurate positions and velocities for the Clouds, and all predicted a nearby Stream with total ionized and neutral masses consistent with observations.

DISCUSSION AND CONCLUSIONS
The new first-passage interaction history of the Magellanic Clouds presented here leads to a dramatically different 3D spatial positioning of the Stream than previous models. Previous first-passage, tidal simulations led to a Stream extending away from the Clouds out to distances upward of 200 kpc, whereas our new model results in the Stream angling up toward the MW reaching as close as 20 kpc to the Sun. While there are many differences between previous models and the model presented here (see Section 2.2.1), two differences in particular lead to this dramatic shift in positioning of the Stream: 1. a qualitative difference in the SMC's orbit around the LMC, and shows the neutral gas originating in the disks of the LMC and SMC compared to the observed data from Nidever et al. (2010) shown in contours (black, gray, and white correspond to 10 21 , 10 20 , and 10 19 cm −2 ). The centroids of the LMC and SMC stellar disks are marked by white circles. Panel (e) shows the local standard of rest (LSR) velocity gradient along the Stream with data shown as contours (Nidever et al. 2010). Panel (f) shows the line-of-sight distance to the gas in the Stream along its length with the centroids of the Clouds marked with circles. Note that the bulk of the simulated Stream is significantly closer to us than the Magellanic Clouds are (∼20 kpc vs. ∼60 kpc). 2. the inclusion of the Galactic and Magellanic Coronae.
First, the orientation of the SMC's orbit around the LMC is qualitatively different than the most recent previous models (B12; P18; L20). To compare these different orbital histories, we need to analyze the orbits in a consistent coordinate system, so we will discuss the relative motion of the SMC around the LMC when viewed projected onto the y-z plane as defined in Figure 1 (relative to the MW disk, which is in the x-y plane with the Sun located at (x, y, z) = (−8.3, 0, 0.027) kpc) 2 . When viewed from this perspective, our model has the SMC on a counterclockwise orbit around the LMC, whereas in the B12 model, the SMC rotates around the LMC clockwise. In the clockwise orbit, the Stream is tidally thrown out in the +y direction with a velocity in the −z direction, leading to it stretching away from the MW disk. Whereas in the counterclockwise orbit, the Stream is still tidally stripped in the +y direction, but its velocity is in the +z direction. This leads to a Stream angled up in the +y and +z directions, resulting in low line-of-sight distances to the Sun. A schematic of these two orbital orientations is shown in Figure 4(a)-(d).
Second, the inclusion of the Galactic and Magellanic Coronae is crucial. While previous parameter space searches of orbits for the Clouds would possibly have explored this orbital configuration in the past (e.g. Růžička et al. 2010;Guglielmo et al. 2014), none included the Galactic and Magellanic Coronae. These are key elements as the ram pressure and friction from the MW's CGM are able to "push" the Stream into its presentday position and velocity. Without the MW's gaseous halo, the Stream would collide with the MW disk before the Clouds reach their present-day positions. Moreover, the Magellanic Corona is required to shield the neutral Stream from the intense forces and pressures as it moves through high-density regions of the Galactic corona L20.
A Stream reaching ∼20 kpc has a number of implications. First, the total observed mass (neutral and ionized) of the Stream would be reduced. Total mass estimates of the Stream have previously assumed the Stream has a similar distance as the Clouds (55 kpc), but they include a scale factor of (d/55 kpc) 2 which equals 0.2 for d = 24.7 kpc (the column density-weighted average distance in the simulation). This leads to new values 2 Note that this is the same perspective as Figure 3 in B12 but is not consistent with the coordinate systems in B12's Figure 2, L20's Extended Data Figure 2, or P18's Figure 2. Because of rotations performed before the Clouds fall into the MW potential in these models, the y-z perspective discussed above is approximately equivalent to mirroring these figures across the y-axis.
of 9.7 × 10 7 M of neutral gas (Brüns et al. 2005) and 4×10 8 M of ionized gas (Fox et al. 2014), although this ionized gas mass should be considered as a lower limit as the spatial extent of the Stream's ionized phase may be significantly greater than the area on the sky that has currently been explored, and the ionized gas may be multiphase (Fox et al. 2014). Summing the masses of all gravitationally unbound particles that fall in the region of the trailing Stream in our simulation gives values of 2.0×10 8 M (neutral) and 3.2×10 9 M (ionized) for our model of the Stream, in good agreement with the observations. Second, as shown in Figure 4f, the stellar component of the Stream is also nearby, with stars predicted at d 20 kpc (with a mean surface brightness of 31 mag arcsec −2 ). Previous works have predicted that any stars associated with the Magellanic System are at large distances (B12; P18), but our new predictions suggest that continued searches for the Stream's stellar component are worthwhile (Zaritsky et al. 2020). Third, the interaction between the Stream and the MW CGM will be enhanced due to its proximity, because of the higher MW corona density. At a closer distance, the Stream would be closer to pressure equilibrium with the MW CGM, which could help explain its multiphase nature (Wolfire et al. 1995) and the high number of "head-tail" clouds seen in high-resolution H I observations . It could also lead to dramatically shorter lifetimes for the Stream in the future (Murali 2000;Bland-Hawthorn et al. 2007). Indeed, due to its angle of approach with respect to the MW, the Stream may even directly collide with the MW disk within the next ∼50 Myr. This moves the timescale of gas accretion onto the disk up by a factor of ∼10 from previous predictions. The new distance to the Stream therefore implies a factor of ∼8 less gas accreting onto the MW ∼10 times earlier than previously thought, resulting in approximately the same infall rate of ∼4-7 M yr −1 as derived before (Fox et al. 2014). Fourth, a closer Stream should be subject to a more intense UV radiation field from the MW and hence significantly brighter in Hα. This could explain the high observed Hα emission from the Stream (Bland-Hawthorn et al. 2013Barger et al. 2017); at a distance of 20 kpc, the Galactic UV background could lead to Hα emission as high as 150 − 300 mR (Tepper-García et al. 2015), in agreement with observed levels along most of the Stream. However, this increased radiation still cannot explain the extremely high Hα emission seen in the region of the Stream under the South Galactic Pole, which may indicate a recent Seyfert flare from the Galactic Center (Bland-Hawthorn et al. 2019). The enhanced radiation field intensity at d=20 kpc would also affect the ionization level inferred from UV metal-line studies of the Stream (Fox et al. 2014; however, the ionized/neutral ratio of ∼3:1 is distance independent since both the neutral and ionized masses scale as d 2 . There are currently no observational distance constraints on the Stream, although parts of the Leading Arm have constraints of < 20 − 30 kpc (McClure-Griffiths et al. 2008;Price-Whelan et al. 2018;Antwi-Danso et al. 2020). To test our prediction of a 20 kpc Stream, UV or optical spectroscopic studies could be performed to look for absorption at Magellanic velocities toward distant MW halo stars, such as blue horizontal branch stars with distances from Gaia. Lehner & Howk (2011) found no UV absorption at Magellanic velocities when looking at 28 halo stars out to z-distances of 12.6 kpc, but searches to larger distances using new stellar catalogs are needed. This offers a pathway for confirming our prediction of a nearby Stream.