Repository: yt-project/unyt
Pull request #3 (2018-04-18 02:15:12)

make objects in unyt.unit_objects Unit instances
This makes `unyt` behave more like astropy.units: the unit objects importable from the top-level `unyt` namespace become `Unit` objects, not `unyt_quantity` instances. The physical constants in `unyt.physical_constants` that are also in the top-level `unyt` namespace remain quantities. 

This change reduces the overhead of doing something like:

```python
from unyt import g
[1, 2, 3]*g
```

by as much as a factor of 4.

Unfortunately this change also triggers a number of other changes to tests (mostly in examples that import something directly from unyt and use it as a quantity or an array). I'm not sure whether that's something we should worry about given that the most common use cases continue to work and most of the examples in the tests are toy examples.


Patch:
diff --git a/docs/usage.rst b/docs/usage.rst
index 1266898..02f1633 100644
--- a/docs/usage.rst
+++ b/docs/usage.rst
@@ -98,8 +98,8 @@ python syntax. This means you must use `**` and not `^`, even when writing a
 unit as a string:
 
   >>> from unyt import kg, m
-  >>> print((kg/m**3).to('g/cm**3'))
-  0.001 g/cm**3
+  >>> print((10*kg/m**3).to('g/cm**3'))
+  0.01 g/cm**3
 
 Arithmetic and units
 --------------------
@@ -111,16 +111,14 @@ mistake in your units in a mathematical formula. To see what I mean by that,
 let's take a look at the following examples::
 
   >>> from unyt import cm, m, ft, yard
-  >>> print("{}, {}, {}, {}".format(cm, m, ft, yard))
-  1.0 cm, 1.0 m, 1.0 ft, 1.0 yd
   >>> print(3*cm + 4*m - 5*ft + 6*yard)
   799.24 cm
 
-Despite the fact that the four unit symbols used in the above example have four
-different units, :mod:`unyt` is able to automatically convert the value of all
-three units into a common unit and return the result in those units. Note that
-for expressions where the return units are ambiguous, :mod:`unyt` always returns
-data in the units of the leftmost object in an expression::
+Despite the fact that the four unit symbols used in the above example correspond
+to four different units, :mod:`unyt` is able to automatically convert the value
+of all three units into a common unit and return the result in those units. Note
+that for expressions where the return units are ambiguous, :mod:`unyt` always
+returns data in the units of the leftmost object in an expression::
 
   >>> print(4*m + 3*cm - 5*ft + 6*yard)  # doctest: +FLOAT_CMP
   7.9924 m
@@ -166,7 +164,7 @@ If you make a mistake by adding two things that have different dimensions,
 code:
 
   >>> from unyt import kg, m
-  >>> kg + m  # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
+  >>> 3*kg + 5*m  # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
   Traceback (most recent call last):
   ...
   unyt.exceptions.UnitOperationError: The <ufunc 'add'> operator for
@@ -249,14 +247,14 @@ you know which units you would like to convert it to, you can make use of the
 :meth:`unyt_array.to <unyt.array.unyt_array.to>` function:
 
   >>> from unyt import mile
-  >>> mile.to('ft')
+  >>> (1.0*mile).to('ft')
   unyt_quantity(5280., 'ft')
 
 If you try to convert to a unit with different dimensions, :mod:`unyt` will
 raise an error:
 
   >>> from unyt import mile
-  >>> mile.to('lb')  # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
+  >>> (1.0*mile).to('lb')  # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
   Traceback (most recent call last):
   ...
   unyt.exceptions.UnitConversionError: Unit dimensionalities do not match.
@@ -408,7 +406,7 @@ global registry of unit systems known to the :mod:`unyt` library. That means you
 will immediately be able to use it just like the built-in unit systems:
 
   >>> from unyt import W
-  >>> W.in_base('atomic')
+  >>> (1.0*W).in_base('atomic')
   unyt_quantity(0.59746607, 'mp*nm**2/fs**3')
 
 If you would like your unit system to include an MKS current unit
@@ -459,7 +457,7 @@ the unit conversion operation. For exmaple, let's use the ``schwarzschild``
 equivalence to calculate the mass of a black hole with a radius of one AU:
 
   >>> from unyt import AU
-  >>> AU.to('Msun', equivalence='schwarzschild')
+  >>> (1.0*AU).to('Msun', equivalence='schwarzschild')
   unyt_quantity(50658673.46804737, 'Msun')
 
 Both the methods that convert data in-place and the ones that return a copy
@@ -561,8 +559,7 @@ data from external sources:
   >>> data = np.random.random((100, 100))
   >>> data_with_units = data*g.units
   >>> data_with_units.base is data
-  False
-  >>> # huh? fixme
+  True
 
 
 Integrating :mod:`unyt` Into a Python Library
diff --git a/unyt/__init__.py b/unyt/__init__.py
index 80af694..854a040 100644
--- a/unyt/__init__.py
+++ b/unyt/__init__.py
@@ -72,7 +72,7 @@ from unyt.unit_systems import UnitSystem  # NOQA
 # constants used to *construct* a physical constant) in this namespace
 def import_quantities(module, global_namespace):
     for key, value in module.__dict__.items():
-        if isinstance(value, unyt_quantity):
+        if isinstance(value, (unyt_quantity, Unit)):
             global_namespace[key] = value
 
 
diff --git a/unyt/array.py b/unyt/array.py
index c4242c3..721831d 100644
--- a/unyt/array.py
+++ b/unyt/array.py
@@ -244,9 +244,9 @@ def _get_inp_u_unary(ufunc, inputs, out_arr=None):
 def _get_inp_u_binary(ufunc, inputs):
     inp1 = _coerce_iterable_units(inputs[0])
     inp2 = _coerce_iterable_units(inputs[1])
-    unit1 = getattr(inp1, 'units', None)
-    unit2 = getattr(inp2, 'units', None)
-    ret_class = _get_binary_op_return_class(type(inp1), type(inp2))
+    unit1 = getattr(inputs[0], 'units', None) or getattr(inp1, 'units', None)
+    unit2 = getattr(inputs[1], 'units', None) or getattr(inp2, 'units', None)
+    ret_class = _get_binary_op_return_class(type(inputs[0]), type(inputs[1]))
     if unit1 is None:
         unit1 = Unit(registry=getattr(unit2, 'registry', None))
     if unit2 is None and ufunc is not power:
@@ -324,6 +324,8 @@ def _coerce_iterable_units(input_object):
             return unyt_array(input_object)
         return input_object
     else:
+        if isinstance(input_object, Unit):
+            return unyt_quantity(1.0, input_object)
         return input_object
 
 
@@ -1046,8 +1048,8 @@ class unyt_array(np.ndarray):
         Example
         -------
         >>> from unyt import Newton, km
-        >>> print((Newton/km).in_cgs())
-        1.0 g/s**2
+        >>> print((10*Newton/km).in_cgs())
+        10.0 g/s**2
         """
         try:
             to_units = self.units.get_cgs_equivalent()
@@ -1069,7 +1071,7 @@ class unyt_array(np.ndarray):
         Example
         -------
         >>> from unyt import mile
-        >>> print(mile.in_mks())
+        >>> print((1*mile).in_mks())
         1609.344 m
         """
         try:
@@ -1176,7 +1178,7 @@ class unyt_array(np.ndarray):
         Example
         -------
         >>> from unyt import km
-        >>> km.list_equivalencies()
+        >>> (1.0*km).list_equivalencies()
         spectral: length <-> frequency <-> energy
         schwarzschild: mass <-> length
         compton: mass <-> length
@@ -1191,11 +1193,11 @@ class unyt_array(np.ndarray):
         Example
         -------
         >>> from unyt import km, keV
-        >>> km.has_equivalent('spectral')
+        >>> (1.0*km).has_equivalent('spectral')
         True
-        >>> print(km.to_equivalent('MHz', equivalence='spectral'))
+        >>> print((1*km).to_equivalent('MHz', equivalence='spectral'))
         0.29979245800000004 MHz
-        >>> print(keV.to_equivalent('angstrom', equivalence='spectral'))
+        >>> print((1*keV).to_equivalent('angstrom', equivalence='spectral'))
         12.398419315219659 angstrom
         """
         return self.units.has_equivalent(equivalence)
@@ -2011,8 +2013,8 @@ class unyt_array(np.ndarray):
                     unit = self._ufunc_registry[ufunc](u)
                 ret_class = type(self)
             elif len(inputs) == 2:
-                unit_operator = self._ufunc_registry[ufunc]
                 inps, units, ret_class = _get_inp_u_binary(ufunc, inputs)
+                unit_operator = self._ufunc_registry[ufunc]
                 if unit_operator in (_comparison_unit, _arctan2_unit):
                     inps, units = _handle_comparison_units(
                         inps, units, ufunc, ret_class)
@@ -2424,10 +2426,10 @@ def ustack(arrs, axis=0):
 def _get_binary_op_return_class(cls1, cls2):
     if cls1 is cls2:
         return cls1
-    if ((cls1 in (np.ndarray, np.matrix, np.ma.masked_array) or
+    if ((cls1 in (Unit, np.ndarray, np.matrix, np.ma.masked_array) or
          issubclass(cls1, (numeric_type, np.number, list, tuple)))):
         return cls2
-    if ((cls2 in (np.ndarray, np.matrix, np.ma.masked_array) or
+    if ((cls2 in (Unit, np.ndarray, np.matrix, np.ma.masked_array) or
          issubclass(cls2, (numeric_type, np.number, list, tuple)))):
         return cls1
     if issubclass(cls1, unyt_quantity):
diff --git a/unyt/equivalencies.py b/unyt/equivalencies.py
index 5d3d7d6..a653f12 100644
--- a/unyt/equivalencies.py
+++ b/unyt/equivalencies.py
@@ -94,9 +94,9 @@ class NumberDensityEquivalence(Equivalence):
     >>> print(NumberDensityEquivalence())
     number density: density <-> number density
     >>> from unyt import Msun, pc
-    >>> rho = Msun/pc**3
+    >>> rho = 3*Msun/pc**3
     >>> rho.to_equivalent('cm**-3', 'number_density', mu=1.4)
-    unyt_quantity(28.88289965, 'cm**(-3)')
+    unyt_quantity(86.64869896, 'cm**(-3)')
     """
     type_name = "number_density"
     _dims = (density, number_density,)
@@ -171,8 +171,8 @@ class MassEnergyEquivalence(Equivalence):
     >>> print(MassEnergyEquivalence())
     mass_energy: mass <-> energy
     >>> from unyt import g
-    >>> g.to_equivalent('J', 'mass_energy')
-    unyt_quantity(8.98755179e+13, 'J')
+    >>> (3.5*g).to_equivalent('J', 'mass_energy')
+    unyt_quantity(3.14564313e+14, 'J')
 
     """
     type_name = "mass_energy"
@@ -207,9 +207,9 @@ class SpectralEquivalence(Equivalence):
     >>> print(SpectralEquivalence())
     spectral: length <-> frequency <-> energy
     >>> from unyt import angstrom, km
-    >>> angstrom.to_equivalent('keV', 'spectral')
-    unyt_quantity(12.39841932, 'keV')
-    >>> km.to_equivalent('MHz', 'spectral')
+    >>> (3*angstrom).to_equivalent('keV', 'spectral')
+    unyt_quantity(4.13280644, 'keV')
+    >>> (1*km).to_equivalent('MHz', 'spectral')
     unyt_quantity(0.29979246, 'MHz')
     """
     type_name = "spectral"
@@ -391,9 +391,9 @@ class SchwarzschildEquivalence(Equivalence):
     >>> print(SchwarzschildEquivalence())
     schwarzschild: mass <-> length
     >>> from unyt import Msun, AU
-    >>> Msun.to_equivalent('km', 'schwarzschild')
-    unyt_quantity(2.95305543, 'km')
-    >>> AU.to_equivalent('Msun', 'schwarzschild')
+    >>> (10*Msun).to_equivalent('km', 'schwarzschild')
+    unyt_quantity(29.5305543, 'km')
+    >>> (1*AU).to_equivalent('Msun', 'schwarzschild')
     unyt_quantity(50658673.46804737, 'Msun')
     """
     type_name = "schwarzschild"
@@ -506,8 +506,8 @@ class ElectromagneticSI(Equivalence):
     Example
     -------
     >>> from unyt import gauss
-    >>> gauss.to_equivalent('T', 'SI')
-    unyt_quantity(0.0001, 'T')
+    >>> (10*gauss).to_equivalent('T', 'SI')
+    unyt_quantity(0.001, 'T')
     """
     type_name = "SI"
     alternate_names = ["si", "MKS", "mks"]
@@ -541,8 +541,8 @@ class ElectromagneticCGS(Equivalence):
     Example
     -------
     >>> from unyt import Tesla
-    >>> Tesla.to_equivalent('G', 'CGS')
-    unyt_quantity(10000., 'G')
+    >>> (10*Tesla).to_equivalent('G', 'CGS')
+    unyt_quantity(100000., 'G')
     """
     type_name = "CGS"
     alternate_names = ["cgs"]
diff --git a/unyt/unit_object.py b/unyt/unit_object.py
index 9d9e99d..56431b7 100644
--- a/unyt/unit_object.py
+++ b/unyt/unit_object.py
@@ -189,6 +189,8 @@ class Unit(Expr):
     __slots__ = ["expr", "is_atomic", "base_value", "base_offset",
                  "dimensions", "registry", "_latex_repr"]
 
+    __array_priority__ = 3.0
+
     def __new__(cls, unit_expr=sympy_one, base_value=None, base_offset=0.0,
                 dimensions=None, registry=None, latex_repr=None,
                 **assumptions):
@@ -341,6 +343,10 @@ class Unit(Expr):
         self._latex_repr = _get_latex_representation(expr, self.registry)
         return self._latex_repr
 
+    @property
+    def units(self):
+        return self
+
     # Some sympy conventions
     def __getnewargs__(self):
         return (self.expr, self.is_atomic, self.base_value, self.dimensions,
@@ -417,10 +423,14 @@ class Unit(Expr):
         if not isinstance(u, Unit):
             if isinstance(u, (numeric_type, list, tuple, np.ndarray)):
                 from unyt.array import unyt_quantity, unyt_array
-                data = np.asarray(u)
+                try:
+                    units = u.units*self
+                except AttributeError:
+                    units = self
+                data = np.asarray(u, dtype='float64')
                 if data.shape == ():
-                    return unyt_quantity(u, self)
-                return unyt_array(u, self)
+                    return unyt_quantity(data, units)
+                return unyt_array(data, units)
             else:
                 raise InvalidUnitOperation(
                     "Tried to multiply a Unit object with '%s' (type %s). "
@@ -485,7 +495,7 @@ class Unit(Expr):
         """ Take Unit to power p (float). """
         try:
             p = Rational(str(p)).limit_denominator()
-        except ValueError:
+        except (ValueError, TypeError):
             raise InvalidUnitOperation("Tried to take a Unit object to the "
                                        "power '%s' (type %s). Failed to cast "
                                        "it to a float." % (p, type(p)))
diff --git a/unyt/unit_symbols.py b/unyt/unit_symbols.py
index 7003832..d8b45ca 100644
--- a/unyt/unit_symbols.py
+++ b/unyt/unit_symbols.py
@@ -19,6 +19,7 @@ For example::
 # The full license is in the LICENSE file, distributed with this software.
 # -----------------------------------------------------------------------------
 
+from unyt.unit_object import Unit
 from unyt.array import unyt_quantity as quan
 
 #
@@ -26,227 +27,227 @@ from unyt.array import unyt_quantity as quan
 #
 
 #: femtometer
-fm = femtometer = quan(1.0, "fm")
+fm = femtometer = Unit("fm")
 #: picometer
-pm = picometer = quan(1.0, "pm")
+pm = picometer = Unit("pm")
 #: nanometer
-nm = nanometer = quan(1.0, "nm")
+nm = nanometer = Unit("nm")
 #: micrometer
-um = micrometer = quan(1.0, "um")
+um = micrometer = Unit("um")
 #: millimeter
-mm = millimeter = quan(1.0, "mm")
+mm = millimeter = Unit("mm")
 #: centimeter
-cm = centimeter = quan(1.0, "cm")
+cm = centimeter = Unit("cm")
 #: meter
-m = meter = quan(1.0, "m")
+m = meter = Unit("m")
 #: kilometer
-km = kilometer = quan(1.0, "km")
+km = kilometer = Unit("km")
 #: Megameter
-Mm = Megameter = megameter = quan(1.0, "Mm")
+Mm = Megameter = megameter = Unit("Mm")
 
 #
 # parsec
 #
 
 #: parsec
-pc = parsec = quan(1.0, "pc")
+pc = parsec = Unit("pc")
 #: kiloparsec
-kpc = kiloparsec = quan(1.0, "kpc")
+kpc = kiloparsec = Unit("kpc")
 #: megaparsec
-Mpc = mpc = megaparsec = quan(1.0, "Mpc")
+Mpc = mpc = megaparsec = Unit("Mpc")
 #: gigaparsec
-Gpc = gpc = Gigaparsec = quan(1.0, "Gpc")
+Gpc = gpc = Gigaparsec = Unit("Gpc")
 
 #
 # gram
 #
 
 #: picogram
-pg = picogram = quan(1.0, "pg")
+pg = picogram = Unit("pg")
 #: nanogram
-ng = nanogram = quan(1.0, "ng")
+ng = nanogram = Unit("ng")
 #: micogram
-ug = microgram = quan(1.0, "ug")
+ug = microgram = Unit("ug")
 #: milligram
-mg = milligram = quan(1.0, "mg")
+mg = milligram = Unit("mg")
 #: gram
-g = gram = quan(1.0, "g")
+g = gram = Unit("g")
 #: kilogram
-kg = kilogram = quan(1.0, "kg")
+kg = kilogram = Unit("kg")
 #: megagram
-Mg = megagramme = tonne = metric_ton = quan(1.0, "kg")
+Mg = megagramme = tonne = metric_ton = Unit("kg")
 
 #
 # second
 #
 
 #: femtosecond
-fs = femtoseconds = quan(1.0, "fs")
+fs = femtoseconds = Unit("fs")
 #: picosecond
-ps = picosecond = quan(1.0, "ps")
+ps = picosecond = Unit("ps")
 #: nanosecond
-ns = nanosecond = quan(1.0, "ns")
+ns = nanosecond = Unit("ns")
 #: millisecond
-ms = millisecond = quan(1.0, "ms")
+ms = millisecond = Unit("ms")
 #: second
-s = second = quan(1.0, "s")
+s = second = Unit("s")
 #: kilosecond
-ks = kilosecond = quan(1.0, "ks")
+ks = kilosecond = Unit("ks")
 #: megasecond
-Ms = megasecond = quan(1.0, "Ms")
+Ms = megasecond = Unit("Ms")
 #: gigasecond
-Gs = gigasecond = quan(1.0, "Gs")
+Gs = gigasecond = Unit("Gs")
 
 #
 # minute
 #
 
 #: minute
-min = minute = quan(1.0, "min")
+min = minute = Unit("min")
 
 #
 # hr
 #
 
 #: hour
-hr = hour = quan(1.0, "hr")
+hr = hour = Unit("hr")
 
 #
 # day
 #
 
 #: day
-day = quan(1.0, "day")
+day = Unit("day")
 
 #
 # year
 #
 
 #: year
-yr = year = quan(1.0, "yr")
+yr = year = Unit("yr")
 #: kiloyear
-kyr = kiloyear = quan(1.0, "kyr")
+kyr = kiloyear = Unit("kyr")
 #: Megayear
-Myr = Megayear = megayear = quan(1.0, "Myr")
+Myr = Megayear = megayear = Unit("Myr")
 #: Gigayear
-Gyr = Gigayear = gigayear = quan(1.0, "Gyr")
+Gyr = Gigayear = gigayear = Unit("Gyr")
 
 #
 # Temperatures
 #
 
 #: Degree kelvin
-degree_kelvin = Kelvin = K = quan(1.0, "K")
+degree_kelvin = Kelvin = K = Unit("K")
 #: Degree fahrenheit
-degree_fahrenheit = degF = quan(1.0, "degF")
+degree_fahrenheit = degF = Unit("degF")
 #: Degree Celsius
-degree_celsius = degC = quan(1.0, "degC")
+degree_celsius = degC = Unit("degC")
 #:
-degree_rankine = R = quan(1.0, "R")
+degree_rankine = R = Unit("R")
 
 #
 # Misc CGS
 #
 
 #: dyne (CGS force)
-dyne = dyn = quan(1.0, "dyne")
+dyne = dyn = Unit("dyne")
 #: erg (CGS energy)
-erg = ergs = quan(1.0, "erg")
+erg = ergs = Unit("erg")
 #:
-barye = Ba = quan(1.0, "Ba")
+barye = Ba = Unit("Ba")
 
 #
 # Misc SI
 #
 
 #: Newton (SI force)
-N = Newton = newton = quan(1.0, "N")
+N = Newton = newton = Unit("N")
 #: Joule (SI energy)
-J = Joule = joule = quan(1.0, "J")
+J = Joule = joule = Unit("J")
 #: Watt (SI power)
-W = Watt = watt = quan(1.0, "W")
+W = Watt = watt = Unit("W")
 #: Hertz
-Hz = Hertz = hertz = quan(1.0, "Hz")
+Hz = Hertz = hertz = Unit("Hz")
 #: Kilohertz
-kHz = khz = kilohertz = quan(1.0, "kHz")
+kHz = khz = kilohertz = Unit("kHz")
 #: Megahertz
-MHz = mhz = megahertz = quan(1.0, "MHz")
+MHz = mhz = megahertz = Unit("MHz")
 #: Gigahertz
-GHz = ghz = gigahertz = quan(1.0, "GHz")
+GHz = ghz = gigahertz = Unit("GHz")
 #: Terahertz
-THz = thz = terahertz = quan(1.0, "THz")
+THz = thz = terahertz = Unit("THz")
 #: Pascal
-Pa = pascal = quan(1.0, "Pa")
+Pa = pascal = Unit("Pa")
 #: sievert
-Sv = sievert = quan(1.0, "Sv")
+Sv = sievert = Unit("Sv")
 
 #
 # Imperial units
 #
 
 #: foot
-ft = foot = quan(1.0, "ft")
+ft = foot = Unit("ft")
 #: mile
-mile = quan(1.0, "mile")
+mile = Unit("mile")
 #: furlong
 furlong = quan(660, "ft")
 #: yard
-yard = yd = quan(1, "yd")
+yard = yd = Unit("yd")
 #: pound
-lb = pound = quan(1.0, "lb")
+lb = pound = Unit("lb")
 #: pound-force
-lbf = pound_force = quan(1.0, "lbf")
+lbf = pound_force = Unit("lbf")
 #: atmosphere
-atm = atmosphere = quan(1.0, "atm")
+atm = atmosphere = Unit("atm")
 #: horsepower
-hp = horsepower = quan(1.0, "hp")
+hp = horsepower = Unit("hp")
 
 #
 # Solar units
 #
 
 #: Mass of the sun
-Msun = quan(1.0, "Msun")
+Msun = Unit("Msun")
 #: Radius of the sun
-Rsun = R_sun = solar_radius = quan(1.0, "Rsun")
+Rsun = R_sun = solar_radius = Unit("Rsun")
 #: Radius of the sun
-rsun = r_sun = quan(1.0, "rsun")
+rsun = r_sun = Unit("rsun")
 #: Luminosity of the sun
-Lsun = lsun = l_sun = solar_luminosity = quan(1.0, "Lsun")
+Lsun = lsun = l_sun = solar_luminosity = Unit("Lsun")
 #: Temperature of the sun
-Tsun = T_sun = solar_temperature = quan(1.0, "Tsun")
+Tsun = T_sun = solar_temperature = Unit("Tsun")
 #: Metallicity of the sun
-Zsun = Z_sun = solar_metallicity = quan(1.0, "Zsun")
+Zsun = Z_sun = solar_metallicity = Unit("Zsun")
 
 #
 # Misc Astronomical units
 #
 
 #: Astronomical unit
-AU = astronomical_unit = quan(1.0, "AU")
+AU = astronomical_unit = Unit("AU")
 #: Astronomical unit
-au = quan(1.0, "au")
+au = Unit("au")
 #: Light year
-ly = light_year = quan(1.0, "ly")
+ly = light_year = Unit("ly")
 #: Radius of the Earth
-Rearth = R_earth = earth_radius = quan(1.0, 'R_earth')
+Rearth = R_earth = earth_radius = Unit('R_earth')
 #: Radius of the Earth
-rearth = r_earth = quan(1.0, 'r_earth')
+rearth = r_earth = Unit('r_earth')
 #: Radius of Jupiter
-Rjup = R_jup = jupiter_radius = quan(1.0, 'R_jup')
+Rjup = R_jup = jupiter_radius = Unit('R_jup')
 #: Radius of Jupiter
-rjup = r_jup = quan(1.0, 'r_jup')
+rjup = r_jup = Unit('r_jup')
 #: Jansky
-Jy = jansky = quan(1.0, "Jy")
+Jy = jansky = Unit("Jy")
 
 #
 # Physical units
 #
 
 #: electronvolt
-eV = electron_volt = electronvolt = quan(1.0, "eV")
+eV = electron_volt = electronvolt = Unit("eV")
 #: kiloelectronvolt
-keV = kilo_electron_volt = kiloelectronvolt = quan(1.0, "keV")
+keV = kilo_electron_volt = kiloelectronvolt = Unit("keV")
 #: Megaelectronvolt
 MeV = mega_electron_volt = Megaelectronvolt = megaelectronvolt = quan(
     1.0, "MeV")
@@ -254,31 +255,31 @@ MeV = mega_electron_volt = Megaelectronvolt = megaelectronvolt = quan(
 GeV = giga_electron_volt = Gigaelectronvolt = gigaelectronvolt = quan(
     1.0, "GeV")
 #: Atomic mass unit
-amu = atomic_mass_unit = quan(1.0, "amu")
-mol = quan(1.0, "mol")
+amu = atomic_mass_unit = Unit("amu")
+mol = Unit("mol")
 #: Angstrom
-angstrom = quan(1.0, "angstrom")
+angstrom = Unit("angstrom")
 #: Electron mass
-me = electron_mass = quan(1.0, "me")
+me = electron_mass = Unit("me")
 
 #
 # Angle units
 #
 
 #: Degree (angle)
-deg = degree = quan(1.0, "degree")
+deg = degree = Unit("degree")
 #: Radian
-rad = radian = quan(1.0, "radian")
+rad = radian = Unit("radian")
 #: Arcsecond
-arcsec = arcsecond = quan(1.0, "arcsec")
+arcsec = arcsecond = Unit("arcsec")
 #: Arcminute
-arcmin = arcminute = quan(1.0, "arcmin")
+arcmin = arcminute = Unit("arcmin")
 #: milliarcsecond
-mas = milliarcsecond = quan(1.0, "mas")
+mas = milliarcsecond = Unit("mas")
 #: steradian
-sr = steradian = quan(1.0, "steradian")
+sr = steradian = Unit("steradian")
 #: hourangle
-HA = hourangle = quan(1.0, "hourangle")
+HA = hourangle = Unit("hourangle")
 
 
 #
@@ -286,47 +287,47 @@ HA = hourangle = quan(1.0, "hourangle")
 #
 
 #: electrostatic unit
-electrostatic_unit = esu = quan(1.0, "esu")
+electrostatic_unit = esu = Unit("esu")
 #: Gauss
-gauss = G = quan(1.0, "gauss")
+gauss = G = Unit("gauss")
 #: Statampere
-statampere = statA = quan(1.0, "statA")
+statampere = statA = Unit("statA")
 #: Statvolt
-statvolt = statV = quan(1.0, "statV")
+statvolt = statV = Unit("statV")
 #: Statohm
-statohm = quan(1.0, "statohm")
+statohm = Unit("statohm")
 
 #
 # SI electromagnetic units
 #
 
 #: Coulomb
-C = coulomb = Coulomb = quan(1.0, "C")
+C = coulomb = Coulomb = Unit("C")
 #: Tesla
-T = tesla = Tesla = quan(1.0, "T")
+T = tesla = Tesla = Unit("T")
 #: Ampere
-A = ampere = Ampere = quan(1.0, "A")
+A = ampere = Ampere = Unit("A")
 #: Volt
-V = volt = Volt = quan(1.0, "V")
+V = volt = Volt = Unit("V")
 #: Ohm
-ohm = Ohm = quan(1.0, "ohm")
+ohm = Ohm = Unit("ohm")
 
 #
 # Geographic units
 #
 
 #: Degree latitude
-latitude = lat = quan(1.0, "lat")
+latitude = lat = Unit("lat")
 #: Degree longitude
-longitude = lon = quan(1.0, "lon")
+longitude = lon = Unit("lon")
 
 #
 # Misc dimensionless units
 #
 
 #: Count of something
-counts = count = quan(1.0, 'counts')
+counts = count = Unit('counts')
 #: Number of photons
-photons = photon = quan(1.0, 'photons')
+photons = photon = Unit('photons')
 #: dimensionless
-_ = dimensionless = quan(1.0, '')
+_ = dimensionless = Unit('')


Test patch:
diff --git a/unyt/tests/test_unit_systems.py b/unyt/tests/test_unit_systems.py
index 139251a..bf778ba 100644
--- a/unyt/tests/test_unit_systems.py
+++ b/unyt/tests/test_unit_systems.py
@@ -60,16 +60,18 @@ def test_unit_system_id():
 
 
 def test_cgs_mks_unit_conversions():
+    t = 1*Tesla
+    g = 1*gauss
     with pytest.raises(EquivalentDimsError):
-        Tesla.to('G')
-    assert Tesla.to_equivalent('G', 'cgs') == 1e4*gauss
-    assert Tesla.to_equivalent('G', 'CGS') == 1e4*gauss
+        t.to('G')
+    assert t.to_equivalent('G', 'cgs') == 1e4*gauss
+    assert t.to_equivalent('G', 'CGS') == 1e4*gauss
     with pytest.raises(EquivalentDimsError):
-        gauss.to('T')
-    assert gauss.to_equivalent("T", "mks") == 1e-4*Tesla
-    assert gauss.to_equivalent("T", "MKS") == 1e-4*Tesla
-    assert gauss.to_equivalent("T", "si") == 1e-4*Tesla
-    assert gauss.to_equivalent("T", "SI") == 1e-4*Tesla
+        g.to('T')
+    assert g.to_equivalent("T", "mks") == 1e-4*Tesla
+    assert g.to_equivalent("T", "MKS") == 1e-4*Tesla
+    assert g.to_equivalent("T", "si") == 1e-4*Tesla
+    assert g.to_equivalent("T", "SI") == 1e-4*Tesla
 
 
 def test_bad_unit_system():
diff --git a/unyt/tests/test_unyt_array.py b/unyt/tests/test_unyt_array.py
index c992d65..da927b7 100644
--- a/unyt/tests/test_unyt_array.py
+++ b/unyt/tests/test_unyt_array.py
@@ -45,6 +45,7 @@ from unyt.array import (
     savetxt
 )
 from unyt.exceptions import (
+    InvalidUnitOperation,
     UnitOperationError,
     UnitParseError,
     UfuncUnitError,
@@ -476,11 +477,13 @@ def test_power():
 
     cm_arr = np.array([1.0, 1.0]) * cm
 
-    assert_equal(cm**3, unyt_quantity(1, 'cm**3'))
-    assert_equal(np.power(cm, 3), unyt_quantity(1, 'cm**3'))
-    assert_equal(cm**unyt_quantity(3), unyt_quantity(1, 'cm**3'))
+    assert_equal((1*cm)**3, unyt_quantity(1, 'cm**3'))
+    assert_equal(np.power((1*cm), 3), unyt_quantity(1, 'cm**3'))
+    assert_equal((1*cm)**unyt_quantity(3), unyt_quantity(1, 'cm**3'))
     with pytest.raises(UnitOperationError):
-        np.power(cm, unyt_quantity(3, 'g'))
+        np.power((1*cm), unyt_quantity(3, 'g'))
+    with pytest.raises(InvalidUnitOperation):
+        np.power(cm, cm)
 
     assert_equal(cm_arr**3, unyt_array([1, 1], 'cm**3'))
     assert_equal(np.power(cm_arr, 3), unyt_array([1, 1], 'cm**3'))
@@ -1239,7 +1242,7 @@ def test_equivalencies():
 
     # frequency to wavelength
 
-    nu = u.MHz.copy()
+    nu = 1*u.MHz
     lam = nu.to('km', 'spectral')
     assert_allclose_units(lam, u.clight/nu)
     nu.convert_to_units('m', 'spectral')
@@ -1249,7 +1252,7 @@ def test_equivalencies():
 
     # frequency to photon energy
 
-    nu = u.MHz.copy()
+    nu = 1*u.MHz
     E = nu.to('erg', 'spectral')
     assert_allclose_units(E, u.hcgs*nu)
     nu.convert_to_units('J', 'spectral')
@@ -1347,7 +1350,7 @@ def test_equivalencies():
 
     # Schwarzschild
 
-    msun = u.mass_sun_cgs.copy()
+    msun = 1*u.Msun
     msun.convert_to_equivalent('km', 'schwarzschild')
     R = u.mass_sun_cgs.in_units("kpc", "schwarzschild")
     assert_allclose_units(msun, R)
@@ -1360,7 +1363,7 @@ def test_equivalencies():
 
     # Compton
 
-    me = u.me.copy()
+    me = 1*u.me
     me.convert_to_units('nm', 'compton')
     length = u.me.in_units("angstrom", "compton")
     assert_allclose_units(length, me)
@@ -1431,7 +1434,7 @@ def test_electromagnetic():
     qp_mks = u.qp.in_units("C", "SI")
     assert_equal(qp_mks.units.dimensions, dimensions.charge_mks)
     assert_almost_equal(qp_mks.v, 10.0*u.qp.v/speed_of_light_cm_per_s)
-    qp = u.qp.copy()
+    qp = 1*u.qp
     qp.convert_to_units("C", "SI")
     assert_equal(qp.units.dimensions, dimensions.charge_mks)
     assert_almost_equal(qp.v, 10*u.qp.v/u.clight.v)
@@ -1450,7 +1453,7 @@ def test_electromagnetic():
     qp_mks_k = u.qp.in_units("kC", "SI")
     assert_array_almost_equal(
         qp_mks_k.v, 1.0e-2*u.qp.v/speed_of_light_cm_per_s)
-    qp = u.qp.copy()
+    qp = 1*u.qp
     qp.convert_to_units('kC', 'SI')
     assert_almost_equal(qp, qp_mks_k)