Rockstar Comparison

I've downloaded and compiled Peter Behroozi's code Rockstar, which is a 6-D phase space subhalo finder.  This code was one of the more robust and well behaved codes in the subhalo comparison paper that I read, and furthermore, was incredibly well behaved from a user's standpoint as well: compiled immediately, ran (on one processor) very quickly, and is decently documented.  Incidentally, this code provides some parent/child relationships when run on multiple time steps, and also will output the appropriately formatted data for running the public merger tree code known as Consistent Trees
( http://code.google.com/p/consistent-trees ).  This seems to be another reason to consider using it.  Mike has some concerns about the parallelization scheme, which may not scale well he worries.  However, Rockstar was run on the Aquarius Level 1 data which has 4.21 billion particles, so it is not completely devoid of scalability (though 4 is obviously different than the 68 that we are working with ).

Here are some preliminary results.  I've looked at one halo (halo 0000004321 in the lcdm_po simulation) of mass ~9e13.  Here are all the subhalos found by vfind/pos (from the .cen files, so I assume these are NOT leaf halos, but indeed the ones after the olapt step), and also the Rockstar subhalos.

Vfind/pos is in red, Rockstar in blue, and the blue particles were plotted second.  A few things are clear here:  Vfind/pos is identifying many more halos as individual halos (1118) than Rockstar (306), especially near the center, despite the fact that Rockstar keeps subhalos >20 particles, and Vfind/Pos subhalos with >30 (I could fix this later, it was a mistake).  I am guessing that toward the center, the rockstar halos are identified has having a much bigger "Mass" than the Vfind halos.  Possibly there is an adjustable parameter in the olapt step that could bring these into better alignment.

Vfind/pos seems to find nearly all of the Rockstar halos, but the reverse is not quite the case.  I decided that Vfind is just plotting way more subhalos, so I first decided to rank them by central density (this is the ONLY quantity I have available in the .cen files, not Npart and not mass either), and plot the first N of them, so that I am plotting the same number of objects for the two finders.

This doesn't look very good, and we can understand why, because the halos toward the center are definitely going to have higher central densities on average, just as a selection effect... we can't find low central density subhalos contrasted against the parent toward the center of the parent.  So in the ranking, I am preferentially selecting these central guys, and you can see them all bunched up there at the center.  So... I wanted to see if I could make this look a little better.  I ranked the Rockstar subhalos by number of particles, and started plotting only the most massive of them, and also plotting some arbitrary number of Vfind/pos halos, tuned to get the best matchup in the outer halos as I could manage.  Here's a plot, but keep in mind this has been fiddled with to give visually not so terrible results.

This is starting to look somewhat more reasonable.  I have plotted the 50 most massive Rockstar halos, and the 75 halos with highest central density from Vfind/pos.  Aside from the weirdness at the center, there are about 35ish outer halos for which there is overlap, and around 10 in each color that the other finder did not identify (above the thresholds I set, they probably did find them overall).

So this tells us that yes, these things are doing approximately the same calculation, but there are non-trivial differences.   One possible next step would be to try and match up subhalo "masses" from the olapt files (assuming the halo ids match those in the .cen files, I am not clear on this point).  This would be smart because it's likely that central density and "mass" have quite a bit of scatter in their relationship, so these plots compare apples to oranges.  However, there is a good bit of evidence in the literature that "mass" is not a good variable for subhalos.  It is not stable or reliable, and does not match up well on an inter-algorithm basis, even in terms of rank ordering.  There is evidence that vmax (max circular rotational velocity) is a much better variable to use.  We do not currently compute this (that I know of) but we do compute a halo center, so it would not be so terrible to add such a calculation it to hproc or hpick.   The magnitude of vmax will depend on 1) the center finding algorithm and 2) whether the mass in the background halo is included or not.  However, it appears that the rank ordering for vmax is pretty robust according to the literature, so using it for both finders we can do a better side-by-side comparison.  Rockstar already outputs this quantity.

Since rockstar appears to run on on this halo with only one processor, it will probably manage to be able (in parallel) to cope with even bigger halos in finite time, and we could consider simply moving forward with the project using Rockstar instead of Vfind/pos, rather than spending a lot of time working out what the origin of the differences may be.  This depends on whether we believe there is something fundamentally interesting about the differences and/or we thing Vfind/pos is more well physically motivated than Rockstar.

Feel free to please please give me some feedback, especially on stuff in Red.

Plots from the talk

I think for continuity, it makes sense to throw the plots I used in my talk at SF12 onto this blog, so that everything is all in one place.  All plots were made with with the PO simulation, which is 4096^3 particles in a 500 Mpc box.   Here's a little table I made describing this suite of sims:

These results are preliminary, in the sense that I did not have a complete sample of halos to work with, but above 10^14, the sample is complete.  Incidentally this required much brute-forcing in the form of restarting hte analysis over 25 times, in order to get enough halos in place to put a talk together.  I need help fixing the pipeline.

First, the relation between the velocity dispersion of the dark matter particles in the po simulation, and the FoF halo masses:

Red points are halos for which subhalo data exists, blue are the ones I haven't got done yet (they overlay the red and saturate).  One very clear thing here, the friends of friends algorithm is linking together some objects that are actually two separate halos, and thereby causing a large scatter upward.  As you can also see, this is artifacty.  Mike thinks this is because I didn't let halopick know that the velocities were peculiar velocities (i.e. do not include the hubble flow*) so it was subtracting off a large number before computing the velocity dispersion.  In principle of course this should not matter (mathematically) but probably does numerically.  Still, it's a funny manifestation and I am skeptical of the diagnosis.  This needs to be re-run and checked.

Next, I'd like to show the mass function and dark matter velocity dispersion function for these halos:

It's impressive to have so many tens of 1e15 solar mass halos in the mass function.  I would like to check the mass function against the Battycharia fit, just as a sanity check on these new simulations.
Again here you can see the impact of the relatively long linking length of b=0.2 in the velocity dispersion, where there are several tens of objects at the high end, whereas in the mass function, there are only a few in top several bins.  Also, I think the relatively high scatter in the mass-vd relation is causing the vd function to roll more slowly, rather than being a fairly obvious powerlaw with a roll off, as in the mass function. This might make it slightly less useful as a cosmological probe, but maybe not.  This particular diagnostic is not an observable anyway.

Next I'd like to show my preliminary work on connecting the subhalo population with various HOD populating techniques.  If we identify the subhalos with galaxies, this next plot shows the HOD, that is, the mean occupation of halos as a function of their mass:

Note that the subhalos used to make this plot are strictly the "leaf" subhalos, that is, subhalos with no further substructure of their own.  The top plot just gives you an idea of how many halos are in each bin (e.g. you should probably ignore the top three bins that have fewer than 10 halos in them).  The four colored power laws in the bottom panel are for different subhalo mass thresholds.  We see that these indeed look like a power law.  One thing on the to-do list might be to fit a power law slope to these simulation.  You should be skeptical of the 10^7 halos because the smallest "subhalos" in this sample only have a few particles.  The 10^10 is robust, but as you can see, we start running out out of substructure at this mass.  This is mainly because at the 10^11 level, most subhalos will have their own substructure, and therefore do not appear in the "leaf" set.  Another thing on the to-do list:  Get occupation from the olapt files instead of the find files.  Olapt halos are bound subhalos that can contain the leaf halos, so for example, we expect many more of the larger objects than in the leaf sample.  I don't know if it's interesting, but it might be cute to see if the slope of the HOD depends strongly on whether we use leaves only, or olapt instead.  We should think about what it is we learn about when we compare these.  Also, I suspect that this will strongly depend on how well resolved the simulation is.  

Next I'd like to talk about the statistics in a particular bin of the plot above.  You'll note that I've placed error bars on the plot, these represent the 25th and 75th percentile of the ranked occupations in the bin.  Of course the masses do not evenly sample the bin because of the steeply sloping mass function, so the distribution of occupations is fairly skew.  This next plot shows the distribution for a few narrower non-overlapping bins in halo mass:

This is for a particular subhalo mass threshold, and the host halo mass bins used to make this plot are 1-2e14, 3-4e14, 5-6e14, 7-8e14, 9-10e14.   This plot shows that all the mass bins have roughly the same width in occupation.  The assumption that is commonly used in the generation of mock catalogs is that the distribution is poisson about the mean.  It would be nice to explicitly check this.  It could be, for example, that the validity of this assumption could depend on the subhalo threshold that is used, since we know that tidal disruption timescales depend both on the infall times and the infall masses.  The next two plots show more explicitly the dependance of this distribution on both host mass, and subhalo cutoff.

These red and green plots show occupation histograms for two different host halo masses (noted in the plot titles for each panel).  The four panels correspond to four different subhalo mass thresholds.  I see no convincing evidence that the occupation distribution depends in a significant way on the subhalo threshold, but this needs to be checked quantitatively.  To do this, I need to figure out how to convolve the the poisson distribution with the underlying mass function to get the predicted distribution of subhalo occupations for a given mass bin, if they are really sampling from Poisson.  This is another item on the to-do list.  It will be interesting to see if these results change when using olapt halos instead of leaf halos.  Incidentally, you might be fooled into thinking that the fourth panel in the red shows some qualitative difference in the distribution, but I think this is a side effect fo the fact that the lower half of the mass bin may be too small to host 10^10 solar mass leaves.  So I don't think that change of shape is significant.



*this has changed in the file format since Feb of this year

Various things are broken

So, things are working ...sortof.  And sortof not.  Item 1: I ran hpick and got a bunch of halos, but either hpick had a bug in it, or I introduced a bug when I attempted to print out only halos with >1000 particles, or the parent.find file is corrupt and missing halos.  Here is the evidence:

If you recall my plot from before, the pn simulation did not have gaps in the hids, so this seems to be wrong.  Unless halopick kept hn instead of new_hn, that is also possible, but I am pretty sure not.

Nonetheless, because I was in the mood to pretend to do science for a change I made a mass function plot, which will look much more sane after this issue gets sorted out.   Here it is:

And, for kicks and giggles, a velocity dispersion function too:

 

These will all improve once I figure out where the missing halos are.

ID catalogs are done

Just a quick post to say, I now have lists of all halos, their masses, and the line number that they start on.

Here is the complete _pe_ simulation, thinnned by a factor of 100:

Sorry about that green corner there.  And here is the complete _pn_ simulation, also thinned by 100.

I talked to mike about the weird hid gaps, he says that is okay, but given that pn doesn't have them, I am not sure whether to believe him or not.  Anyway.

Testing Latex

$\int_0^\infty x^2 dx = $


This is made possible by this person's blog entry.

Preliminary halo mass plot

For several weeks now, I have been working on generating a mass limited sample of halos.  None of the .find files are sorted by halo mass, the option are halo id, or particle id.   So I have modified haloproc.c to generate a list of halos as it processes them for velocity dispersions and other statistics.  This list contains the halo ID, the number of particles in that halo, and the startline in the .find file.  My plan is to sort this list, and threshold it, and generate individual halo files for all halos above the threshold.  I have already done the individual halo files for just the first 500 halos (and run a subhalo finder on each of them), but I think the subhalo velocity function statistics will not be valid unless we are using ALL halos in the sim above a given mass.  Here is a preliminary plot from the lcdm_pe simulation:

It shows a couple of things.  1) simply thresholding on halo id is terrible, because neighboring id numbers can vary by something like an order of magnitude in mass.  2) there seem to be certain halo id numbers that are just missing.  


Mike, you should tell me if this is what you were expecting to see. 

Incidentally, looking at this plot equivalent for one of the smaller simulations showed that there was a second population of very high mass halos, with very large halo id numbers, so that this relation wasn't even approximately monatonic.  Unfortunately I cannot find this plot, I wish I'd thrown it on the blog!! The same MIGHT be true in this simulation too, you'll note that there are only about 4 million halos plotted here, whereas this sim has something closer to 70 million.

I am unsure what caused haloproc.c to stop after 18 billion particles (possibly I simply didn't ask for quite enough time in my batch submission), but I am working on trying to restart it (not from the very beginning) so that we can get a completed halo list for this simulation.

When I have something better I will post again.

velocity dispersion function

mass versus subhalo velocity dispersion

The zeroth order question is the following:  If you observe the velocity dispersion of galaxies in a halo, how tightly does that correlate to the halo mass?  To answer this I've taken halos in the lcdm_pe simulation, and run subhalo finders on each of them.  Using the mean velocities of each of the subhalos in a halo, I've computed a velocity dispersion.  I have plotted this against the host halo mass.  Here is a plot, using all identified subhalos.

The top panel is the total velocity dispersion, and the bottom panel is the radial velocity dispersion, the quantity that can potentially be observed.   There is a noticeable correlation, with a lot of scatter.

Next, I thought it would be interesting to place mass thresholds on the subhalos used.  This is interesting because we believe that different morphological types of galaxies live in dark matter halos of different masses, e.g. big old red dead ellipticals can be expected to have a higher characteristic dark halo mass than young star-forming spirals.

Below is a sequence of subhalo thresholding, ranging from 1e10 to 1e12.  Note that the y axis on the 5e11 and 1e12 plots have been significantly expanded.

What we see is that as we go to larger subhalos, the velocity dispersion is a much less accurate tracer of host halo mass.  Assuming the brightest galaxies live in the largest subhalos, this means it's unlikely that choosing the brightest few galaxies to compute the velocity dispersion will yeild a result that correlates at all tightly with the mass.

Originally we were thinking of generating light cylinders and computing the velocity dispersion among all galaxies along the line of sight, but to me, the scatter for the bright guys looks large enough that there may not be a point.

Notes to own self:  Procedure for getting these plots:
1) modified tree19/apps/etc/haloproc.c to dump a separate file for every halo in the .find file of lcdm_pe simulation
2) ran isodensity subhalo finder tree16/vfind64/pos/vfind.x86_64-mpi on each halo from step 1
3) ran original tree19/apps/etc/haloproc.corei7 on the .find files from step 2
4) to get subhalo masses ran ~mswarren/tree16/apps/hmassr/hmassr.x86_64 to generate .halos files
5) code to make these plots lives in projects/velfunc/hproc/shvdvsmass

Halo Gallery: Corrected

The bug in Mike's halo files has been corrected, and these are the updated plots of the highest velocity particles in the halo, at early times. 

Halo Gallery II

The previous halo gallery was annoying because it was impossible to see structure, because the plots wer saturated.  I thinned the particles out by a factor of 50...  Here are the same set of halos, with slightly better contrast.