[yt-dev] [continued] RAMSES units bug in yt-3.0

Wed May 14 07:38:07 PDT 2014

Hi Matt,

Sorry for the confusion - I've cross-checked against another dataset and 
the cell mass should be more like 2.5e8Msun, as you find; I thought the 
density values seemed very low (as well as the mass fraction of gas to 
dark matter) and hence added a factor of boxlen**3, which seemed more 
reasonable, but checking against another run it seems that the other run 
just has very low densities, and so the density unit doesn't need 
converting by boxlen**3 after all, so the only thing to be careful of is 
to make sure that cell-based spatial units are multiplied by boxlen 
while particle-based spatial units aren't; this way both the lengths and 
masses should be consistent. The particle mass should still be 1.8e10Msun.

Thanks for your work on this!

Sam

On 14/05/14 14:22, Matthew Turk wrote:
> Hi all,
>
> I've got Sam's dataset here.  As a quick note, with current yt-3.0
> tip, all results are correct for cosmology.
>
> With the outstanding pull request, here are the results:
>
> Cell Mass: 2.474e+08 Msun (5.300e+13 desired, 4.669e-06 ratio)
> Part Mass: 3.865e+15 Msun (1.800e+10 desired, 2.147e+05 ratio)
>
> Without it:
>
> Cell Mass: 1.146e+03 Msun (5.300e+13 desired, 2.161e-11 ratio)
> Part Mass: 1.789e+10 Msun (1.800e+10 desired, 9.941e-01 ratio)
>
> So the particle mass is correct in current yt-3.0 tip, and the cell
> mass isn't right in either the tip or the pull request.  Now, if I use
> yt-3.0 tip and compare cell_mass*boxlen**6 against the known, I get
> the right answer.
>
> To solve this asymmetry between the particles and the cell density,
> I've added a "code_density" unit and updated the pull request.  I now
> get correct results for non-cosmo and cosmo runs.  This is in
> changeset 8847f04588b1:
>
> https://bitbucket.org/yt_analysis/yt/pull-request/895/fixing-ramses-density-and-mass-units
>
> -Matt
>
> On Tue, May 13, 2014 at 8:59 AM, Sam Geen <samgeen at astro.ox.ac.uk> wrote:
>> I'll try to hunt down and upload a sample dataset shortly.
>>
>>
>> On 13/05/2014 14:52, Matthew Turk wrote:
>>> On Tue, May 13, 2014 at 8:48 AM, Sam Geen <samgeen at astro.ox.ac.uk> wrote:
>>>> Ah, very good, thanks for the link to that old PR! I was worried I was
>>>> going
>>>> crazy... in any case, whichever is the easiest way to make the
>>>> position/length units consistent.
>>> I change the particle positions in the ingestion phase, so at no point
>>> should they appear incorrect.  I'd really like to figure out the
>>> boxlen != 1.0 ...
>>>
>>>> By the way, I took a quick look at reading the grav_ files recently,
>>>> although the output seems kind of strange, so I haven't issued a PR yet
>>>> until I understand what it's doing (plus I still don't 100% understand
>>>> the
>>>> units). I suspect it might be reading in the wrong order or something. It
>>>> also complicates the code, so it might be worth refactoring a bit to make
>>>> a
>>>> universal "AMR-like" data_structures/io interface for the files called
>>>> hydro_, grav_, rt_, etc, rather than what I did which was to scatter a
>>>> bunch
>>>> of conditionals into the existing functions.
>>> That could be -- I think there was some discussion that they might
>>> also include ghost zones?
>>>
>>> I've also run some checks on cosmology runs.  Here's what I've found.
>>> With no changes to the current tip of yt-3.0, I have run this script:
>>>
>>> --
>>> ds = yt.load("output_00010/info_00010.txt")
>>> dd = ds.all_data()
>>>
>>> # Let's find the total mass
>>> pmass = dd["particle_mass"].sum().in_units("Msun")
>>> cmass = dd["cell_mass"].sum().in_units("Msun")
>>> rho = ds.cosmology.critical_density(ds.current_redshift)
>>> crit_mass = (rho * dd["cell_volume"]).sum().in_units("Msun")
>>> print pmass
>>> print cmass
>>> print pmass / (pmass + cmass)
>>> print ((pmass + cmass) / crit_mass)
>>> --
>>>
>>> The results are good.  I get roughly 0.84 for pmass / (pmass + cmass),
>>> and I get 0.998 for the final check, which seems within the parameters
>>> of the simulation.  So the cosmological units should be correct, and
>>> more to the point, they match up with what I'm expecting having read
>>> over the code carefully.
>>>
>>> So it looks to me like there is in fact a boxlen missing, but that as
>>> it stands, the particle and density units are correct for boxlen =
>>> 1.0.  Nick's message is worrisome because the factor isn't exactly 1e3
>>> between his max density.  I have a boxlen != 1.0 dataset, but I have
>>> no reference values to compare it again.  Any chance there's one with
>>> *reference* results for total mass in particles, total mass in gas,
>>> etc?
>>>
>>> -Matt
>>>
>>>> On 13/05/2014 14:21, Matthew Turk wrote:
>>>>> Hi Sam,
>>>>>
>>>>> So, here's an old pull request I found:
>>>>>
>>>>>
>>>>>
>>>>> http://hg.yt-project.org/yt-3.0/pull-request/62/adding-boxlen-to-ramses-units-for-mass-and
>>>>>
>>>>> and then there's the units.f90 from RAMSES itself:
>>>>>
>>>>>
>>>>>
>>>>> https://bitbucket.org/rteyssie/ramses/src/81ab29b8a405e7ccf6bf30fd67582b4edacddd6e/trunk/ramses/amr/units.f90?at=master
>>>>>
>>>>> I'm going to take care of this, update the PR, and ping you both.
>>>>> (And maybe we should move to yt-dev? :)
>>>>>
>>>>> -Matt
>>>>>
>>>>> On Tue, May 13, 2014 at 8:17 AM, Sam Geen <samgeen at astro.ox.ac.uk>
>>>>> wrote:
>>>>>> Huh, that difference in the max density is very strange. I can't think
>>>>>> why
>>>>>> you'd have that difference; it's not an obvious multiple.
>>>>>>
>>>>>> Actually, I'm looking at the code and it's possible that the particle
>>>>>> positions are scaled from 0 to boxlen and the grid positions from 0 to
>>>>>> 1,
>>>>>> though I'd need to open up the raw values in a RAMSES file to confirm
>>>>>> this.
>>>>>> If this is true, I suppose one "fix" could be to hard-multiply the cell
>>>>>> sizes by pf["boxlen"] when you load them so that the length units are
>>>>>> consistent for all data. I *think* that the older version of yt-3.0
>>>>>> gave
>>>>>> correct answers, though, so that should at least be a way to calibrate
>>>>>> the
>>>>>> new version.
>>>>>>
>>>>>>
>>>>>> On 13/05/2014 14:01, nick moeckel wrote:
>>>>>>
>>>>>> I had a chance to test out an older and a newer version of yt-3 on a
>>>>>> small
>>>>>> dataset that has boxlen=10. This seems to confirm a boxlen**3 factor
>>>>>> somewhere, although there's a further difference in the max density.
>>>>>>
>>>>>> more recent version (hg id -i gives f4838a2165c0):
>>>>>>
>>>>>> after dd = ds.all_data():
>>>>>>
>>>>>> In [11]: dd.quantities.extrema('density')
>>>>>>
>>>>>> Out[11]: (2.00739832795e-28 g/cm**3, 4.54771338163e-24 g/cm**3)
>>>>>>
>>>>>>
>>>>>> older version (hg id -i gives 3e8b733c9ee9):
>>>>>>
>>>>>> after dd = ds.h.all_data():
>>>>>>
>>>>>> In [10]: dd.quantities.extrema('Density')
>>>>>>
>>>>>> Out[10]: [(2.0073983279492632e-25, 4.0079147772594363e-21)]
>>>>>>
>>>>>>
>>>>>>
>>>>>> On Fri, May 9, 2014 at 5:50 PM, Sam Geen <samgeen at astro.ox.ac.uk>
>>>>>> wrote:
>>>>>>> Yep, that looks like it should work. I'll try to run it on some
>>>>>>> particle
>>>>>>> data when I get the time, but like I said I'm 99% sure the mass units
>>>>>>> should
>>>>>>> be identical for both grid hydro and particle data.
>>>>>>>
>>>>>>>
>>>>>>> On 09/05/14 17:46, Matthew Turk wrote:
>>>>>>>> Hi Sam,
>>>>>>>>
>>>>>>>> Okay, I've looked over a bit, and I think the correct change would
>>>>>>>> be:
>>>>>>>>
>>>>>>>> mass_unit = rho_u * length_unit**3.
>>>>>>>>
>>>>>>>> That should include the correct length unit, and I think will reduce
>>>>>>>> "density" back to the "unit_d" that's in the parameter file.  If this
>>>>>>>> looks okay to you, I will push it, but I really do want to make sure
>>>>>>>> the particle masses are correct.  Can Nick or Romain provide a bit of
>>>>>>>> guidance here?
>>>>>>>>
>>>>>>>> -Matt
>>>>>>>>
>>>>>>>> On Fri, May 9, 2014 at 11:25 AM, Sam Geen <samgeen at astro.ox.ac.uk>
>>>>>>>> wrote:
>>>>>>>>> OK, interesting. In theory RAMSES should have identical units for
>>>>>>>>> both
>>>>>>>>> particles and gas. I can hunt down a run with particles to test if
>>>>>>>>> you
>>>>>>>>> like.
>>>>>>>>>
>>>>>>>>> Thanks!
>>>>>>>>> Sam
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> On 09/05/14 17:22, Matthew Turk wrote:
>>>>>>>>>> Hi Sam,
>>>>>>>>>>
>>>>>>>>>> On reflection, I think this might be related to getting the
>>>>>>>>>> *particle*
>>>>>>>>>> masses correct.  I will take a look at it as soon as I can.
>>>>>>>>>>
>>>>>>>>>> -Matt
>>>>>>>>>>
>>>>>>>>>> On Fri, May 9, 2014 at 11:00 AM, Sam Geen <samgeen at astro.ox.ac.uk>
>>>>>>>>>> wrote:
>>>>>>>>>>> Hmm, that looks like it should be "mass_unit = rho_u *
>>>>>>>>>>> length_unit**3"
>>>>>>>>>>> in
>>>>>>>>>>> line 492. You're right that it mentions the boxlength issue,
>>>>>>>>>>> though.
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> On 09/05/14 16:50, Matthew Turk wrote:
>>>>>>>>>>>> Hi Sam,
>>>>>>>>>>>>
>>>>>>>>>>>> Okay, sounds good.  Looking at how code unit attributes are set
>>>>>>>>>>>> up:
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>> https://bitbucket.org/yt_analysis/yt/src/a14a150c7c81850df81346162bdaff271e77eb50/yt/frontends/ramses/data_structures.py?at=yt-3.0#cl-482
>>>>>>>>>>>>
>>>>>>>>>>>> suggests to me that length_unit takes into account boxlen, and
>>>>>>>>>>>> mass_unit does not.  The comments have some indication why this
>>>>>>>>>>>> might
>>>>>>>>>>>> be.
>>>>>>>>>>>>
>>>>>>>>>>>> -Matt
>>>>>>>>>>>>
>>>>>>>>>>>> On Fri, May 9, 2014 at 10:46 AM, Sam Geen
>>>>>>>>>>>> <samgeen at astro.ox.ac.uk>
>>>>>>>>>>>> wrote:
>>>>>>>>>>>>> Yep; both this and dd["density"] give cgs values that are too
>>>>>>>>>>>>> small
>>>>>>>>>>>>> by
>>>>>>>>>>>>> roughly a factor of boxlen**3.
>>>>>>>>>>>>>
>>>>>>>>>>>>> One other thing I need to try is to make sure I'm using the very
>>>>>>>>>>>>> latest
>>>>>>>>>>>>> version of YT; I've been playing around with the Ramses frontend
>>>>>>>>>>>>> so
>>>>>>>>>>>>> it's
>>>>>>>>>>>>> possible my version is somehow out of sync. Will let you know if
>>>>>>>>>>>>> this
>>>>>>>>>>>>> fixes
>>>>>>>>>>>>> things.
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> On 09/05/14 16:42, Matthew Turk wrote:
>>>>>>>>>>>>>> Hi Sam,
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Can you verify the units are in fact incorrect in *cgs*?
>>>>>>>>>>>>>> Something
>>>>>>>>>>>>>> like this would work:
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> ds = load(...)
>>>>>>>>>>>>>> dd = ds.all_data()
>>>>>>>>>>>>>> print dd.quantities.total_quantitiy("cell_mass").in_cgs()
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Your second message makes me wonder if there's just a slipup in
>>>>>>>>>>>>>> how
>>>>>>>>>>>>>> the units are returned.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> -Matt
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> On Fri, May 9, 2014 at 10:37 AM, Sam Geen
>>>>>>>>>>>>>> <samgeen at astro.ox.ac.uk>
>>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>>> Sorry for the spam; a second bug I've seen is that the density
>>>>>>>>>>>>>>> and
>>>>>>>>>>>>>>> pressure
>>>>>>>>>>>>>>> unit labels on figures appears to be broken; it seems to print
>>>>>>>>>>>>>>> a
>>>>>>>>>>>>>>> latex-mangled code names for the units rather than the cgs
>>>>>>>>>>>>>>> name
>>>>>>>>>>>>>>> of
>>>>>>>>>>>>>>> the
>>>>>>>>>>>>>>> units; see attached example. Temperature is fine.
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> On 09/05/14 16:28, Sam Geen wrote:
>>>>>>>>>>>>>>>> Hi,
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Just to say that in the latest version of yt-3.0 (i.e. since
>>>>>>>>>>>>>>>> various
>>>>>>>>>>>>>>>> fields were renamed or re-implemented), I've found a bug in
>>>>>>>>>>>>>>>> the
>>>>>>>>>>>>>>>> implementation of "cell_mass", which is giving results that
>>>>>>>>>>>>>>>> are
>>>>>>>>>>>>>>>> too
>>>>>>>>>>>>>>>> low
>>>>>>>>>>>>>>>> in
>>>>>>>>>>>>>>>> my runs; I believe the issue is that it's missing a factor of
>>>>>>>>>>>>>>>> pf["boxlen"]**3 (which is of course only a problem if boxlen
>>>>>>>>>>>>>>>> is
>>>>>>>>>>>>>>>> not
>>>>>>>>>>>>>>>> 1).
>>>>>>>>>>>>>>>> The
>>>>>>>>>>>>>>>> previous "CellMass_Msun" worked fine. If I get time I might
>>>>>>>>>>>>>>>> take
>>>>>>>>>>>>>>>> a
>>>>>>>>>>>>>>>> look
>>>>>>>>>>>>>>>> and
>>>>>>>>>>>>>>>> issue a pull request, but otherwise I'm just flagging this in
>>>>>>>>>>>>>>>> case
>>>>>>>>>>>>>>>> someone
>>>>>>>>>>>>>>>> else runs into problems; you can just manually multiply the
>>>>>>>>>>>>>>>> result
>>>>>>>>>>>>>>>> by
>>>>>>>>>>>>>>>> pf["boxlen"]**3 until it's fixed.
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Cheers,
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Sam
>>>>>>>>>>>>>>>> _______________________________________________
>>>>>>>>>>>>>>>> yt-users mailing list
>>>>>>>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> _______________________________________________
>>>>>>>>>>>>>>> yt-users mailing list
>>>>>>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>> _______________________________________________
>>>>>>>>>>>>>> yt-users mailing list
>>>>>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>>>>>
>>>>>>>>>>>>> _______________________________________________
>>>>>>>>>>>>> yt-users mailing list
>>>>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>>>> _______________________________________________
>>>>>>>>>>>> yt-users mailing list
>>>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> yt-users mailing list
>>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>> _______________________________________________
>>>>>>>>>> yt-users mailing list
>>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>>>
>>>>>>>>> _______________________________________________
>>>>>>>>> yt-users mailing list
>>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>> _______________________________________________
>>>>>>>> yt-users mailing list
>>>>>>>> yt-users at lists.spacepope.org
>>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> yt-users mailing list
>>>>>>> yt-users at lists.spacepope.org
>>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>
>>>>>>
>>>>>>
>>>>>> _______________________________________________
>>>>>> yt-users mailing list
>>>>>> yt-users at lists.spacepope.org
>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>
>>>>>>
>>>>>>
>>>>>> _______________________________________________
>>>>>> yt-users mailing list
>>>>>> yt-users at lists.spacepope.org
>>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>>>
>>>>> _______________________________________________
>>>>> yt-users mailing list
>>>>> yt-users at lists.spacepope.org
>>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>
>>>> _______________________________________________
>>>> yt-dev mailing list
>>>> yt-dev at lists.spacepope.org
>>>> http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
>>> _______________________________________________
>>> yt-dev mailing list
>>> yt-dev at lists.spacepope.org
>>> http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
>>
>> _______________________________________________
>> yt-dev mailing list
>> yt-dev at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
> _______________________________________________
> yt-dev mailing list
> yt-dev at lists.spacepope.org
> http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org