BaPSF / bapsflib

A Python Toolkit for BaPSF
12 stars 8 forks source link

Use astropy.units subpackage to represent quantities with units? #26

Open namurphy opened 5 years ago

namurphy commented 5 years ago

It's exciting to find an open source Python software analysis package for big plasma physics experiments! We as a community need to do more of this.

It looks like many of the physical values in this package are stored as floats or arrays. There exist a few Python packages to represent quantities with units. For PlasmaPy, we use astropy.units for compatibility with Astropy, SunPy, and other packages in the Python ecosystems for astronomy and heliophysics. The use of a units package has helped immensely and makes it a lot easier for people to be able to choose between SI or cgs units. Here's an example:

>>> from astropy import units as u
>>> distance = 5 * u.m
>>> distance = 88 * u.imperial.mi
>>> distance = 88 * u.imperial.mile
>>> time = 1 * u.hour
>>> velocity = distance / time
>>> print(velocity)
88.0 mi / h
>>> velocity.to('m/s')
<Quantity 39.33952 m / s>
>>> distance + time  # using astropy.units acts as an additional error check
astropy.units.core.UnitConversionError: Can only apply 'add' function to quantities with compatible dimensions

The astropy.constants subpackage has most of the fundamental physical constants that we need for plasma physics too.

There's also a really useful decorator, astropy.units.quantity_input, that uses annotations in functions to check that the units are physically compatible. One caveat is that last I checked (sometime in 2018), it didn't work well if you have some arguments that are supposed to be astropy.units.Quantity instances and some that aren't.

from astropy import units as u

@u.quantity_input
def get_velocity(distance: u.m, time: u.s) -> u.m / u.s:
    return distance / time
namurphy commented 5 years ago

Oh, I forgot to mention something really important for us plasma physicists. We can use eV as temperature if and only if we use their equivalencies feature.

>>> kT = 1 * u.eV
>>> kT.to('K')
astropy.units.core.UnitConversionError: 'eV' (energy) and 'K' (temperature) are not convertible
>>> kT.to('K', equivalencies=u.temperature_energy())
<Quantity 11604.5220604 K>
>>> 
rocco8773 commented 5 years ago

Hi @namurphy! Thanks for adding your thoughts and advice on the package. I completely agree that we need more of these packages developed for the big machines across the field. The work you and others have done on PlasmaPy has been a great step and addition to the field.

We are planning to integrate astropy.units into the bapsflib.plasma package and the bapsflib package as a whole. I discovered astropy.units a little late in the development, so it has only been implement in the bapsflib.lapd package.