GalacticDynamics-Oxford / Agama

Action-based galaxy modeling framework
Other
74 stars 37 forks source link

Add a getUnits() function in py_wrapper? #15

Closed adrn closed 3 years ago

adrn commented 3 years ago

I'm working on some code that uses both Gala and Agama, but needs to construct an Agama potential on-the-fly and compare with unit-ful output from Gala, and so it would be useful if there was a way to retrieve the current units that are used by Agama. Is there a way to do this already that I missed? Since there is agama.setUnits() it might make sense to provide a agama.getUnits() to retrieve the current unit state?

eugvas commented 3 years ago

perhaps, although I doubt this is really useful in practice. the setUnits() directive should be called at the beginning of the script (or never, if the convention G=1 is used throughout), before creating any physically relevant objects such as potentials, DFs, etc. - otherwise the values returned by their methods will be hopelessly screwed. So the only reason for existence of getUnits() could be that the user has forgotten what s/he wrote at the beginning of the script! do you think this is really needed?

adrn commented 3 years ago

I see another case: I write some code that uses both Gala and Agama internally, so it needs to know what unit system Agama is working in to properly compare with Gala outputs. That code could internally call setUnits(), but that changes the global state -- what if the user had previously called setUnits() and is not expecting the values to have changed outside of the function I provide? getUnits() would allow my code to just grab the previous unit state.

eugvas commented 3 years ago

ah ok, so this could be used in a function called from the user script, but the function itself doesn't know how the units were set up. okay, so the most recent commit has getUnits(), which returns a dictionary with four dimensional units (only three are independent): mass, velocity, length, time. If astropy is installed in the system, they will be provided as astropy.Quantity (maybe makes more sense to give them as astropy.UnitBase? I am not too fluent in this to appreciate the difference..), otherwise just as plain numbers expressing the units in terms of Msun, km/s, kpc and Myr. Btw, the setUnits() function can now accept astropy quantities as well (instead of plain numbers in the same base units as listed above). But the remaining functions and classes in the library still do not accept or return astropy-compatible arrays of Quantity or whatever is in vogue these days. Agama is optimized for speed and astropy sometimes introduces too much overhead, that's why I actually created a lightweight custom coordinate conversion subsystem in agama.

adrn commented 3 years ago

Thanks @eugvas ! Returning the units as Quantity objects is fine.