Solving the puzzle of subhalo spins

Investigating the spin parameter distribution of subhaloes in two high resolution isolated halo simulations, re- cent work by Onions et al. suggested that typical subhalo spins are consistently lower than the spin distribution found for field haloes. To further examine this puzzle, we have analyzed simulations of a cosmological volume with sufficient resolution to resolve a significant subhalo population. We confirm the result of Onions et al. and show that the typical spin of a subhalo decreases with decreasing mass and increasing proximity to the host halo center. We interpret this as the growing influence of tidal stripping in removing the outer layers, and hence the higher angular momentum particles, of the subhaloes as they move within the host potential. Investigating the redshift dependence of this effect, we find that the typical subhalo spin is smaller with decreasing redshift. This indicates a temporal evolution as expected in the tidal stripping scenario.

Due to a lack of resolution in previous generations of large cosmological simulations, subhalo spins have not been thoroughly investigated so far, despite their application within current semi-analytic models (Guo et al. 2011). Initial work by Lee & Lemson (2013) analyzed the spins of the two most massive substructures of Local Group like systems in the Millennium-II simulation (Boylan-Kolchin et al. 2009) and revealed possible consequences for the application of subhalo spins to near-field cosmology. Onions et al. (2013) investigated the spin distribution of subhaloes in two high resolution simulations of a Milky Waylike halo (the Aquarius simulation (Springel et al. 2008) and the GHALO simulation (Stadel et al. 2009)) analyzed by a variety of subhalo finders. They suggested that subhalo spins are significantly offset to lower values than those seen in typical distribution functions fitted to haloes (Bullock et al. 2001;Bett et al. 2007). As this result is independent of the subhalo finder used, it suggests that this is a true physical effect. This could not be investigated further because their resimulation did not contain a large field halo population. Excluding subhaloes, Colín et al. (2004) found the spin parameter distribution of isolated dwarf dark matter haloes to be perfectly consistent with that of larger haloes. This suggests that the consistently lower spin of substructure is not due to the generally smaller mass of subhaloes, but is more likely related to tidal stripping of high angular momentum material. On the other hand, this offset could also be due to differences between the Aquarius and GHALO simulations and those used by Bullock to define the field relation. To answer this question we require a single simulation that simultaneously includes both a significant subhalo and field halo population.
In this work we use purpose built simulations, specifically designed to contain both a field and subhalo population, to investigate the difference in spin distribution functions between subhaloes and haloes. In Section 2, we present these simulations and the corresponding (sub-)halo catalogues. The different theoretical models of dimensionless spin parameters are described in Section 3, while our results are summarized in Section 4. We discuss our work and conclude in Section 5.

SIMULATION DATA
As we require our halo and subhalo masses to span a wide dynamic range (10 8 M 10 15 M ) , we have run four dark matter only comoving cosmological boxes containing 512 3 particles, with linear sizes of 8 h −1 Mpc, 20 h −1 Mpc, 50 h −1 Mpc, 100 h −1 Mpc respectively (hereafter BoxA, BoxB, BoxC, BoxD). The softening lengths are chosen to be 4% of the mean separation between particles. This set of simulations can both sufficiently resolve subhalo spins (at least 300 particles per subhalo; Bett et al. 2007) and have significant statistics for haloes (c.f. Table 1). We also ran two simulations with the same parameters as those used for the BoxA simulation except for the gravitational softening parameter. BoxA S1 has a smaller softening length whereas BoxA S2 has a larger softening length. We also ran a low resolution simulation containing 256 3 particles and the same linear size, 8 h −1 Mpc as BoxA which we designate BoxLo. The mass resolution of BoxA (2.6 × 10 5 h −1 M per particle) is very close to that of the Aquarius-A simulation at level 4 (Springel et al. 2008) which had a particle mass of 2.7 × 10 5 h −1 M in the high resolution region. This is roughly three times better mass resolution than the Millennium-II simulation (m p = 6.9 × 10 6 h −1 M ) used by Lee & Lemson (2013). The cosmology was chosen to be the same as in the Aquarius simulation, i.e. ΛCDM with Ω M = 0.25, Ω Λ = 0.75, σ 8 = 0.9, n s = 1, h = 0.73. Initial conditions were generated at z = 127 by the code N-GenIC using the Zel'dovich approximation (written by Volker Springel) to linearly evolve positions from an initially glasslike state. This was then evolved to the present day using   (Springel 2005).
All simulations except BoxLo were analyzed with the (sub-)halo finding code SUBFIND (Springel et al. 2001), AHF (Gill et al. 2004;Knollmann & Knebe 2009) and HBT (Han et al. 2012). BoxLo was only analyzed with SUBFIND . A summary of our simulations is given in Table 1 and details of all the halo finding algorithms we have used and a discussion of their relative merits can be found in Knebe et al. (2013).

THEORY
The dimensionless spin parameter indicates how much a collection of particles is supported by the angular momentum against gravitational collapse assuming gravitational equilibrium, where a negligible spin parameter represents minimal support, while the value for a completely supported system depends on the chosen parametrisation. There are two standard parametrisations defined by Peebles (1969) and Bullock et al. (2001), respectively, that we describe in the following two sections. Hetznecker & Burkert (2006) showed that Bullock's parametrisation is less dependent on redshift evolution than Peebles' parametrisation. This is due to it being more robust to variations in the position of the structure's outer radius and therefore not as strongly affected by the many minor mergers over a halo's merging history. Therefore, the two descriptions are not readily interchangeable and results need to be compared using the same parameter. Peebles (1969) proposed to parameterise the (sub-)halo spin in the following way:

Peebles Spin Parameter
where J is total angular momentum, E the energy and M the mass of the (sub-)structure. With this choice a value of λ 0.4 represents a purely rotationally supported object (Frenk & White 2012). Applying this parametrisation Bett et al. (2007) determined the spin distribution of haloes in the Millennium Simulation . The Millennium Simulation has a mass resolution of m p = 8.6 × 10 8 h −1 M and therefore contains very few subhaloes. The vast majority of the objects in the TREEclean catalogue of Bett et al. (2007) are haloes rather than subhaloes. The probability density function of log λ they found to fit the distribution used the following parametrisation; where A is given by, and Γ(α) is the gamma function. They found λ 0 = 0.04326 and α = 2.509 best fit the distribution of halo spins. Bullock et al. (2001) proposed a dimensionless spin parameter of the form:

Bullock Spin Parameter
where J is the angular momentum within a virilized sphere with radius R and mass M , and V is the circular velocity at the virial radius (V 2 = GM/R). They also proposed a parametrisation of the probability density function based on Barnes & Efstathiou (1987), They found the best fit for haloes is given by λ 0 = 0.035 and σ = 0.5.

RESULTS
Note that in this section we will show results for the discrete, normalized derivative of the spin distribution function ∆N (< log λ)/∆ log λ/N tot and ∆N (< λ )/∆λ /N total , while the fitted functions are for the continuous probability density function P (log λ) = dN (< log λ)/d log λ/N total and P (λ ) = dN (< λ )/dλ /N total , respectively. For Peebles spin, we set the bin width to be ∆ log λ = (log λ max − log λ min )/100. For Bullock spin, the bin width is ∆λ = (λ max − λ min )/100.

Halo finding code, softening, and resolution test
In this section, we first test whether simulation data set and the specific choice of substructure finding code will affect the derived spin of (sub-)haloes. Macciò et al. (2008) has tested halo spins with different cosmological parameters. They found that the spin distributions of haloes is essentially independent of cosmology, at least for changes between WMAP1,WMAP3 and WMAP5. We choose not to confirm this result here. BoxA, BoxA S1 and BoxA S2 are used to compare different force resolutions. BoxA and BoxLo are used to compare different mass resolutions. To ensure reliable properties are recovered only (sub-)haloes with more than 300 particles are selected throughout this work. While calculating the spin of haloes, all their substructures are removed from them. Spin distributions are fitted by Equation 2 and Equation 5, and the fitting parameters are listed in Table 2. BoxA * 1132 --215 --* This higher particle number threshold is used to compare with BoxLo over the same halo mass range. Bullock Spin haloes subhaloes λ 0 σ λ 0 σ SUBFIND BoxA S1 0.0308 ± 0.00027 0.655 ± 0.007 0.0123 ± 0.00017 0.827 ± 0.009 BoxA 0.0310 ± 0.00025 0.629 ± 0.007 0.0167 ± 0.00019 0.664 ± 0.009 BoxA S2 0.0320 ± 0.00022 0.615 ± 0.006 0.0224 ± 0.00028 0.629 ± 0.010 HBT BoxA S1 0.0303 ± 0.00026 0.638 ± 0.007 0.0150 ± 0.00013 0.754 ± 0.007 BoxA 0.0309 ± 0.00026 0.626 ± 0.007 0.0198 ± 0.00018 0.674 ± 0.008 BoxA S2 0.0320 ± 0.00025 0.605 ± 0.007 0.0262 ± 0.00031 0.644 ± 0.010 AHF BoxA S1 0.0278 ± 0.00025 0.677 ± 0.008 0.0178 ± 0.00028 0.768 ± 0.013 BoxA 0.0283 ± 0.00019 0.655 ± 0.006 0.0248 ± 0.00046 0.761 ± 0.016 BoxA S2 0.0290 ± 0.00022 0.613 ± 0.006 0.0323 ± 0.00059 0.779 ± 0.015 As Table 2 shows, the recovered spin properties of haloes are largely independent of the choice of gravitational softening. For subhaloes there is a slight trend for the Peebles spin parameter to increase as the softening is increased but this effect is only barely resolved. Such a trend would be expected as a larger softening will produce a shallower core potential, lowering slightly the central kinetic energy and altering the energetics and angular momentum profile, thus affecting the spin parameter. For the Peebles measure (Equation 1), the change in energy is outweighed by the change in angular momentum but the two effects counteract each other. For the Bullock spin parameter (Equation 4), only angular momentum has an affect on the spins. So we find that the Bullock spin of subhaloes is more sensitive to the softening, as shown in Table 2.
Table 2 also shows that, contrary to Onions et al. (2013), the three halo finding methods do not recover consistent spin parameters. While they all agree on the halo spins, AHF recovers significantly larger spins on average for the subhalo population than either SUBFIND or HBT which are consistent with each other. The subhalo spins for AHF are broadly consistent with the field population, particularly for larger values of the gravitational softening. This is also discrepant with Onions et al. (2013) who found lower spin parameter values for their subhaloes. This result is due to the failure of AHF to resolve a significant fraction of subhaloes within the simulations. The subhalo numbers given in Table 1 indicate that around 36% of the subhaloes containing 300 or more particles in BoxA are missed by AHF. Difficulties for AHF in resolving substructures where the density contrast between the subhalo and the main halo is expected to be small have also been reported elsewhere (Avila et al. 2014). Further evidence for this issue is the rising incidence of missing substructures as the box size is increased evidenced in Table 1: for the largest box (BoxD), AHF misses 71% of the subhaloes found by SUBFIND. AHF is missing small subhaloes in the outskirts of the host halo and (as we shall demonstrate later) these small subhaloes are precisely the ones with the lowest spin parameters.
This naturally produces a population of subhaloes with higher average spin parameter for AHF. As HBT and SUBFIND produce consistent results and AHF fails to recover the complete subhalo population we choose to concentrate our analy-sis on SUBFIND for the remainder of this paper.
In Table 3 we test the influence of mass resolution. In the brackets after simulation name, we note the particle number threshold chosen. This is set in order to match (sub)halo masses between BoxA and BoxLo, the lower resolution version of this simulation. This ensures that the halo and subhalo catalogues for BoxLo(N p ≥ 300), BoxA(N p ≥ 2400) and BoxA S2(N p ≥ 2400) have the same mass range. The gravitational softening lengths for both BoxLo and BoxA were set to 4% of mean particle separation. BoxLo and BoxA S2 have the same absolute softening length (1.25h −1 kpc). The results show that mass resolution has almost no effect on spin distribution. The first three rows of each sub-table also give a hint about the influence of softening: as suggested above softening mainly affects the spin of small subhaloes. This is not surprising for the reasons already indicated.

Haloes vs. Subhaloes
So far we have seen that the lower spin parameter distribution observed for subhaloes appears to be a robust result that does not depend upon the choice of halo finder, gravitational softening or mass resolution. Here we explore a possible physical origin for the lower subhalo spins.
The left and right panel of Figure 1 show the Peebles spin and Bullock spin distributions of subhaloes as well as haloes in all the simulation boxes and the respective fitting functions from Bett et al. (2007), Bullock et al. (2001) and Onions et al. (2013). To assure robustness, only structures resolved with at least 300 particles are included. For haloes we remove particles contained within substructures when we calculate the spin. We make fits to the histograms in Figure 1 using Equation 2 and Equation 5. The respective fitting parameters are given in Table 4. Figure 1, Table 2, Table 3 and Table 4 all show that the subhalo spin distribution is different from the halo spin distribution. This further confirms earlier results by Onions et al. (2013) who found that the spin distribution of subhaloes in the Level 4 resolution Aquarius simulation (Springel et al. 2008) is significantly different to the one derived by Bett et al. (2007) for haloes in the Millennium Simulation.
Our results reveal new information about the spin of subhaloes. As the box size grows from BoxA to BoxD, the discrepancy between λ 0 of haloes and subhaloes decreases gradually, i.e. the scale of the effect is mass dependent with larger subhaloes tending to have higher spin. One possibility is that in small simulations such as BoxA or isolated halo models such as Aquarius-A to E studied by Onions et al. (2013), large substructures are generally absent. In the next section we will demonstrate that subhalo spins increase with subhalo mass while halo spins do not have a significant mass dependence.
As an aside, it should be noted that while our fits do not match those given by Bullock et al. (2001), Bett et al. (2007) and Onions et al. (2013) exactly they are within the range of results covered by these works. In practice, previous work does not arrive at an agreement on the exact value of haloes' spin. Most of these studies fix λ 0 in the range 0.031 − 0.045, with σ between 0.48−0.64 (see section 4.1 and Fig.7 in Shaw et al. 2006). In Bett et al. (2007), they found median values of λ med = 0.0367 − 0.0429 for the entire population of haloes, depending on the definition of halo, and λ med = 0.043 for the catalogue of haloes they refined. Different sets of simulations and halo selection criteria may result in this variation of the recovered spin parameter. As Subsection 4.1 demonstrated, such factors as mass resolution and gravitational soft-  Table 1. ening influence the spin. On the other hand, the discrepancy between the spin of haloes and subhaloes within our simulations is much larger than the bias among simulations. It should therefore be regarded as an intrinsic physical property rather than a result of different data sets.
4.3. Mass Dependence As we suggest above, mass dependence can explain the discrepancy between the spins of haloes and subhaloes. To validate this, we further explore the mass dependence of (sub-)halo spin. Figure 2 and Figure 3 show the two dimensional histogram of spin against (sub-)halo mass. They present a  Bett et al. (2007), the parameters have much smaller uncertainty because there are much lager populations of halos in the work (17709121 raw FOF halos including 1332239 "clean" ones). ** In Bullock et al. (2001), the parameters have larger uncertainty because they use less haloes for fitting (only 500 haloes). *** Note that these parameters differ from the ones originally stated in Onions et al. (2013). The original values were derived using an incorrect fitting routine. The values stated here are the correct values fitted to the Aquarius L4 data set. straightforward picture of how the spin distribution changes with mass. We use four simulations to expand the mass range. Contours for each simulation at the same redshift are normalized and stacked together into one subplot. We then divide the sample into 40 bins by log(sub-)halo mass, fit the distribution by Equation 2 and Equation 5 (if the sample volume in that bin is large enough). λ 0 (λ 0 ) in each mass bin is calculated and marked on the plots with a cross, and a linear fit to λ 0 (λ 0 ) against mass is indicated by the red solid line. From the subplots corresponding to redshift 0 (top left and top right), we can see that, for subhaloes, λ 0 (λ 0 ) clearly increases with increasing subhalo mass. In contrast, the spin distribution of haloes is almost independent of mass. The increasing subhalo spin with mass is even more pronounced relationship for the Bullock spin parameter shown in Figure 3 because for the Bullock spin the (sub)halo mass has a higher weight (c.f. Equation 1,Equation 4). This results in a larger discrepancy between haloes and subhaloes at the low mass end. In Onions et al. (2013), their samples are from a Milky Way like re-simulation, which contains subhaloes similar to those found in BoxA. So the subhalo spin distribution in their work is closer to that from BoxA. Onions et al. (2013) suggest that the physical mechanism that drives this difference is mass stripping. Subhalo particles with high angular momentum are stripped preferentially which leads to a decrease in the spin parameter. Subhaloes with low mass are usually the ones stripped most severely. Our results strongly support the claims Onions et al. (2013) made. The slightly positive slope of λ 0,halo (M halo ) and λ 0,halo (M halo ) is inconsistent with some previous work. They found that λ 0,halo is constant or has a slightly negative slope ( (Bett et al. 2007;Macciò et al. 2007, c.f.)). However we should not forget this trend includes the effect of systematic bias between simulations. Figure 10 in Macciò et al. (2008) shows that λ 0 (M halo ) has a slope of 0.005 in a simulation using WMAP1, which is the same cosmology as that used here.

Radial Dependence
To understand further whether tidal stripping of high angular momentum material could cause the lower subhalo spin distribution, we investigate the radial dependence of the subhalo spins. Subhaloes located closer to the center of their host halo are likely to have undergone stronger tidal stripping than those nearer the virial radius. Onions et al. (2013) has done some tests to support their argument, e.g. they analyze the average spin parameter of subhaloes at different distances from the center of the host halo. Here we perform a more detailed test. We stack subhalo samples from four simulations together and then make the two dimensional histogram of spin against their centric distance. Then in each radial bin we fit the subhalo sample using Equation 2 and Equation 5. Finally we make a linear fit as λ 0 (r) = cr/R vir + d or λ 0 (r) = cr/R vir + d. The results are displayed in Figure 4 and show that the spin of subhaloes is suppressed close to the center of the host halo. This is consistent with the argument that subhaloes loose their high angular momentum particles as they are stripped of their outer layers after infall into a main halo.
4.5. Redshift Dependence So far all our analysis was conducted on the z = 0 snapshot. However, the spin distribution is known to change with redshift (Hetznecker & Burkert 2006). Hence, we investigated the redshift dependence of the halo vs. subhalo spin distribution offset. In Figure 2 and Figure 3,we give the Peebles and Bullock spin distributions of (sub-)haloes at z = 0, z = 0.84, z = 2.38 and z = 5.21, respectively. The spin of subhaloes at the lower mass end decreases significantly with time, while the spin of massive haloes increases slightly with time. We calculate the λ 0 (λ 0 ) in each mass bin and then use a linear function logλ 0 = a * log(M/10 10 h −1 M )+b (λ 0 for Bullock spin) to fit λ 0 (λ 0 ) against mass. We list the parameters of each fitting line in Table 5. It is clear the difference between the halo and subhalo spin distribution increases with time. This is consistent with the argument that tidal stripping causes the difference. Affected by stripping, subhaloes loose more and more high angular momentum particles as time passes.
To confirm this result we checked how the spin of a single halo changes with time. We randomly select 6 subhaloes and plot their spin against redshift. We constrain the samples selected so that they are likely to be the subhaloes heavily stripped. As expected, heavily stripped subhaloes should have low mass (here we choose masses less than 10 9 h −1 M ) and have long histories, forming prior to redshift 8. As shown in Figure 5, the spin of subhaloes declines at low redshift. We have checked many more subhaloes not displayed in Figure 5 and we find that most of them display the same trend. We also plot a dashed line for a halo as reference. The halo is in the same mass range as those subhaloes selected. Its spin almost does not change at low redshift. This piece of evidence strongly supports the claim that subhaloes suffer from stripping, loosing their spin over time.
The information at the high redshift end in Figure 5 is not reliable since the progenitors don't contain very many particles. This results in the large fluctuations seen here.

DISCUSSION & CONCLUSIONS
In this work we compared the spin distribution function of haloes and subhaloes in sets of cosmological box simulations. We found that the halo spin distribution function is well fitted by the parametrisations given by Bullock et al. (2001) and Bett et al. (2007) for the Bullock and Peebles (Peebles 1969) spin parameter, respectively. For the subhalo spin distribution function, however, the typical spin of a small subhalo is significantly lower. This was previously suggested by Onions  Peebles (upper) and Bullock (lower) spin distribution of SUBFIND subhaloes against their relative distance from the host halo center. r is the distance of subhaloes from the center of their host halo. R vir is the virial radius of host halo. The crosses indicate the peak of every radial bin and the red line is a linear fitting to these marks. Contours with different colors represent different level of number density ( calculated as dN/dlogλ(orλ )/dlog(r/R vir )/10000 ) as indicated in the color bars to the right. et al. (2013) for the subhaloes within Milky Way-like haloes, but is confirmed here for a full cosmological volume. We investigated the origin of the difference between the halo and subhalo spin distributions. We examined the influence of (sub-)halo finder, spin parametrisation and resolution to confirm that these factors are not the origin of the difference between the halo and subhalo spin distributions. In this process we confirmed the difficulties that the AHF finder has in recovering substructures reliably where the density contrast between the main halo and the subhalo is low. We recom-  mend that the AHF halo finder should be treated with caution in situations where a complete unbiased sample of the subhalo population is required. In this paper, we have argued that the difference between the spin distributions is physical and it is caused by tidal stripping of subhaloes removing high angular momentum material. This argument is strongly supported by three pieces of evidence presented here. Firstly, subhaloes tend to have lower spin when compared to haloes of the same mass. This discrepancy gets larger towards the low mass end. Secondly, the spin distribution of subhaloes is radially dependent within a host halo. Subhaloes closer to the host halo center, which are expected to have been more tidally stripped, have lower spin than those closer to the virial radius. Thirdly the difference between halo and subhalo spin increases with time and hence is being caused by a dynamical effect, such as tidal stripping.
In summary, we have demonstrated that subhaloes typically have lower spin than haloes because tidal stripping removes their highest angular momentum material. This can have an important consequence for galaxy properties that require spin parameter information. Galaxy properties that are related to spin are more likely to be correlated to the spin of the subhalo before infall and not necessarily to its present value.