[yt-dev] Issue #1055: RAMSES units error with boxlen != 1 (yt_analysis/yt)

[Ji-hoon Kim]

New issue 1055: RAMSES units error with boxlen != 1
https://bitbucket.org/yt_analysis/yt/issues/1055/ramses-units-error-with-boxlen-1

Ji-hoon Kim:

Dear All,

The current way of assigning lengh_unit and density_unit doesn't seem to properly work for a RAMSES dataset with boxlen != 1.  FYI, I am testing yt-3.2 on the Ramses entry for the AGORA isolated disk comparison (no RT).  As I understand, yt-dev has been working fine for this particular dataset until PR #1610 is merged. I think part of the problem -- though maybe not the entirety of the problem -- is that "unit_l" gets multiplied by boxlen twice (frontend/ramses/data_structures.py; first in line 498, then in 521).  I think this alone may warrant a reinvestigation of the RAMSES units in yt.  

```
#!python
       length_unit = self.parameters['unit_l'] * boxlen
       density_unit = self.parameters['unit_d']/ boxlen**3
       mass_unit = density_unit * length_unit**3 
...
       self.density_unit = self.quan(density_unit, 'g/cm**3')
       self.length_unit = self.quan(length_unit * boxlen, "cm")
       self.mass_unit = self.quan(mass_unit, "g")
```

I tried seemingly the most obvious fix of removing the second boxlen multiplication, but it fails again with too low density (off by boxlen**3).  

```
#!python
>>> pf.all_data()["density"]()
YTArray([  2.87245413e-38,   2.87426855e-38,   2.87426181e-38, ...,
        3.07532664e-38,   3.07330200e-38,   3.07943311e-38]) g/cm**3
>>> pf.all_data()["density"].median()
3.09606619982e-33 g/cm**3
```

FWIW, I did learn about the history behind the data_structures.py by reading through previous conversations dating back last year -- among Ricarda, Nathan, Desika, Ben, Sam, Matt, and others.  Especially after going through its lengthy history and multiple changes and reverts, I wasn't certain that I could *find and verify* the correct way of setting units for all Ramses users and datasets out there.  But, once Ricarda and a broader yt-Ramses user base come up with the **most general** way of assigning units, I will be more than happy to test the fix on the dataset I have.    

Many thanks as always to the yt community!

Best wishes,

Ji-hoon

* * *

PS.  Below are just my two cents.  

IMHO, *as far as length_unit and density_unit are concerned*, I think pre-PR#1610-era solution is more consistent with pymses (snippet in https://bitbucket.org/yt_analysis/yt/issues/939/units-error-in-ramses-frontend).  Please see below for comparison.   

However, there is one important outstanding issue.  Judging from the comments in the pymses snippet, the “unit_mass” in pymses seems *only* for DM particles, and boxlen renormalization is *only* applied to “unit_mass".  (Am I correct here?)   In contrast, “mass_unit” in yt could be for any entity including star particles.  This might explain the issue Desika originally reported with stellar masses in another AGORA-disk run (http://lists.spacepope.org/pipermail/yt-dev-spacepope.org/2015-June/019345.html).  Unfortunately, however, this is a pure speculation based *not* on extensive exposures to various Ramses datasets, *but* only on my reading of the pymses snippet.  

[A]  Pre-PR#1610-era length_unit and density_unit (e.g. changeset 242b9f739613):  

```
#!python
       length_unit = self.parameters['unit_l']
       density_unit = self.parameters['unit_d']
       mass_unit = density_unit * (length_unit * boxlen)**3
...
       self.density_unit = self.quan(density_unit, 'g/cm**3')
       self.length_unit = self.quan(length_unit * boxlen, "cm")
       self.mass_unit = self.quan(mass_unit, "g")
```

[B]  Compare with the pymses code snippet:

```
#!python
   info_dict["unit_length"] = float(par_dict["unit_l"])*C.cm
   info_dict["unit_density"] = float(par_dict["unit_d"])*C.g/C.cm**3
...
   # Boxlen renormalisation : used only for dark matter particles,
   # so we do this just before computing the mass unit
   info_dict["unit_length"] = info_dict["unit_length"] * info_dict["boxlen"]
   info_dict["unit_mass"] = info_dict["unit_density"] \
                   * info_dict["unit_length"]**3
```


I am sure that Ricarda and other experienced Ramses users will be better positioned to suggest the *most general* solution.  I look forward to hearing their insights.  (Questions I couldn't quite answer for myself: Would any suggested changes cause troubles in other units such as velocity, pressure, or B-field?  If my reading of the pymses code is correct, how could yt possibly assign two mass_units?)