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

Matthew Turk matthewturk at gmail.com
Wed May 14 05:22:02 PDT 2014 

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

More information about the yt-dev mailing list