The Magellanic System: the puzzle of the leading gas stream

TL;DR

This study examines the formation of the Magellanic Stream and Leading Arm, revealing that the hot Galactic corona inhibits the Leading Arm's formation, challenging previous tidal models and suggesting alternative origins or properties of the corona.

Contribution

The paper introduces a Galactic model with a weakly magnetized, spinning hot corona, showing its impact on the formation of gaseous features around the Magellanic Clouds.

Findings

The Magellanic Stream's location and extension are consistent with previous models.

The hot corona inhibits the formation of the Leading Arm in the models.

Two alternative explanations for the Leading Arm are proposed, involving corona density or a different origin.

Abstract

The Magellanic Clouds (MCs) are the most massive gas-bearing systems falling into the Galaxy at the present epoch. They show clear signs of interaction, manifested in particular by the Magellanic Stream, a spectacular gaseous wake that trails from the MCs extending more than 150 degree across the sky. Ahead of the MCs is the "Leading Arm" usually interpreted as the tidal counterpart of the Magellanic Stream, an assumption we now call into question. We revisit the formation of these gaseous features in a first-infall scenario, including for the first time a Galactic model with a weakly magnetised, spinning hot corona. In agreement with previous studies, we broadly recover the location and the extension of the Stream on the sky. In contrast, we find that the formation of the Leading Arm -- that is otherwise present in models without a corona -- is inhibited by the hydrodynamic interaction with the hot component. These results hold with or without coronal rotation or a weak, ambient magnetic field. Since the existence of the hot corona is well established, we are led to two possible interpretations: (i) the Leading Arm survives because the coronal density beyond 20 kpc is a factor of 10 lower than required by conventional spheroidal coronal x-ray models, consistent with recent claims of rapid coronal rotation; or (ii) the `Leading Arm' is cool gas trailing from a frontrunner, a satellite moving ahead of the MCs, consistent with its higher metallicity compared to the trailing stream. Both scenarios raise issues that we discuss.

Figures

Tables

Table 1. Table 1: Structural parameters of the Magellanic Clouds. The parameters values, with exception of r tr subscript 𝑟 tr r_{\rm tr} and Z 𝑍 Z , are adapted from Pardy et al. ( 2018 , their model 9:1) . Column headers are as follows: M t subscript 𝑀 t M_{\rm t} , total mass ( 10 9 superscript 10 9 10^{9}  M ⊙ ); r s subscript 𝑟 s r_{\rm s} , scalelength (kpc); r tr subscript 𝑟 tr r_{\rm tr} , truncation radius (kpc); N p subscript 𝑁 p N_{\rm p} , particle number ( 10 4 superscript 10 4 10^{4} ); Z 𝑍 Z , gas metallicity (Z ⊙ ).

	Profile	$M_t$	$r_s$	$r_{\mathrm{tr}}$	$N_p$	$Z$
LMC
DM halo${}^a$	NFW	175	13	75f	100	–
Stellar disc${}^b$	Exp${}^e$	2.5	2.4 d	9.6 g	60	–
Gas disc	Exp${}^c$	2.0	9.5	38g	100h	0.2
SMC
DM halo${}^a$	NFW	15	3.8	28f	11	–
Stellar disc${}^b$	Exp${}^e$	0.8	1.2 d	4.8g	6	–
Gas disc	Exp${}^c$	2.0	4.7	19g	36h	0.1

Table 2. Table 2: Overview of Galaxy models.

Name	DM halo	Corona	Peak $v_c$a	$B_c$b
Static DM halo	Static	–	–	–
Live DM halo	Live	–	–	–
Non-rotating corona	Live	Yes	0	0
Slow corona	Live	Yes	70	0
Fast corona	Live	Yes	135	0
Magnetised coronac	Live	Yes	0	0.1

Equations

$${\sigma }_r^2(r)=\frac{1}{{\rho }_{\mathrm{c}}(r)}{\int }_r^{\infty }{\rho }_{\mathrm{c}}(r^{\prime })\frac{\partial \Phi (r^{\prime })}{\partial r^{\prime }}d{r}^{\prime },$$

$$e\equiv \frac{1}{\gamma -1}\frac{k{T}_{\mathrm{c}}}{\mu m_p}=\frac{1}{(\gamma -1)}({\sigma }_r^2-v_{\mathrm{c}}^2).$$

$$P_{\mathrm{ram}}=\frac{1}{2}{C}_{\mathrm{D}}\text{leavevmode}{\rho }_{\mathrm{c}}{\left({\vec{v}}_{\mathrm{g}}-{\vec{v}}_{\mathrm{c}}\right)}^2.$$

$$\left|{\vec{v}}_{\mathrm{g}}-{\vec{v}}_{\mathrm{c}}\right|\ge \left|v_{\mathrm{g}}-v_{\mathrm{c}}\right|.$$

$$n_{\mathrm{d}}(R,z)=n_{\mathrm{d},0}\exp \left[-R/R_{\mathrm{d}}\right]\exp \left[-z/z_{\mathrm{d}}\right],$$

$$n_c(r)=n_{\mathrm{c},0}{\left[r/r_c\right]}^{-3\beta },$$

Full text

††pagerange: The Magellanic System: the puzzle of the leading gas stream–B

1 Introduction

The Magellanic System is the most important source of cold gas accretion within the Galaxy at the present epoch. The total gas mass of the Magellanic Stream (or simply, the Stream Wannier & Wrixon, 1972; Mathewson et al., 1974), the Magellanic Bridge (Brüns et al., 2005), the Magellanic Clouds and the Leading Arm (Lu et al., 1998; Putman et al., 1998), is estimated at $2\times 10^{9}\leavevmode\nobreak\ {\>\!{\rm M}_{\odot}}\,[r/55\,{\>\!{\rm kpc}}]^{2}$ for a median cloud distance $r$ (Fox et al., 2014), significantly higher than the mass of all high-velocity clouds (HVCs) around the Galaxy together ($\sim 10^{8}\,{\>\!{\rm M}_{\odot}}$; Putman et al., 2012). If the gas accretion time onto the Galaxy is fixed by the MC merger time ($\sim 2.5$ Gyr; Cautun et al., 2018), the gas inflow rate is roughly $1\leavevmode\nobreak\ {\>\!{\rm M}_{\odot}}\leavevmode\nobreak\ {\>\!{\rm yr}^{-1}}$, of the same order as the present-day star formation rate of the Galaxy (e.g. Robitaille & Whitney, 2010).

It is however unclear whether and over what timescale the gas will settle onto the Galactic disc. Full simulations of gas-stripped satellites within a realistic Galaxy model that includes a hot corona show that stripped gas may rain onto the Galactic disc on timescales well below 1 Gyr (Tepper-García & Bland-Hawthorn, 2018b). On the other hand, observations of recombination emission (H$\alpha$; Weiner & Williams, 1996; Reynolds et al., 1998; Putman et al., 2003b; Madsen, 2012) and ultra-violet (UV) absorption line measurements (Fox et al., 2013; Richter et al., 2013; Barger et al., 2017) indicate that mass budget of the Stream is dominated by ionised gas. A high ionisation fraction may be common to most HVCs (Lehner & Howk, 2011) in part due to the Galactic UV emission (Bland-Hawthorn & Maloney, 1999, 2001), particularly for HVCs close to the disc. For clouds far from the disc, in the absence of alternative ionisation mechanisms (e.g. stellar radiation; Ostheimer et al., 1997), the observations are best explained by constrained hydrodynamic models of infalling gas shock-ionised by the interaction with the Galactic corona (‘shock cascade’ model; q.v. Tepper-García et al., 2015). These models indicate that the gas is ionised on a timescale of order a few 100 Myr and will likely evaporate in the corona before reaching the star-forming disc (cf. Heitsch & Putman, 2009).

The interpretation of the data is complicated by the observation that the Galactic corona is magnetised, with mean field values spanning the range $\sim 0.1-10$ $\mu$G (Sun & Reich, 2010) depending on Galactocentric radius, as is common for spiral galaxies (Beck, 2016). An ambient magnetic field – neglected in all models of the Stream to date – has been long known to provide a mechanism that suppresses the onset of instabilities that ultimately drive the destruction of the infalling gas (Mac Low et al., 1994), thereby increasing its survival chances (for a discussion, see Bland-Hawthorn, 2009).

The Galactic corona, which has been accreted and has evolved over billions of years (van de Voort et al., 2011), plays a crucial role in the evolution of gas falling onto the Galaxy. Thus, a full account of its properties is required if we are to arrive at a solid understanding of the Magellanic System. A simple model of the formation of the Magellanic System is based on the paradigm referred to as ‘first-infall scenario’ (Besla et al., 2012).111For a recent review see D’Onghia & Fox (2016). However, this model relies on a crude description of the Galaxy; in particular, it neglects the presence of the circumgalactic hot corona.

Earlier models addressed the interaction between gas stripped away from the MCs and the Galactic corona, although using simpler assumptions. Following a purely analytic approach, Meurer et al. (1985, see also ) calculated the effect of drag from the Galactic corona on the gas in the Magellanic Bridge, and concluded that this could explain the formation of the Stream. This work was followed by Gardiner (1999), and extended by Diaz & Bekki (2011), using pure $N$-body models. Both groups modelled the gas drag exerted on the infalling gas by a diffuse, uniform background by adding a drag term to the equation of motion of the otherwise collisionless simulation particles. The first self-consistent $N$-body, hydrodynamical model was presented by Mastropietro et al. (2005), who modelled the interaction between a hot hydrostatic corona around the Galaxy and the LMC, but ignored the presence of the SMC.

Most recently, Hammer et al. (2015) revisit an idea first expressed by Sawa & Fujimoto (2005), and argue that the trailing stream and the ‘Leading Arm’ arises from of a convoy of baryon-dominated tidal dwarf systems, which includes the MCs, that lost gas along their orbit as a result of ram pressure stripping by the intergalactic medium and the Galactic corona. Here, we reconsider this idea as one of the plausible scenarios explaining the leading stream gas, i.e. the gas leads the MCs but does not emerge from it dynamically.

Even with their adopted simplifications, the previous studies were able to show that the properties of the Magellanic System are sensitive to the presence of an external corona along its orbit. Here we provide a more sophisticated perspective on the evolution of the Magellanic System, in particular of the Stream and the Leading Arm, using an improved model for the Galaxy, which includes – for the first time – a magnetised Galactic corona with a range of properties consistent with observations. This is the core of our paper.

We are not attempting to reproduce the properties of the Magellanic System in detail, e.g. the precise location and kinematic properties of the MCs, or the properties of the gas structures such as its total mass, its ionisation state, or the kinematic and structural properties of either the trailing stream or the leading arm. For all the partial success of others in doing so (q.v. D’Onghia & Fox, 2016), we maintain that it is futile to pursue such an ambitious task given that the present-day position and kinematics of the MCs, the Stream and the Leading Arm do not unambiguously constrain the required initial conditions, such as the total masses of the MCs and their structure, the inclination of the SMC with respect to its orbit relative to the LMC, their spin, or their full three-dimensional (3D) position and velocity at infall which depend on the unknown Galactic potential. The initial conditions rely on a large number of model-dependent parameters (for a discussion see Hammer et al., 2015). Since numerical experiments are costly, it is prohibitive to explore in full the available parameter space looking for an optimal solution. Techniques such as genetic algorithms can in fact be used to perform the search more efficiently (Guglielmo et al., 2014), but even here the available parameter space needs to be significantly constrained beforehand.

Our approach rests on adapting a recent model (Pardy et al., 2018, dubbed ‘9:1’) – based on Besla et al. (2012, their model 2) – that broadly reproduces some of the observed properties of the Stream and the Leading Arm. We systematically modify this model to investigate the effect of several factors (dynamical friction, gas drag, confinement by dark matter and magnetic fields) on the evolution of the gas structures associated to the Magellanic System, with particular emphasis on the Leading Arm, in a series of controlled and self-consistent numerical experiments. We argue that the relevant physical processes at play should be captured correctly in our models. Thus, our results should be applicable to other similar systems.

This paper is organised as follows. In Section 2, we describe our adopted base model and variations thereof. In Section 3, we present the results for each of our models. We discuss the implications of our findings in Section 4, and finish with some concluding remarks in Section 5.

2 The First-Infall Scenario Revisited

2.1 The base model

Our fiducial model follows the 9:1 model by Pardy et al. (2018). This model is constructed essentially in two steps: i) The LMC and SMC are evolved for $\sim 6$ Gyr as a binary pair in isolation starting from some prescribed initial conditions; ii) The evolved binary system is placed into a Galaxy-like potential at one Galactic virial radius ($\sim 220$ kpc), and evolved further for roughly 1 Gyr (see Appendix A for further details). The gravitational effect of the Galaxy on the MCs is approximated by a static dark matter (DM) host halo parametrised by a spherical Navarro et al. (1997, NFW) profile with a virial mass of $\sim 1.5\times 10^{12}$ M*⊙*, and a scale radius of $\sim 20$ kpc (corresponding to a concentration of 12). In contrast, the MCs are modelled as fully live, multi-component galaxies, each consisting of a stellar disc and a gas disc embedded in a DM host subhalo. The structural parameters of the MCs are summarised in Table 1.

In this model, the MCs have – by design – large initial gas fractions, implying that higher gas masses are deposited onto the structures that later form the Stream and the Leading Arm, compared to models with lower initial gas fractions. Yet the column densities along the Stream and the Leading Arm in this model are lower than observed (see Pardy et al., 2018, for a discussion). We note that our choice of simulation code (grid-based rather than SPH; see Section 3) forces upon us the need to define an arbitrary truncation radius for each model component in order for them fit within the simulation volume.222The truncation radius needs to be chosen with care, given that orbital history of the SMC around the LMC is very sensitive to their masses (dynamical friction), and the truncation radius affects the actual mass of the corresponding component. We have chosen values that most closely reproduce Pardy et al. (2018)’s 9:1 model. We refer to the base model as ‘Static DM halo’ model.

2.2 Extended models

We consider, in addition, the following stepwise improvements of the Galaxy model with respect to the starting point described above.333These modifications only affect the evolution of the Magellanic System at infall; in other words, the evolution of the binary pair in isolation is identical in all the models considered here. These ‘extended’ models are summarised in Table 2.

First, we remove the restriction of a static DM host halo and approximate the gravitational effect of the Galaxy onto the MCs using a fully live DM host halo with identical structural properties. This modification allows to study the effect of dynamical friction on the orbital history of the MCs and its impact on the formation of the Stream and the Leading Arm. This model is dubbed ‘Live DM halo’ model.

As an extension of the latter model, we include a static (i.e. non-rotating), spheroidal Galactic corona. The inclusion of this component in our model makes it possible to explore the impact of gas drag along the orbit of the MCs on the infalling gas within a realistic Galaxy model in a self-consistent manner.

Our model corona follows initially the same profile and extension as the DM host halo, adopting identical structural parameters. It has a total initial mass of $\sim 3\times 10^{10}$ M*⊙* within 250 kpc, and properties (density profile, temperature profile) in broad agreement with the properties of the Galactic corona as inferred from a number of observations (Figures 1 - 3; q.v. Bland-Hawthorn & Gerhard, 2016). It is worth noting that our coronal model is consistent with the most recent model of the Galactic hot halo by Bregman et al. (2018). We refer to this model as ‘Non-rotating corona’.

Recent kinematic models based on x-ray spectra suggest that the Galactic corona is spinning fast ($181\pm 41$ km s*-1* out to $R\approx 50$ kpc; Hodges-Kluck et al., 2016). Such a rapidly rotating corona can be accommodated within barotropic (Pezzulli et al., 2017) and baroclinic (Sormani et al., 2018) equilibrium models within the context of cosmologically motivated angular momentum distributions. Therefore, we improve upon our ‘Non-rotating corona’ model by allowing the corona to have net rotation over a range of radii, $v_{\rm c}(r)$. Thus we are able to explore for the first time the effect of the torque exerted by the corona on the gas streams associated with the Magellanic System.

We set the rotation speed of the corona according to the specific angular momentum of the Galactic DM halo, which follows roughly its cumulative mass profile (Bullock et al., 2001), adopting a fixed value for the spin parameter $\lambda$. We consider two variants: a ‘Slow corona’ model and a ‘Fast corona’ model, with peak rotation speeds in the plane of the Galaxy444In our models, the rotation speed decreases with distance from the plane Sormani et al. (cf. 2018). ($z=0$) of $v_{\rm c}\approx 70$ km s*-1* and $v_{\rm c}\approx 135$ km s*-1*, corresponding to $\lambda=0.08$ and $\lambda=0.16$, respectively These values are at the high end of the typical range for $\lambda$ found at $z=0$ in $N$-body, hydrodynamic cosmological simulations of structure formation (Zjupa & Springel, 2017). The peak rotation speed in the ‘Fast corona’ model is marginally consistent with the low end of the Galactic coronal rotation speed inferred from observations (Figure 3). A higher rotation speed would lead to a coronal temperature that is too low to be consistent with the temperature inferred from x-ray data. This is a consequence of the condition of thermal equilibrium imposed on the model corona, which naturally leads to lower mean temperatures for higher net rotation speeds (Figure 2).

The radial velocity dispersion of the gas, implied by the spherically symmetric Jeans (1915)’s equation with isotropic velocity distribution, is given by

[TABLE]

where $\rho_{\rm c}(r)=\mu\leavevmode\nobreak\ m_{p}\leavevmode\nobreak\ n_{\rm c}(r)$ describes the mass density profile of the corona, $n_{\rm c}$ is the corresponding particle density, and $\Phi$ is the total potential of the Galaxy. The specific internal energy of the gas, corrected for its net rotation speed $v_{\rm c}(r)$, is given by

[TABLE]

Thus, for a fixed $\rho_{c}$ and $\Phi$, there is a maximum rotation speed such that $T_{\rm c}$ is consistent with observations. A higher rotation speed requires either: i) abandoning the assumption of thermal equilibrium, which makes it significantly harder to model the corona; or ii) adopting a more massive DM host halo for the Galaxy, ruled out by a number of observations (e.g. Piffl et al., 2014; Kafle et al., 2014; Dierickx & Loeb, 2017).

The last extended model we consider corresponds to a weakly magnetised version of our ‘Non-rotating corona’ model, referred to as the ‘Magnetised corona’ model. We adopt a magnetic field initially uniform along the $z$ axis set initially to 0.1 $\mu$G everywhere. We are aware that this is not entirely consistent with the structure of the magnetic field around the Galaxy (e.g. Jansson & Farrar, 2012). In particular, the field is stronger (weaker) compared to the observed field far away (close) to the Galactic Centre. However, our goal is not to model the evolution of the field (see e.g. Rieder & Teyssier, 2016) but to explore the effect of an initially weak ambient magnetic field on the infalling gas with the simplest possible configuration.555We note that the magnetic pressure, $P_{B}=B^{2}/8\pi$, is negligible compared to the thermal pressure in the corona, $P_{\rm c}$, i.e. $\beta=P_{\rm c}/P_{B}\gg 1$, such that the presence of a magnetic field does not affect the thermal properties of the model corona. The use of a more realistic field structure is left for future work.

It is worth emphasising that while the gas in and around the MCs is known to be magnetised (q.v. Kaczmarek et al., 2017), we have chosen to ignore this circumstance in order to cleanly isolate the effect of a magnetised medium (the Galactic corona) on non-magnetised infalling gas clouds (the gas associated to the MCs).

For the current work, it was unnecessary to adopt more sophisticated Galaxy models (e.g. models that include a stellar bulge, a stellar disc, or a gas disc) as we have done in previous work (Tepper-García & Bland-Hawthorn, 2018a, b). While these baryonic components may generally dominate the potential close to the Galactic Centre, their gravitational effect is marginal at the mean distance of the MCs today ($r_{\rm MC}\approx 55$ kpc; Walker, 2012; Graczyk et al., 2014). Nonetheless, the gas disc may have an important effect: it may negatively affect the lifetime of gas structures that interact with it (e.g. Galyardt & Shelton, 2016). In view of our results (yet to be discussed), we opt to ignore this component for now, and leave more complex Galaxy models for future work.

2.3 Numerical setup

The initial conditions (i.e. the particle positions and their velocities; and additionally the internal energy for the gas components) for each of our models described above are generated with the dice code (Perret et al., 2014).666The sequence of geometric and kinematic transformations used to map the evolved, isolated MC binary pair onto the Galaxy is described in Appendix A. See also footnote 9.

The evolution of the MCs in isolation and its subsequent infall into the Galaxy are calculated from the corresponding initial conditions using Ramses, an $N$-body, magneto-hydrodynamic (MHD) grid code with adaptive mesh refinement (AMR; version 3.0 of the code last described by Teyssier, 2002). For the run adopting a magnetised corona, we exploit the MHD capabilities of Ramses, which implements the constrained transport, divergence-preserving (to numerical precision) MHD solver of Teyssier et al. (2006, see also ). The contribution to the overall gravitational field of every component (gaseous or collisionless) is taken into account at all times by solving the Poisson equation777The source term is given by the sum of the individual component densities at each point. on the AMR grid using the multi-grid scheme of Guillet & Teyssier (2011). The simulations of the isolated MCs and the composite system (MCs and Galaxy) are run in a cubic volume of 500 kpc and 1 Mpc length per side, respectively.

The base grid level is set at 7, corresponding to a physical spatial resolution of roughly 4 kpc (8 kpc) in our 500 kpc (1 Mpc) box. We adopt a maximum refinement level of 13 when calculating the evolution of the MCs in isolation, and 14 for the composite system. The physical spatial resolution is $\gtrsim 60$ pc in either case. The minimum cell size we adopt is thus well below the typical linear sizes of clouds along the Stream inferred from radio observations (Hsu et al., 2011).

The grid is refined (coarsened) at runtime whenever the gas mass in a cell exceeds (falls short of) $\sim 2.5\times 10^{6}\leavevmode\nobreak\ {\>\!{\rm M}_{\odot}}$ or the number of collisionless particles is larger (smaller) than 40. Also, cell refinement is applied in order for the Jeans (1915) length to be resolved at all times with at least four cells (Truelove et al., 1997).

In addition to the above refinement criteria, we have used a set of fixed, nested, spherical grids centred on the model Galaxy of increasing size (with radii ranging from 0.4 kpc to 450 kpc) and decreasing resolution (by a factor 2 from grid to grid) moving from the centre outwards.888This approach is similar to the so-called Enhanced Halo Resolution (EHR) used by Hummels et al. (2018). The innermost level has a resolution equal to our maximum adopted refinement level (14). This approach was necessary in order to roughly maintain the same resolution across the simulation volume in models with a live DM halo compared to the model with a static DM halo. For the latter, the AMR grid is refined using a density-based approach using the analytic DM density profile. In simulations with a live DM halo, the refinement strategy associated with collisionless components is governed by the number of particles within a given voxel, which depends on the actual $N$-body realisation of the DM density profile and the adopted mass resolution.999The setup files used to create initial conditions as well as those containing the details such as the hydro solver, the refinement technique, etc. used to run our simulation are freely available upon request to the corresponding author (TTG).

We include radiative gas cooling by hydrogen, helium, and heavier elements following Theuns et al. (1998). The gaseous components are all taken to be monoatomic ($\gamma=5/3$), ideal gases characterised by an hydrogen fraction ${\>X_{\textnormal{\sc h}}}=0.76$ and an initially uniform metallicity $Z$ that varies from component to component (Table 1). The metallicity of the Galactic corona is set to $Z=0.3$ Z*⊙* (Miller & Bregman, 2015). The average gas metallicity of the MCs today (as given by the abundance of $\alpha$ elements) is $Z\approx 0.3$ Z*⊙* (SMC) and $Z\approx 0.4$ Z*⊙*  (LMC; Russell & Dopita, 1992). Chemical evolution models predict that the gas phase abundances in these galaxies 1 - 2 Gyr ago were roughly 0.2 – 0.3 dex lower than today (Harris & Zaritsky, 2004, 2009). Therefore, we set the initial gas metallicity of the SMC (LMC) to $Z=0.1\leavevmode\nobreak\ {\>\!{\rm Z}_{\odot}}$ (0.2 Z*⊙*).101010We ignore the value of the gas metallicity adopted by Pardy et al. (2018). Our approach neglects the chemical enrichment of the gas via stellar feedback. A higher metallicity would enhance the cooling of the bound gas allowing the gas to sink further into the potential well of the dwarfs, making stripping more difficult. Energetic stellar feedback – which we ignore; see below – would have the opposite effect. For the time being we simply assume that these effects roughly balance out. We note that this assumption does not affect the outcome of our experiments as we shall discuss later.

In contrast to Pardy et al. (2018), we do not include star formation (SF) in our models. However, this process does not appear to have a significant impact on the overall properties of the Magellanic System in their model. Indeed, as we show in Section 3 and in Appendix A, we are able to broadly replicate their results in spite of our ignoring this process. Since we do not include SF, we neglect the associated feedback, as do Pardy et al. (2018) and Besla et al. (2012). While stellar feedback and their associated outflows in the LMC subject to the ram pressure of the Galactic corona may contribute to the formation of the Stream, its relative importance in terms of mass is minor (Salem et al., 2015; Bustard et al., 2018).

We neglect the ionising effect of UV radiation field of the Galaxy or the MCs on the gas. This is justified since we are not interested here in the ionisation state of the gas, which has negligible impact on the magnitude of gas drag (see Section 3.2). Also, the heating of the ionising radiation on the gas is expected to be small (e.g. Nichols & Bland-Hawthorn, 2011).

In a nutshell, our experiments consist of full 3D, $N$-body, magneto-hydrodynamic, grid-based (AMR), non-adiabatic, numerical calculations of the infall of the MCs into the Galaxy that take into account the self-gravity of each component.

3 The evolution of the Magellanic System

Quite generally, we find in all the models considered here that the MCs have present-day positions and velocities broadly consistent with their observed values, in agreement with Pardy et al. (2018). Similarly, we find that the Stream roughly matches its observed extension and position across the sky. However, we are led to some striking conclusions with respect to the formation of the Leading Arm.

Our main findings111111In order to keep the image file sizes manageable, all the physical space figures have been created using up to 8 (out of 14) AMR refinement levels; this implies that all the quantities in the simulation cube (gas density, magnetic field, etc.) are in fact slightly higher than shown the corresponding figures. are summarised in Figures 4, 7, and 9. There, we show the distribution of the gas associated with the Magellanic System, i.e. initially bound to the MCs, at the present epoch for each our models. Each row corresponds to a different different Galaxy model, and consists of three panels, which display, respectively, the gas distribution both in physical space projected along two orthogonal directions, and in the observed space. It is worth emphasising that the column densities shown presented here correspond to total gas densities, and thus put a strict upper limit on the expected H i column density of the gas (which we do not consider here). The simulation results in physical space are mapped onto observed space by assuming the Sun is located at $(x,y,z)=(-8,0,0)$ kpc at all times, and the North Galactic Pole points along the $z$ direction. We do not adopt the exact solar coordinates (cf. Bland-Hawthorn & Gerhard, 2016) since we are not interested here in reproducing the precise location of the Magellanic System and its components, as stated above.

In all models, we find that the total gas mass of the trailing gas stream alone is $\lesssim 10^{9}$ M*⊙*, at the low end of the range of masses of the Magellanic Stream estimated using ultra-violet absorption spectroscopy of bright background sources (Fox et al., 2014). Our result is agreement with the models by Besla et al. (2012) and Pardy et al. (2018).

We now discuss in detail the results for each of our models and their differences where relevant.

3.1 The effect of dynamical friction

As noted above, our model ‘Static DM halo’ is, by design, identical to the model 9:1 by Pardy et al. (2018). As we show here and in Appendix A, our results agree remarkably well with theirs. It is worth stressing that this agreement is not at all trivial, given the significant difference in the codes and methods used by each group both to generate the initial conditions and to calculate the evolution of the system.

The effect of dynamical friction on the orbit of the MCs can be studied by comparing the ‘Static DM halo’ and ‘Live DM halo’ models (Figure 4). Overall, we find no significant difference in the large-scale distribution of the Magellanic System between these two models. However, there are some noticeable differences in the orbital history of the MCs in the presence of a live DM host halo. In this case, the pericentric distance of the MCs are smaller, and their corresponding space velocities lower, compared to their infall into a static DM host halo (Figure 14).

In addition to generating dynamical friction, a live DM host halo responds as a whole to the gravitational interaction with a massive system such as the LMC-SMC binary pair, i.e. it is free to ‘recoil’. The recoil of the Galactic DM halo may have important consequences on the orbital history of infalling systems (Gómez et al., 2015). In the presence of a live Galaxy, as the MCs fall in, they and the Galaxy move around their common centre of mass. A fraction of the orbital energy of the MCs is absorbed by the motion of the Galaxy as a whole, and therefore the MCs move slower and sink further. We find that the recoil of the Galactic DM host halo in our models is non-negligible, on the order of 10 kpc when the MCs roughly reach their present-day position. This has a noticeable effect on the positions and kinematics of the MCs (Appendix A), but it does not affect substantially the overall position and extension of the Stream.

There is a subtlety inherent in our choice of simulation code. Like any grid-based code, Ramses does not allow the existence of a true vacuum within the simulation volume. At initialisation, each cell is ‘filled’ with a fluid defined by its density, pressure, and kinetic energy, however small. In our ‘Static DM halo’ and ‘Live DM halo’ models, we have chosen the initial density of the cells not associated to the MCs gas to a very low value, and chosen the temperature to match the mean temperature of the gas associated to the MCs at infall. The fact that we are able to broadly reproduce Pardy et al. (2018)’s results indicates that these choices have not affected our results in any significant way.

3.2 The effect of gas drag

Our model series allows us for the first time to study the impact of the Galactic corona on the infalling gas in a systematic and self-consistent manner. A key requirement of our models is that the properties of the model corona – notably its mean density distribution – be consistent with the observed properties Galactic corona during the infall of the MCs.

Given the strong differences between each of the coronal models in terms of their internal kinematics, their structure may evolve in time in a noticeable way. This is illustrated in Figure 5. When $v_{\rm c}>0$ km s*-1*, a corona that is initially spherically symmetric by design, eventually flattens along its spin axis. The effect is more severe for higher rotation speeds. Such a flattening is expected in equilibrium models (Pezzulli et al., 2017; Sormani et al., 2018).121212We caution against a misinterpretation of the word ‘flattening’. At times, the term is used to indicate a shallower, albeit spherically symmetric, density profile of the hot corona compared to some reference profile (e.g. Bregman et al., 2018). Here, the term is used to indicate that the hot corona is ‘squashed’ along its spin axis, with the result that its 3D structure is no longer spherically symmetric.

In spite of the evolution and the consequent difference in structure between the corona models, their mean density profile is comparable at all times, and does not change significantly over the course of our simulations. In Figure 6, we compare the initial (i.e. at infall) and final (i.e at the present epoch) density profiles of the corona in each of our relevant models well beyond the present-day distance of the MCs. We find that the model corona evolves away from its initial state but remains well within the available observational constraints. In particular, the coronal density drops slightly at $r>50$ kpc such that the gas drag is slightly reduced along the MCs’ orbit.

The inclusion of a Galactic corona in our simulations introduces a technical difficulty when using AMR codes. In contrast to the models that lack a corona, it is no longer straightforward to isolate (and visualise) the gas initially associated with the MCs. In order to overcome this difficulty, we attach a passive scalar (colour tracer) to the gas initially associated with the MCs prior to placing the evolved binary system at the Galactic virial radius. Specifically, we set the value of the gas tracer associated with the MCs to 20, and to 0 everywhere else. Note that these choices are arbitrary and do not affect in any way our results.

The presence of a hot corona embedded within a live DM host halo does not affect in any significant way the orbital history of the MCs compared to the model that features a live DM host halo but no corona (see Figure 15). But it does have a noticeable impact on the evolution of the infalling gas.

A visual comparison between models ‘Live DM halo’ (second row in Figure 4) and ‘Non-rotating corona’ (first row in Figure 7) readily reveals our more striking result: In the presence of a Galactic corona, the Leading Arm is absent from the final configuration that otherwise resembles the observed Magellanic System. In contrast, the Stream appears slightly more extended (by roughly 30*∘*) compared to models that lack a Galactic corona.

The absence of a leading gaseous stream is common to all the models that include a Galactic corona, regardless of the corona’s rotation speed, as can be seen by comparing models ‘Non-rotating corona’, ‘Slow corona’, and ‘Fast corona’. Note that the ‘Fast corona’ model features a gas structure that, projected on the sky, could be reminiscent of the Leading Arm. Such a feature results from the gas clouds that are decelerated by the drag force, thus falling onto lower orbits and eventually leading the MCs (e.g. Moore & Davis, 1994). However, the distribution of gas in physical space in our ‘Fast corona’ model reveals that it is not a truly leading gas stream (cf. bottom row of Figures 7 and 11).

The absence of a leading gas stream in all models that include a Galactic corona can be understood in terms of drag forces. The magnitude of the drag force, or equivalently the ram pressure $P_{\rm ram}$, of the corona on the infalling gas associated to the Magellanic System depends on the mass density of the corona $\rho_{\rm c}$ and the infall velocity of the gas $\vec{v}_{\rm g}$ relative to the velocity of the corona $\vec{v}_{\rm c}$ as

[TABLE]

The drag coefficient $C_{\rm D}$ is a parametrisation of the effective momentum transfer between the infalling gas and the corona (Benjamin & Danly, 1997).

Both $\rho_{\rm c}$ and $v_{\rm g}\equiv|\vec{v}_{\rm g}|$ increase as the MCs approach the Galaxy, and so does the destructive effect of the Galactic corona on the infalling gas. Since the MCs move on a slightly prograde orbit with respect to corona’s rotation,131313The orbit is in fact near-polar. the effective ram pressure exerted by the corona on the infalling gas is expected to be somewhat smaller in models with a faster spinning corona compared to models with a non-rotating corona.141414Both these statements are easily proven for any $\vec{v}_{\rm g}$ and $\vec{v}_{\rm c}$ such that $v_{\rm c}<v_{\rm g}$, a condition that is satisfied in our models, with aid of the reverse triangle inequality,

$\left|\vec{v}_{\rm g}-\vec{v}_{\rm c}\right|\geq\left|v_{\rm g}-v_{\rm c}\right|\,.$

In the first case, consider a fixed $v_{\rm c}$. Trivially, the ram pressure increases because $v_{\rm g}$ increases. In the second case, consider a fixed $v_{\rm g}$. Here the ram pressure decreases because $v_{\rm c}$ increases. These processes are all illustrated in Figure 8. There, we compare the evolution of the infalling gas in a model without a Galactic corona (left column), a model with a non-rotating corona (centre), and a model with a fast spinning corona (right). Each row corresponds to a different snapshot from the infall of the MCs up to $\sim 100$ Myr before the present day. The pressure of the gas is indicated by the grey scale. Overlaid are a series of contours corresponding to total gas column density of the gas associated to the Magellanic System with values in the range $\log_{10}\left[{\rm N}_{\rm gas}({\>\!{\rm cm}^{-2}})\right]=(18,20.5)$ in steps of 0.5 dex.

In the absence of a Galactic corona (left column), the gas the MCs have stripped from one another while evolving in isolation, and which is surrounded by design by a gas structure consisting of a leading and a trailing components, simply becomes tidally stretched as the system falls into the Galaxy. In the presence of a Galactic corona (central and right columns), the infalling gas is affected early on by the corona’s drag. The build-up of ram pressure ahead of the MCs and the resulting gradient along the gas as it falls in is apparent. Moving at roughly 200 km s*-1* through a stationary corona with a density ${\rm n}_{\rm c}\sim 10^{-5}$ cm*-3*  (at $\sim 200$ kpc), the leading edge of the gas experiences a ram pressure a few times $10$ cm*-3* K. The ram pressure is significantly higher at smaller distances because of the increasing ambient density. It appears somewhat lower if the corona is rotating, which is a consequence of two factors: 1) the reduced velocity difference (slightly prograde orbit); ii) the flattened density profile, which features a lower density at a given distance compared to a non-rotating corona (see Fig. 5).

At infall, the leading gas stream that otherwise gives rise to the Leading Arm in corona-free models is pushed against itself by drag forces and ‘sweeps’ past the denser gas around the MCs. The formation of the Leading Arm is hindered early on. The effect of the gas being compressed along its leading edge and ablated is consistent with the shock cascade found in our previous constrained simulations (Bland-Hawthorn et al., 2007; Tepper-García et al., 2015). In consequence, the total mass in the trailing gas stream is slightly higher in these models compared to models that lack a corona.

It could be argued that the gas structures around the MCs that lead to the formation of the leading gas stream in corona-free models may be negatively impacted when the evolved isolated binary system is placed into the model Galaxy that features a corona. To test our results against a potential numerical artefact, we have run a control model similar to the ‘Non-rotating corona’ model, but with the initial extension of the Galactic corona limited to 150 kpc (rather than $\sim 250$ kpc). In this model, the MCs and their associated gas do not interact with the corona until their reach a Galactocentric distance of $\sim 150$ kpc. Reassuringly, we find no significant difference in the final configuration in these control models compared to our models presented here. In particular, none of the control models displays a leading gas stream (not shown).

We note finally that, because the drag depends on the mass density, rather than the total particle density, of the ambient medium, models that take into account the ionisation state of the gas (which we have ignored in our models) should not yield significantly different results.

3.3 Confinement by an ambient magnetic field

The Galactic corona and disc gas have long been known to be magnetised (Beck et al., 1996), with a complex structure and a strength that generally declines rapidly away from the Galactic plane (Jansson & Farrar, 2012). Moreover, there is evidence for the presence of an enhanced, coherent magnetic field around the Leading Arm151515In the region identified as LA II north of the Galactic plane (Venzmer et al., 2012). with a magnitude along the line of sight $\gtrsim 6$ $\mu$G (McClure-Griffiths et al., 2010).

Early studies based on two-dimensional (2D) simulations of the interaction between gas clouds and magnetised media (Mac Low et al., 1994; Jones et al., 1996; Konz et al., 2002) found that the presence of even a weak ambient magnetic field can extend a cloud’s lifetime by draping around the gas and dampening the onset of boundary instabilities that eventually lead to the destruction of the cloud. This result has been confirmed by recent, more advanced 3D numerical experiments (McCourt et al., 2015; Banda-Barragán et al., 2018), which also have shown that the draped field around the infalling gas can be amplified by factors of $\sim 100$ relative to the background (Grønnow et al., 2017).

But magnetic fields may have a negative impact on infalling clouds as well. Sophisticated 3D simulations (e.g. Kwak et al., 2009; Grønnow et al., 2018) demonstrate that clouds moving through a magnetised medium experience a somewhat enhanced drag, compared to initially identical clouds moving in a non-magnetised medium. In other words, infalling gas clouds may be significantly delayed from accreting onto the Galaxy as a consequence of the Galactic global magnetic field (see also Birnboim, 2009). Thus, whether and on which timescale the infalling gas reaches the disc will depend on the competition between two opposing effects: An enhanced drag force on the one hand, and an extended cloud lifetime on the other.

Our model dubbed ‘Magnetised corona’ is designed to investigate the outcome of this competition. As detailed in Section 2.2, we set up a model where the Galactic corona is initially identical to the ‘Non-rotating corona’ model, but threaded by an initially uniform magnetic field at a level of 0.1 $\mu$G, and use it to recompute the evolution of the Magellanic System.161616Given the absence of any significant dependence on the kinematic properties of the corona discussed in Section 3.2, for simplicity we do not consider here a spinning, magnetised corona model. The result is shown in Figure 9. In agreement with previous studies, we find that the motion of the infalling gas trough the magnetised corona sweeps-up the ambient magnetic field, yielding a draped, enhanced field with a magnitude $\sim 0.1-10$ $\mu$G all around the infalling gas (see bottom panels), a range of values consistent with observations. Although we do not see evidence for an enhanced drag compared to the results of our previous models, the leading gas stream apparently does not survive in this model either.

It is fair to ask whether this result is sensitive to initial configuration and magnitude of the magnetic field we have chosen, both of which are somewhat arbitrary. As by shown previous work (Kwak et al., 2009; Grønnow et al., 2018), an initially uniform ambient magnetic field transverse to a gas cloud’s motion enhances its survival. Also, the shielding effect increases with the field strength. Therefore, our adopted initial field configuration, with an strength far away from the Galaxy that is high relative to the strength of the Galactic magnetic field at the same location, and a significant transverse component along the orbit of the MCs, should favour the survival of the infalling gas.

We thus conclude that the presence of an ambient magnetic field does not help the leading gas stream survive.

3.4 Confinement by dark matter

Finally, we explore whether the survival of the Leading Arm can be attributed to the confinement of the infalling gas by DM. As suggested by previous work (Quilis & Moore, 2001; Nichols et al., 2014; Plöckinger & Hensler, 2012; Tepper-García & Bland-Hawthorn, 2018a), a DM halo around the infalling gas may reduce the stripping effect of ram pressure by confining the baryons within its potential well. In this scenario, DM tidally stripped from the MCs ahead of their orbit could have provided the required confinement to the Leading Arm, through the Galactic corona (and eventually across the Galactic gas disc).

In order to assess whether DM confinement may have played a role in the survival of the Leading Arm, we analyse the distribution of DM associated to the MCs at the present epoch in two of our models: ‘Live DM halo’ and ‘Non-rotating corona’. The results are shown, respectively, in the top and bottom rows of Figure 10. There, we show for each of these models the gas distribution in physical space along two orthogonal projections. Overlaid we show the projected DM density initially associated to either the LMC or the SMC. The result for the other corona models, in particular the ‘Magnetised corona’, are comparable to the result of the ‘Non-rotating corona’ model, and are therefore not discussed here.

Clearly visible are, for each dwarf, a leading and a trailing DM stream. Since the LMC is initially larger, the streams associated to its DM subhalo are more extended. It is worth emphasising that the DM distribution of either dwarf galaxy is similar in both models. This is expected since DM, being collisionless by assumption, should be sensitive only to tidal forces. Since the potential is dominated globally by the Galactic DM halo and locally by the DM subhalos of the dwarfs, both of which are initially identical in both models, there is no reason for the DM distribution to be too dissimilar between models.

The DM streams around the MCs in the ‘Live DM halo’ model (top row) appear to be confining the leading gas stream, and also – although to a lesser degree – part of the trailing gas stream. However, we just mentioned that the DM streams are similar in both models. Therefore, if DM confinement were at play here, the same effect should be observed in the ‘Magnetised corona’ model (bottom row). Clearly, this is not the case.

4 Discussion

4.1 The need for ram-pressure models

The difficulties arising in ‘ram pressure models’ when attempting to reproduce a leading gas stream, with either an analytic description of the gas drag (Meurer et al., 1985; Heller & Rohlfs, 1994; Moore & Davis, 1994) or with a self-consistent hydrodynamic treatment (Mastropietro et al., 2005), were recognised early on. This circumstance led to the view that the competing ‘tidal models’ (e.g. Fujimoto & Sofue, 1976; Murai & Fujimoto, 1980; Gardiner et al., 1994) appeared more successful in explaining the formation of the Magellanic System (e.g. Putman et al., 1998). But new, accurate proper motion measurements of the MCs (Kallivayalil et al., 2006) would lead to a radical change of view. Specifically, these measurements implied that the MCs are likely on their first approach to the Galaxy (‘first-infall scenario’; Besla et al., 2007), as had been suggested much earlier (q.v. Lin et al., 1995). This implication was only reinforced by the latest measurements of the MCs’ proper motions (Kallivayalil et al., 2013).

If the MCs are truly on their first approach to the Galaxy, there has not been enough time for the tidal interaction with the Galaxy to allow for the formation of the Magellanic System. Thus it was argued (Besla et al., 2010) that the key to explaining the origin of the Magellanic System lies in the interaction between the LMC and the SMC prior to being accreted by the Galaxy. Besla et al. (2012) showed that many (but not all) of the properties of the Magellanic System can be reproduced by a model constructed based on this paradigm. Nevertheless, the result that the measurements of the MCs’ proper motions favour a first approach of the MCs onto the Galaxy over alternative orbital histories depends strongly on the Galactic potential, and it is also subject to the condition that the MCs have been a bound binary pair for a significant fraction of the Hubble time; a strong assumption indeed. In fact, others (Diaz & Bekki, 2012) have shown that the long-term connection of the MCs and their first approach to the Galaxy are unwanted assumptions in order to reproduce some key properties of the Magellanic System not reproduced by first-infall models, such as the observed spatial and kinematic bifurcation of the Stream (Putman et al., 2003a; Nidever et al., 2008).

On the other hand, circumstantial evidence for the long-term association of the LMC and the SMC is provided by a study of the global star-formation history of the LMC, which suggests that the LMC and SMC may have had coupled episodes of star formation, one occurring roughly 2 Gyr ago and the other some 500 Myr ago (Harris & Zaritsky, 2009), indicating that the MCs have been a strongly interacting pair in the past (but do not need to be bound; see Guglielmo et al., 2014). Moreover, the existence of a large number ($\sim 70$) of Galactic satellites associated with (and moving ahead of) the LMC as suggested by data of the Dark Energy Survey provides strong and independent evidence for a first passage of the LMC around the Galaxy (Deason et al., 2015; Jethwa et al., 2016). We thus consider it highly likely that the MCs are on their first approach as a binary pair.

A noteworthy caveat of Besla et al. (2012)’s model is that it ignores the presence of the Galactic corona. First predicted by Spitzer (1956) to explain the long-term survival of H i high-velocity clouds in the halo of the Galaxy, the presence of a diffuse (${\rm n}\sim 10^{-5}-10^{-2}\leavevmode\nobreak\ {\>\!{\rm cm}^{-3}}$), hot ($\sim 10^{6}$ K), gaseous component encompassing the Galaxy is now well established (q.v. Bregman et al. 2018; Sormani et al., 2018). The impact of a diffuse background on the evolution of infalling gas clouds has been demonstrated repeatedly by a number of independent studies (see below). Therefore, ram pressure models are no longer a matter of choice; they are a physical necessity.

We have shown here that a new ram pressure model based on a realistic description of the Galactic corona, coupled to the first-infall scenario, does not allow for the formation of a leading gas stream like the Leading Arm.

4.2 Potential limitations

The results discussed above are robust to variations in the kinematic properties of the corona, or to the absence / presence of an ambient magnetic field. However, one may wonder whether our main result – the absence of a leading gas stream in models that feature a Galactic corona – is affected by a lack of adequate numerical resolution.

It has long been known that the hydrodynamic interactions governing the evolution of gas clouds moving through a diffuse medium depend on resolution. Recently, Hummels et al. (2018) have shown that the presence and distribution of cool gas in the circumgalactic medium of galaxies can be dramatically affected by resolution. At a more fundamental level, a progressively increasing resolution has been shown to reduce the efficiency of cloud acceleration (i.e. drag; Schneider & Robertson, 2017). However, these authors also point out that a cloud’s survival is appreciably reduced when the resolution is increased. In other words, at a higher resolution the drag experienced by a cloud may be less efficient, but at the same time the cloud is destroyed more rapidly.

A number of studies based on wind-cloud simulations show that infalling clouds with HVC-like densities typically break apart within a few 100 Myr, regardless of the presence of magnetic fields (e.g. Gregori et al., 2000; Heitsch & Putman, 2009; Kwak et al., 2011; Banda-Barragán et al., 2016). Hence, while the typical survival timescale is much shorter than the infall timescale of the MCs onto the Galaxy in our models ($\sim 1$ Gyr), it is not obvious that our adopted resolution is adequate to properly account for the hydrodynamic interaction between the infalling gas and the hot corona.

In order to address this potential shortcoming, we have run two additional simulations, which are in essence identical to our ‘Slow corona’ model, but each at a progressively increasing resolution. The details and result of this experiment are discussed in Appendix B. In brief, we find no appreciable difference between these runs; i.e. in none of them does a leading gas stream survive. This experiment strongly suggests that, while gas drag is less efficient at higher resolution, the leading gas clouds are destroyed early on regardless.

Based on these tests, we conclude that our results – in particular the destruction of the leading gas stream – have not been affected in any significant way by our adopted standard limiting resolution. Our belief is strongly supported by similar results found by previous, albeit less sophisticated ram pressure models of the Magellanic System (see above).

Similarly, we do not believe that our main result is overly sensitive to the physical processes – notably star formation and its associated feedback – we have ignored in our models. We do so even in spite of the early suggestion that the Leading Arm may be a feature associated with stellar feedback within the LMC (the ‘blowout hypothesis’; Nidever et al., 2008, see also Venzmer et al. 2012). According to this scenario, the Leading Arm was created out of gas outflows from a strong star-forming region in the LMC. However, recent numerical work (Salem et al., 2015; Bustard et al., 2018) shows that outflows coupled to the effect of ram pressure may indeed aid to the formation of a trailing gas stream, but they do not contribute in any significant way to the formation of a leading gas stream. We stress that if we insisted in including stellar feedback in our models, a fully self-consistent treatment would require us to take into account the feedback from the Galaxy (i.e. Galactic winds), which may only further contribute to the destruction of the infalling gas (e.g. Cooper et al., 2009).

So, taking the results presented here at face value, how can we understand the existence of the gas referred to as ‘Leading Arm’? In our view, there are only two plausible interpretations, both with far-reaching consequences and issues, which we discuss in the following.

4.3 An alternative Galactic corona model

Although the presence of diffuse, hot gas around the Galaxy seems secure, its structure is still fiercely debated. Two competing models are generally invoked to explain the local x-ray observations: i) a disc-like structure with a scalelength and a scaleheight both of a few kpc (e.g. Yao et al., 2009); or ii) an extended ($\gtrsim 100$ kpc), spheroidal gaseous halo (e.g. Miller & Bregman, 2013). The models of the evolution of the Magellanic System we have presented here so far are based on the latter, more conventional interpretation. A recent analysis of Suzaku archival x-ray data suggests, however, that neither of these interpretations may be appropriate (Nakashima et al., 2018).

The latter study shows that the x-ray coronal emission can be explained with a hybrid of the two competing models introduced above, i.e. a model that consists of both a dense, disc-like component that dominates the x-ray emission; and a diffuse, spheroidal, extended component that dominates the mass of the corona, with a mean temperature of $\sim 3\times 10^{6}$ K. As noted by Nakashima et al. (2018), the x-ray data is highly anisotropic, with large scatter across the sky, making simple spherically symmetric models rather unlikely. Given this result, it is interesting to explore whether a model of the evolution of the Magellanic System within a Galaxy that includes a ‘corona’ based on this hybrid interpretation could allow the formation and survival of the Leading Arm .

To this end, we construct a model of a hot gas distribution around the Galaxy following Nakashima et al. (2018)’s parametrisation, consisting of an axisymmetric disc described by

[TABLE]

with an extension of 30 kpc in both directions, and total mass of $\sim 5\times 10^{7}$ M*⊙*; and a spherically symmetric halo described by

[TABLE]

with a mass of $\sim 2\times 10^{9}$ M*⊙* enclosed within 250 kpc. The corresponding parameter values are $n_{\rm d,0}\approx 4\times 10^{-3}$ cm*-3*, $R_{\rm d}=7.0$ kpc, $z_{\rm d}=2.0$ kpc, and $n_{\rm d,0}\approx 10^{-3}$ cm*-3*, $r_{\rm c}=2.4$ kpc, $\beta=0.51$.171717We neglect the metallicity dependence introduced by Nakashima et al. (2018) in their calculations of each component’s mass. The composite radial density profile corresponding to this model is compared to our ‘Non-rotating corona’ model in Figure 1. We refer to this as the ‘Hot disc-halo’ model. It is worth emphasising that at $R\gtrsim 10$ kpc the density profile resulting from this model roughly corresponds to the density profile of our ‘Non-rotating corona’ model but scaled by a factor $\sim 1/10$.

Using the ‘Hot disc-halo’ model, we calculate the evolution of the Magellanic System as we did before with our conventional corona models. The result of this exercise is shown in Figure 11. Clearly, a leading gas stream survives in this model up to the present-day. Note that the leading stream is, however, less extended compared to the leading gas stream in models that lack a hot corona. These results are not fully unexpected in view of the role ram pressure plays in the survival of a leading gas stream discussed earlier and the lower density of the hot gas in this model compared to our standard corona models.

Our main result is in line with Murali (2000) who found using pure analytic methods that an ambient density $<10^{-5}$ cm*-3* at 50 kpc, i.e. an order of magnitude lower than our standard coronal model, is required for the leading gas stream of the Magellanic System to survive. The density profile of our ‘Hot disc-halo’ model is consistent with this upper limit (Figure 1). Thus we conclude that any coronal model with a mean density lower by a factor of $\gtrsim 10$ compared to our standard coronal model allows the leading gas stream to survive.

4.3.1 Spheroidal vs. flattened coronae

The survival of the leading gas stream in our ‘Hot disc-halo’ model is an argument in favour of flattened over spheroidal coronal models. Here we present additional evidence in support of flattened hot halo models.

Diffuse, hot coronae around massive ($\gtrsim 10^{12}$ M*⊙*) spiral galaxies extending out the virial radius of their host halo in near hydrostatic equilibrium are expected in theories of galaxy formation (Fall & Efstathiou, 1980). X-ray images of quiescent spiral galaxies (e.g. Bogdán et al., 2013b; Bogdán et al., 2013a) seem to support this picture, although no x-ray haloes extending beyond $\sim 70$ kpc have been observed. Nonetheless, extended, spheroidal coronae are still seriously considered in the context of the so-called ‘missing baryons’ problem (Fukugita et al., 1998). Indeed, some claim (q.v. Faerman et al., 2017) that the Galactic corona may be massive enough to bring the baryon fraction of the Galaxy close to the expected cosmic value.

Simulations of galaxy formation within a cosmological framework (Crain et al., 2010; Nuza et al., 2014; Nelson et al., 2016) self-consistently lead to the formation of extended (i.e. beyond the optical disc) structures of hot gas around spiral galaxies; however their shapes are rarely spheroidal, as some models of the Galactic corona assume (e.g. Miller & Bregman, 2015). If the Galactic corona is spinning fast as inferred from x-ray observations of the Galaxy (Hodges-Kluck et al., 2016), then it cannot be spherical, i.e. an initially spheroidal spinning corona will eventually flatten as it spins up (Pezzulli et al., 2017; Sormani et al., 2018). Conversely, a non-rotating corona will eventually acquire enough angular momentum as a result of the momentum transfer from the gas disc ejected in the form of a Galactic fountain (Fraternali & Binney, 2008; Marinacci et al., 2011). In other words, both a significant rotation speed and a spherical shape – both of which are inferred from the same x-ray data, it should be noted–, are mutually exclusive properties of the Galactic corona.

Therefore, we believe that the hot gas around the Galaxy is likely to be a moderately flattened, spinning structure although the degree of flattening is an open question. An equilibrium, baroclinic model of the Galactic corona that is consistent with a number of observational constraints (but that ignores x-ray observations) also shares these properties (Sormani et al., 2018). We note with interest that some evidence for a flattened hot gaseous structure around the system NGC 891 has been recently presented (Hodges-Kluck et al., 2018). This galaxy is widely considered to be the best edge-on surrogate or analogue to the Milky Way.

The specific model of this kind considered here implies a total mass of hot gas ($<10^{10}$) M*⊙* that is low compared to conventional models, and clearly cannot account for the missing baryons in the Galaxy, which are expected to have a total mass $\sim 10^{11}{\>\!{\rm M}_{\odot}}$. This is a consequence of the density profile of this model, which is low compared to some of the available estimates of the gas density at $r\gtrsim 50$ kpc all around the Galaxy (Figure 1). However, it is consistent with the upper limits inferred by Stanimirović et al. (2002), and Murali (2000), as mentioned above.

In summary, taking both the observational evidence and theoretical arguments into account, the distribution of hot gas around the Galaxy likely follows a flattened distribution as a result of its relatively high rotation speed. Such a distribution may permit formation of the Leading Arm within the framework of a ram-pressure model coupled to the first-infall scenario. The caveat is that such a distribution does not leave room to accommodate the large reservoir of baryons required to raise the total baryonic mass fraction of the Galaxy to the cosmic value, as argued by others (Faerman et al., 2017).

4.4 An alternative view on the nature of the Leading Arm

The physical need for the Galactic corona, as well as our current ignorance about its detailed structure, and the difficulty in the formation of a leading gas stream in its presence, leads us to consider an alternative scenario: that the Magellanic Leading Arm is not in fact directly associated to either the LMC or the SMC.

The interpretation of the Leading Arm as the tidal counterpart of the Magellanic Stream is based on two observations: i) The apparent continuity in the kinematic properties of the H i gas from the tip of the Leading Arm (regions LA II and LA IV)181818See Venzmer et al. (2012) for a definition of the different components associated with the Leading Arm. across the Galactic plane to LA I and to the MCs (Putman et al., 1998); and ii) the chemical signature of the H i gas across the Leading Arm region, which is comparable to the chemical signature of the MCs gas (e.g. Sembach et al., 2001).

However, while recent kinematic H i data suggests an LMC origin (Nidever et al., 2008; Venzmer et al., 2012), accurate measurements of the chemical composition of the Leading Arm clouds favour a SMC origin (Fox et al., 2018; Richter et al., 2018). There is, however, an inherent danger in the use of the metallicity as a constraint on the origin of the gas, in that the interaction between two gas phases results in mixing that inevitably washes out the original chemical signature (Henley et al., 2017, see also Tepper-García & Bland-Hawthorn 2018a).

But there is an additional issue specific to the components LA II and LA III north of the Galactic plane ($b>0$∘), henceforth collectively referred to as LA-North. A number of models of gas cloud-gas disc interactions indicate that gas clouds will not survive a collision with the Galactic gas disc (Tenorio-Tagle et al., 1986), regardless of the presence of an ambient magnetic field (Santillán et al., 1999; Galyardt & Shelton, 2016), unless their density is high enough compared to the disc or surrounded by a DM subhalo (Nichols et al., 2014; Tepper-García & Bland-Hawthorn, 2018a). Assuming the Galactic gas disc follows the radially declining H i surface density profile as estimated by Kalberla & Dedes (2008) at the time of the interaction, the total gas column density of the Leading Arm prior to collision had to be $\gtrsim 10^{21}$ cm*-2*, substantially higher than total gas column density of LA I (which has not yet crossed the disc) as inferred from observations (Fox et al., 2018).

We conclude that the evidence so far for the connection of the Leading Arm with the MCs is somewhat inconclusive. This opens the door to, or perhaps even demands, an alternative interpretation about the nature of the Leading Arm.

Mastropietro et al. (2005) showed that the Leading Arm could be explained as a trailing stream that arises from the LMC as it encircles the Galaxy. This interpretation requires the LMC to have orbited the Galaxy once already. Although unlikely, it should be noted that such an orbit is not entirely ruled out by recent measurements of the MCs’ proper motion (Kallivayalil et al., 2013). Nonetheless, it would appear difficult for this model to explain the significant spread in elemental abundances observed among the different Leading Arm components ($\sim 0.1-0.3$ Z*⊙*; Fox et al., 2013; Richter et al., 2018).

A different paradigm has recently been laid out by Hammer et al. (2015), based on an scenario originally proposed by Sawa & Fujimoto (2005), according to which the Magellanic gas structures are associated with a large number of baryon-dominated tidal dwarf systems, whose largest members are the MCs (see also Pawlowski & McGaugh, 2014). In this picture, a group of essentially DM-free, tidal dwarf galaxies lost their gas along their orbit all the way from the system M31 to the Galaxy as a result of ram pressure by the intergalactic medium and the Galactic corona.191919The idea that the MCs underwent a close encounter with M31 before being captured by the Galaxy dates back to Shuter (1992); Byrd et al. (1994). As for Mastropietro et al. (2005)’s model, in Hammer et al. (2015)’s model the Leading Arm is interpreted as a trailing gas stream of the tidal dwarfs moving ahead of the MCs. We consider this idea by appealing to newly discovered satellites (including streams) around the Galaxy as the possible origin of the Leading Arm (e.g. Helmi et al., 2018).

4.4.1 Is the Leading Arm linked to the infall of a satellite moving ahead of the MCs?

The paradigm we are proposing is that the Leading Arm components south of the Galactic plane (i.e. regions LA I, LA III at $b<0$∘, henceforth referred to as LA-South) can still be directly associated to the MCs, but LA-North most likely are not. Our conviction comes from the fact that no pure gas structure could have survived the interaction with the gas disc, as discussed above, unless confined by a deep DM potential, i.e. a dwarf subhalo possibly (but by no means necessarily) associated to the MCs, a system collectively referred to as the Magellanic Group (D’Onghia & Lake, 2008; Nichols et al., 2011).

Evidence for the existence of a Magellanic Group is widespread.202020The idea that other Galactic satellites may be associated to the Magellanic System dates back to Lynden-Bell (1976). Data taken as part of the Dark Energy Survey (DES) suggest that there are $\sim 70$ satellites associated to the LMC (Jethwa et al., 2016). A number of potential candidates of satellites associated to the LMC have also been identified in the Gaia DR2 (Kallivayalil et al., 2018). Simulations of structure formation within a cosmological context predict that LMC-SMC-like dwarf systems are surrounded by their own satellites at all epochs (Wheeler et al., 2015).

Our conjecture is that the LA-North could be a mixture of gas that was stripped from a ‘frontrunner’, a satellite galaxy moving ahead of the MCs, and the interstellar medium (ISM) of the Galaxy at the crossing point. LA-South would then correspond to the gas stripped away from the dwarf by ram-pressure by the corona prior to the frontrunner’s transit of the disc, much alike what has been observed in an external system (Cramer et al., 2018).

We have illustrated these mechanisms in our work on the Smith Cloud (Tepper-García & Bland-Hawthorn, 2018a, see also Galyardt & Shelton 2016). In essence, in that work we showed that a $\sim 10^{9}$ M*⊙* gas-bearing DM subhalo traveling at $\sim 350$ km s*-1* is partially stripped of its gas while interacting with the Galactic corona and mostly stripped when transiting the Galactic gas disc, lifting up gas from the ISM, and leaving behind a gas streamer of mixed composition. Crucially, since the streamer moves around with the Galactic disc, it will in general not be aligned with the orbit of the infalling subhalo, unless it is very young. In other words, the stripped gas will not generally track its parent satellite orbit.

If the Leading Arm components are in fact fragments of a gas stream created by a gas-bearing satellite that has lost its gas via ram-pressure by the Galactic corona and the Galactic gas disc, the only way we can tightly constrain its progenitor is by having a handle on the distance and the age of the Leading Arm. A recently discovered star cluster in the Galactic halo, believed to be associated to the LA II, would place it at a Galactocentric position of $(R,z)\approx(23,15)$ kpc (or an heliocentric distance of $D\approx 29$ kpc; Price-Whelan et al., 2018). The LA II component extends over roughly 30*∘* on the sky. At an heliocentric distance of $\sim 29$ kpc, this corresponds to a physical extension of $\sim 16$ kpc. Based on our work on the Smith Cloud, we estimate an age on the order of 200 Myr for the LA II, and a distance of the progenitor of $\lesssim 100$ kpc. It is worth noting that our age estimate is consistent with the age of the star cluster believed to be associated to the LA II ($\sim 130$ Myr; Price-Whelan et al., 2018).

To sum up, our key constraints on the frontrunner are: i) moving ahead of the MCs on a prograde orbit with respect to theirs; ii) an orbital history that brings it within the optical radius of the Galaxy; iii) a present-day distance consistent with having travelled $\lesssim 100$ kpc along its orbit. It would be desirable to have evidence that the candidate may have lost all of its gas only until recently; and perhaps less important: an orbit roughly consistent with the orbit of the MCs. The reason for the last constraint to be a weak one is two-fold: First, the frontrunner need not be directly associated to the MCs. Secondly, its orbit may have been significantly affected by the interaction with the Galactic gas disc (cf. Tepper-García & Bland-Hawthorn, 2018a, b).

We consider plausible candidates from a subset of the Galactic satellites discovered so far. In the left panel of Figure 12, we show the positions of the satellites in Stream (MS) coordinates (Nidever et al., 2008), where the LMC is at $l_{\rm MS}=0$∘ and $b_{\rm MS}=0$∘, by definition. The pole of this MS coordinate system is at $(l;b)=\left(188.5^{\circ};7.5^{\circ}\right)$, the MS normal $\hat{n}_{\mathrm{MS}}$. On the right of Figure 12, we show the angle between a satellite’s orbital pole $\hat{n}_{\mathrm{sat}}$ and the MS normal, $\theta_{\mathrm{MSnormal}}=\arccos\left(\hat{n}_{\mathrm{sat}}\cdot\hat{n}_{\mathrm{MS}}\right)$, as a function of $l_{\rm MS}$. The data has been taken from Fritz et al. (2018a, b). Note that the gap at around $l_{\rm MS}=30-50$∘ is due to obscuration by the Galactic disc. Systems with $\theta_{\rm MS}<90$∘($\theta_{\rm MS}>90$∘) are moving on prograde (retrograde) orbits with respect to the MCs. The smaller $\theta$, the more consistent the orbit with that of the MCs. The region of interest here is thus identified by $l_{\rm MS}>50$∘and $\theta_{\rm MS}\ll 90$∘.

Clearly, there is already a large number of dwarf galaxies observed to be moving ahead of the Magellanic Clouds, with orbits that are roughly consistent with the orbit of the MCs. And in the impending era of the Large Synoptic Sky Survey (LSST), it is likely that many more are yet to be discovered.

Of all, the most interesting candidates are perhaps Draco I and Ursa Minor. These are bound to the Galaxy (McConnachie, 2012; Kafle et al., 2014). Both Draco I’s (Mashchenko et al., 2006) and Ursa Minor’s (Piatek et al., 2005) proper motion measurements are consistent with orbital histories that bring the dwarfs within the Galactic optical disc, with orbital periods of $\sim 1.5$ Gyr. Their orbits are well aligned with the orbit of the LMC (Figure 12; see also Pawlowski et al., 2015). Both Draco I and Ursa Minor have dynamical masses in excess of $10^{7}$ M*⊙*, large enough to lift up ISM gas during a disc crossing, but it is currently unknown whether they initially had gas that was recently lost (McConnachie, 2012). Interestingly, Caproni et al. (2015) model the gas loss in Ursa Minor and conclude that stellar feedback is not enough, but additional mechanism are required to fully explain its gas loss, perhaps ram pressure stripping by the Galactic gas disc. Nevertheless, based on our current knowledge of their orbits, it is unlikely that any of these transited the Galactic gas disc within our constraint of $\sim 200$ Myr.

The next most interesting candidate is Crater II with an estimated dynamical mass of $\sim 5\times 10^{6}$ M*⊙* (Caldwell et al., 2017), at a distance of $\sim 117$ kpc, an orbit well aligned with that of the MCs, and a pericentric distance within the optical radius of the Galaxy (Fritz et al., 2018b). However, the epoch of its last disc crossing is likely too long ago to be consistent with our time constraint.

Thus, neither of these systems represents really a satisfactory candidate for our sought frontrunner. In particular, none of them shows evidence for young stars. A suitable candidate may have yet to be unveiled. The prospect of doing so, however, seems low, given the expectation of a young age, characterised by bright A-type stars, making it unlikely that such a system has escaped detection so far.

5 Summary

Our new study has highlighted the difficulty of understanding the Leading Arm as a gas stream that is tidally stripped and runs ahead of the Magellanic Clouds. The latter has become the conventional interpretation since the seminal H i study of Putman et al. (1998). Simple dynamical models support this interpretation assuming the Clouds are on first infall. But we have shown that realistic models for the Galaxy that include a pervasive, hot coronal halo challenge this basic assertion because gaseous ‘leading arms’ do not form.

In brief:

1. It may be possible to rescue the conventional ‘Leading Arm’ interpretation if the hot halo is truncated well within the distance of the Clouds, e.g. through flattening by rapid rotation. The appeal of this scenario is that it is consistent with the recent claims based on x-ray data of the Galaxy, and with the recent observation of a flattened hot halo in an external disc galaxy. The disadvantage is that a truncated hot halo cannot accommodate the ‘missing baryons’ claimed by many authors, both observers and simulators, to be locked up in the hot halo of the Galaxy.

2. A belief shared by some is that a population of H i clouds in the Northern Galactic hemisphere appears to be a continuation of the Leading Arm and may be associated with it. We consider this connection to be highly unlikely because the Leading Arm would then need to survive its transit of the Galaxy’s gaseous disc. It is difficult to conceive how the latter can be true given the results of numerous studies that strongly argue against it.

3. We consider an alternative explanation involving a gas-bearing dwarf galaxy running ahead of the Clouds (‘frontrunner’) that is being ram-pressure stripped of its gas. The advantage of this scenario is that it requires the existence of an extended hot halo in line with the conventional picture, and may simultaneously explain: i) the enhanced metallicity of the Leading Arm gas (relative to the trailing Stream); ii) recent claims of star formation in the Leading Arm as a result of ram-pressure acting upon the gas (Cramer et al., 2018; Alig et al., 2018); and iii) the survival of the gas as it transits the disc as a result of the shielding effect of the host potential. The disadvantage is that a suitable candidate dwarf has yet to be identified.

Overall, it seems difficult to escape the conclusion that there is something fundamental about the Magellanic System and the gaseous structure of the Galaxy that we do not yet understand.

Acknowledgments

We thank S. Pardy for providing the geometric transformations used to setup the initial positions and velocities of the MCs for our infall models. We thank the anonymous referee for carefully reading our manuscript and providing insightful comments that helped improving the presentation of our results. TTG acknowledges financial support from the Australian Research Council (ARC) through an Australian Laureate Fellowship awarded to JBH. MSP acknowledges funding by NASA through Hubble Fellowship grant #HST-HF2-51379.001-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS5-26555. We acknowledge the facilities, and the scientific and technical assistance of the Sydney Informatics Hub (SIH) at the University of Sydney and, in particular, access to the high-performance computing facility Artemis and additional resources on the National Computational Infrastructure (NCI), which is supported by the Australian Government, through the University of Sydney’s Grand Challenge Program the Astrophysics Grand Challenge: From Large to Small (CIs: Geraint F. Lewis and JBH). All figures created with gnuplot,212121 http://www.gnuplot.info originally written by Thomas Williams and Colin Kelley, and matplotlib (Hunter, 2007). This research has made use of NASA’s Astrophysics Data System (ADS) Bibliographic Services222222 http://adsabs.harvard.edu , and of Astropy232323 http://www.astropy.org , a community-developed core Python package for Astronomy (Astropy Collaboration 2013).

Appendix A Model verification

Here we present our results intended to replicate some of Pardy et al. (2018)’s results corresponding to their 9:1 model. This step is relevant to make sure that we are using the correct initial conditions for all our models. We stress that the agreement between their model and our replica is non-trivial, given the significant difference in codes (i.e. numerical techniques, and implementation) used by them and us to create and evolve the initial conditions.

A glance at our Figure 14 and comparison to Pardy et al. (2018, their figure 2) readily reveals the similarity in the evolution of the MCs in isolation in their model and our realisation. We note that the 9:1 model by Pardy et al. (2018) consists of extremely gas-rich dwarf galaxies with gas discs as extended as the Galactic gas disc. We have chosen not to use gas disc as extended as theirs, but we do adopt the same mass.

We find the best match to the final configuration for the system evolved in isolation in terms of number of passages of the SMC around the LMC and its relative distance presented by Pardy et al. (2018) after 5.76 Gyr, a perfect match to their quoted time scale of 5.8 Gyr.

In order to model the infall of the MCs, the evolved binary system is mapped onto the Galactic reference frame. In brief, the MCs are placed around the Galaxy’s centre by applying a series of geometric transformations to the position and velocities of the MCs’ constituents. The following are the exact steps we took to set the initial Galactocentric positions and the Galactocentric velocities of the MCs in all of our models, after they have been evolved in isolation, kindly provided by S. Pardy (private communication):242424Note that the rotation as defined here follows the mathematical convention that a positive angle corresponds to a counter-clockwise rotation.

Calculate the centre of mass (CoM) and its velocity (CoV) of the LMC. We get $\vec{r}_{\textrm{\sc lmc}}\approx(-70,3.5,0.5)$ kpc and $\vec{v}_{\textrm{\sc lmc}}\approx(-19,-17,-0.2)$ km s*-1*, respectively; 2. 2.

Shift both the LMC and SMC to a coordinate frame where the CoM of the LMC is at $\vec{r}=(0,0,0)$ and its CoV at $\vec{v}=(0,0,0)$; 3. 3.

Apply a rotation of:

66.00 deg around the $x$ axis; then

• -

231.2 deg around the $z$ axis; 4. 4.

Apply a constant shift to the SMC and the LMC position and velocities of $\vec{r}_{\rm shift}=(48.0,198,-85.0)$ kpc, and $\vec{v}_{\rm shift}=(-17.0,-160,-29.0)$ km s*-1*, respectively; 5. 5.

Apply a rotation of:

31.9 deg around the $x$ axis; then

• -

18.5 deg around the $z$ axis.

The position and the velocity of the SMC we obtain after these transformations is $\vec{r}\approx(-31.4,217,20.8)$ kpc and $\vec{v}\approx(41.3,-274,-73.3)$ km s*-1*, respectively, to be compared to $\vec{r}\approx(-29.9,217,22.0)$ kpc and $\vec{v}\approx(42.9,-284,-80.0)$ km s*-1*from Pardy et al. (2018, see their table 3). Note that the position and the velocity of the LMC we get are identical to theirs, by construction.

The agreement between Pardy et al. (2018)’s and our results is reinforced by a comparison of the orbital history presented in Figure 13 and Pardy et al. (2018, their figure 1).252525Compare also to Besla et al. (2012)’s model 2 (their figures 2b, 4b and 4d). This demonstrates that not only our initial orbital parameters, but also the Galactic potential as well as the masses and structure of the MCs agree well with theirs.

Appendix B Convergence tests

We considered in in Sec. 4.2 the possibility that our results have been affected by a lack of adequate resolution. This is motivated by recent numerical work which demonstrates that a progressively increasing resolution leads to a less efficient drag, but at the same time decreases a cloud’s lifetime. We have run two additional simulations, which are in essence identical to our ‘Slow corona’ model, but at a progressively higher resolution. The model chosen is irrelevant, as any other model behaves roughly in the same way. In one of these simulations we adopt a refinement strategy whereby a cell is refined if the gas mass exceeds $\sim 5\times 10^{5}\leavevmode\nobreak\ {\>\!{\rm M}_{\odot}}$ (rather than $\sim 2.5\times 10^{6}\leavevmode\nobreak\ {\>\!{\rm M}_{\odot}}$ as we have adopted throughout). In the other simulation we adopt a gas mass threshold of $\sim 2.5\times 10^{5}\leavevmode\nobreak\ {\>\!{\rm M}_{\odot}}$. These correspond to an improvement in the gas mass resolution by a factor 5 and 10, respectively, with respect to our initial adopted setting. The spatial resolution is progressively increased by an additional factor of 2 (per dimension), thus yielding a minimum cell size of 30 pc and 15 pc, respectively.

To assess the effect of the increase in resolution on our results, we compare the final distribution of gas in each of these simulations, paying especial attention to whether a leading gas stream is able to survive in any of them. The result is shown in Fig. 16. Clearly, the gas distribution is not appreciably different among these three simulations. In particular, none of the simulations features a leading gas stream, despite an order of magnitude increase in the refinement density threshold, and an increase in spatial resolution of 4$\times$ per linear dimension.