[yt-users] possible bugs in analytical mass function code and halo profiler code

[Michael Kuhlen]

I think I've identified a bug in the analytical mass function code,
and a possible bug or improvement in the halo profiling code:


1) In yt/analysis_modules/halo_mass_function/halo_mass_function.py
   line 327, shouldn't dn_M_z be multiplied by h^3 instead of h^4?

2) In yt/analysis_modules/halo_profiler/halo_filters.py, I've found
   that the way the virial radius is determined often doesn't work
   very well for me.

   a) I commonly find that more than the first few bins in the profile
      have overDensity <= virial_overdensity, even though there's
      nothing wrong with the halo profile. It's just a transient,
      perhaps caused by too small bin sizes. Just a little bit farther
      from the center the profile becomes reasonable again, with
      overDensity > virial_overdensity and decreasing outwards. In
      these cases, this condition (line 100)

      if (overDensity[1] <= virial_overdensity):
          index = 0

      results in this halo getting an incorrect virial radius and
      getting tossed out as "not virialized".

   b) Related to the previous point, I've found that determining the
      virial radius (the bin where overDensity drops below
      virial_overdensity) from the inside-out, as in (line 105)

      for q in (na.arange(len(overDensity)-2))+2:
            if (overDensity[q] < virial_overdensity):
                index = q - 1
                break

      doesn't work as well as an outside-in approach.

   I've gotten much better results, i.e. fewer rejected halos that
   were perfectly fine and better virial radii, by doing this:

   if (overDensity[-1] >= virial_overdensity):
       index = -2
   else:
       for q in (na.arange(len(overDensity),0,-1)-1):
           if (overDensity[q] < virial_overdensity) and
(overDensity[q-1] > virial_overdensity):
               index = q - 1
               break

3) Lastly, one question about the radial profiles. I've noticed that
   the 'CellVolume' field is rarely equal to 4.0/3.0 * pi *
   ('RadiusMpc' * 3.09e24)**3. Typically 'CellVolume' is about 35 -
   40% larger. Maybe this is because the grid cells partly "overshoot"
   the edge of the sphere? But I'm surprised that this is a 35%
   effect... This potentially affects the determination of the virial
   radius.


Thanks,

Mike


-- 
*********************************************************************
*                                                                   *
*  Dr. Michael Kuhlen              Theoretical Astrophysics Center  *
*  email: mqk at astro.berkeley.edu   UC Berkeley                      *
*  cell phone: (831) 588-1468      601 Campbell Hall                *
*  skype username: mikekuhlen      Berkeley, CA 94720               *
*                                                                   *
*********************************************************************