Repository: yt-project/unyt
Pull request #6 (2018-05-11 20:19:29)

almost 100% test coverage, remove MKS/CGS equivalences in favor of direct conversions
I realize this is a huge pull request. I can definitely break this up if that would help.

This gets us to 3 lines away from 100% test coverage. The remaining untested lines are either too difficult to test (e.g. I'd need to implement a dummy numpy ufunc) or are in code i've copy/pasted from sympy and I'd prefer to leave unchanged for simplicity.

There's a number of places where I've deleted code. It may look like those places were adding functionality but according to the test coverage those lines were never actually executed. In all cases I made a good-faith effort to try to trigger the deleted lines and usually I found that logic elsewhere made it impossible to reach that line.

I removed `UnitSystemConstants` since as far as I could see it wasn't actually used anywhere. If there was an intended purpose for that I'm happy to add it back along with a test and some docs.

The `Unit` class is no longer a subclass of sympy's `Expr` class. This makes the code simpler to understand by making it clearer where the sympy "magic" happen (in the `expr` attribute of the `Unit` class, not in the `Unit` class itself).


Patch:
diff --git a/docs/usage.rst b/docs/usage.rst
index db0a142..4374426 100644
--- a/docs/usage.rst
+++ b/docs/usage.rst
@@ -257,8 +257,8 @@ raise an error:
   >>> (1.0*mile).to('lb')  # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
   Traceback (most recent call last):
   ...
-  unyt.exceptions.UnitConversionError: Unit dimensionalities do not match.
-  Tried to convert between mile (dim (length)) and lb (dim (mass)).
+  unyt.exceptions.UnitConversionError: Cannot convert between mile (dim
+  (length)) and lb (dim (mass)).
 
 While we recommend using :meth:`unyt_array.to <unyt.array.unyt_array.to>` in
 most cases to convert arrays or quantities to different units, if you would like
@@ -351,13 +351,30 @@ system. In CGS units the electromagnetic units like Gauss and statA are
 decomposable in terms of the base mass, length, and time units in the unit
 system. For this reason quantities defined in E&M units in CGS units are not
 readily convertible to MKS units and vice verse since the units are not
-dimensionally equivalent. To resolve this, :mod:`unyt` provides a unit
-equivalency system, discussed below, to convert data between semantically
-equivalent but not dimensionally equal units.
+dimensionally equivalent. The :mod:`unyt` library does have limited support for converting electromagnetic units between MKS and CGS, however only simple conversions of data with a single specific unit are supported and no conversions are allowed for complex combinations of units. For example converting between Gauss and Tesla is supported:
 
-The names ``"SI"``, ``"si"``, and ``"MKS"`` are accepted as alternative names
-by :mod:`unyt` for the MKS unit system. Similarly, ``"CGS"`` is acceptable as
-a name for the CGS unit system.
+  >>> from unyt import T
+  >>> (1.0*T).to('G')
+  unyt_quantity(10000., 'G')
+
+But converting a more complicated compound unit will raise an error:
+
+  >>> from unyt import C, T, V
+  >>> (1.0*C*T*V).in_cgs()  # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
+  Traceback (most recent call last):
+  ...
+  unyt.exceptions.UnitsNotReducible: The unit "C*T*V" (dimensions
+  "(length)**2*(mass)**2/((current_mks)*(time)**4)" cannot be reduced to
+  an expression within the cgs system of units.
+
+If you need to work with complex expressions involving electromagnetic units, we
+suggest sticking to either CGS or SI units for the full calculation. There is no
+general way to convert an arbitrary quantity between CGS and SI units if the
+quantity involves electromagnetic units. Instead, it is necessary to do the
+conversion on the equations under consideration, and then recompute the
+necessary quantity in the transformed set of equations. This requires
+understanding the context for a calculation, which unfortunately is beyond the
+scope of a library like :mod:`unyt`.
 
 You can convert data to a unit system :mod:`unyt` knows about using the
 :meth:`unyt_array.in_base <unyt.array.unyt_array.in_base>` and
@@ -401,7 +418,8 @@ atoms and molecules:
      luminous_intensity: cd
     Other Units:
      energy: eV
-
+   >>> print(atomic_unit_system)
+   atomic
    >>> atomic_unit_system['number_density']
    nm**(-3)
    >>> atomic_unit_system['angular_momentum']
diff --git a/tox.ini b/tox.ini
index f2fab89..30f2ce1 100644
--- a/tox.ini
+++ b/tox.ini
@@ -1,5 +1,5 @@
 [tox]
-envlist = docs,begin,py27-min,py{27,34,35,36},end
+envlist = docs,begin,py27-min,py27-mindeps,py{27,34,35,36},end
 
 [travis]
 os =
@@ -9,7 +9,7 @@ python =
     3.6: docs, py36
     3.5: py35
     3.4: py34
-    2.7: py27, py27-min
+    2.7: py27, py27-min, py27-mindeps
 
 [testenv]
 setenv =
@@ -55,6 +55,21 @@ commands =
     pytest --cov=unyt --cov-append --basetemp={envtmpdir}
     coverage report --omit='.tox/*'
 
+[testenv:py27-mindeps]
+deps =
+    pytest
+    sympy
+    numpy
+    six
+    coverage
+    pytest-cov
+    pytest-doctestplus
+    backports.functools_lru_cache
+commands =
+    # don't do doctests on old numpy versions
+    pytest --cov=unyt --cov-append --basetemp={envtmpdir}
+    coverage report --omit='.tox/*'
+
 [testenv:docs]
 changedir = docs
 deps = sphinx
diff --git a/unyt/array.py b/unyt/array.py
index a978ea2..669ca86 100644
--- a/unyt/array.py
+++ b/unyt/array.py
@@ -14,9 +14,13 @@ unyt_array class.
 # -----------------------------------------------------------------------------
 
 import copy
-import numpy as np
-
 from functools import wraps
+try:
+    from functools import lru_cache
+except ImportError:
+    from backports.functools_lru_cache import lru_cache
+from numbers import Number as numeric_type
+import numpy as np
 from numpy import (
     add,
     subtract,
@@ -101,41 +105,36 @@ from numpy import (
     divmod as divmod_,
     isnat,
     heaviside,
+    ones_like,
 )
 from numpy.core.umath import _ones_like
+from sympy import Rational
 
-from unyt.unit_object import (
-    Unit,
-    UnitParseError
-)
-from unyt.unit_registry import UnitRegistry
 from unyt.dimensions import (
     angle,
-    current_mks,
-    em_dimensions,
     temperature
 )
 from unyt.exceptions import (
-    EquivalentDimsError,
     IterableUnitCoercionError,
     InvalidUnitEquivalence,
     InvalidUnitOperation,
+    MKSCGSConversionError,
     UnitConversionError,
     UnitOperationError,
+    UnitsNotReducible,
 )
-from unyt.equivalencies import (
-    equivalence_registry,
-    _check_em_conversion,
+from unyt.equivalencies import equivalence_registry
+from unyt._on_demand_imports import (
+    _astropy,
+    _pint,
 )
-try:
-    from functools import lru_cache
-except ImportError:
-    from backports.functools_lru_cache import lru_cache
-from numbers import Number as numeric_type
-from unyt._on_demand_imports import _astropy
-from sympy import Rational
-from unyt._unit_lookup_table import default_unit_symbol_lut
 from unyt._pint_conversions import convert_pint_units
+from unyt.unit_object import (
+    _check_em_conversion,
+    _em_conversion,
+    Unit,
+)
+from unyt.unit_registry import UnitRegistry
 
 NULL_UNIT = Unit()
 POWER_SIGN_MAPPING = {multiply: 1, divide: -1}
@@ -232,11 +231,7 @@ def _coerce_iterable_units(input_object):
                 raise IterableUnitCoercionError(input_object)
             # This will create a copy of the data in the iterable.
             return unyt_array(input_object)
-        return np.asarray(input_object)
-    else:
-        if isinstance(input_object, Unit):
-            return unyt_quantity(1.0, input_object)
-        return np.asarray(input_object)
+    return np.asarray(input_object)
 
 
 @lru_cache(maxsize=128, typed=False)
@@ -246,14 +241,6 @@ def _unit_repr_check_same(my_units, other_units):
     is compatible with this quantity. Returns Unit object.
 
     """
-    equiv_dims = em_dimensions.get(my_units.dimensions, None)
-    if equiv_dims == other_units.dimensions:
-        if current_mks in equiv_dims.free_symbols:
-            base = "SI"
-        else:
-            base = "CGS"
-        raise EquivalentDimsError(my_units, other_units, base)
-
     if not my_units.same_dimensions_as(other_units):
         raise UnitConversionError(
             my_units, my_units.dimensions, other_units, other_units.dimensions)
@@ -318,6 +305,7 @@ unary_operators = (
     spacing,
     positive,
     isnat,
+    ones_like,
 )
 
 binary_operators = (
@@ -516,15 +504,6 @@ class unyt_array(np.ndarray):
             if registry is not None:
                 obj.units.registry = registry
             return obj
-        if input_array is NotImplemented:
-            return input_array.view(cls)
-        if registry is None and isinstance(input_units, (str, bytes)):
-            if input_units.startswith('code_'):
-                raise UnitParseError(
-                    "Code units used without referring to a dataset. \n"
-                    "Perhaps you meant to do something like this instead: \n"
-                    "ds.arr(%s, \"%s\")" % (input_array, input_units)
-                    )
         if isinstance(input_array, unyt_array):
             ret = input_array.view(cls)
             if input_units is None:
@@ -618,17 +597,20 @@ class unyt_array(np.ndarray):
         >>> print(length)
         [30. 20. 10.] m
         """
+        units = _sanitize_units_convert(units, self.units.registry)
         if equivalence is None:
-            units = _sanitize_units_convert(units, self.units.registry)
-            em_units, em_us = _check_em_conversion(self.units)
-            if em_units is not None:
-                em_units = _sanitize_units_convert(em_units,
-                                                   self.units.registry)
-                if em_units.dimensions == units.dimensions:
-                    return self.convert_to_equivalent(units, em_us)
-            new_units = _unit_repr_check_same(self.units, units)
-            (conversion_factor, offset) = self.units.get_conversion_factor(
-                new_units)
+            try:
+                conv_data = _check_em_conversion(self.units, units)
+            except MKSCGSConversionError:
+                raise UnitConversionError(self.units, self.units.dimensions,
+                                          units, units.dimensions)
+            if any(conv_data):
+                new_units, (conversion_factor, offset) = _em_conversion(
+                    self.units, conv_data, units)
+            else:
+                new_units = _unit_repr_check_same(self.units, units)
+                (conversion_factor, offset) = self.units.get_conversion_factor(
+                    new_units)
 
             self.units = new_units
             values = self.d
@@ -779,15 +761,19 @@ class unyt_array(np.ndarray):
         """
         units = _sanitize_units_convert(units, self.units.registry)
         if equivalence is None:
-            em_units, em_us = _check_em_conversion(self.units)
-            if em_units is not None:
-                em_units = _sanitize_units_convert(em_units,
-                                                   self.units.registry)
-                if em_units.dimensions == units.dimensions:
-                    return self.to_equivalent(units, em_us)
-            new_units = _unit_repr_check_same(self.units, units)
-            (conversion_factor, offset) = self.units.get_conversion_factor(
-                new_units)
+            try:
+                conv_data = _check_em_conversion(self.units, units)
+            except MKSCGSConversionError:
+                raise UnitConversionError(self.units, self.units.dimensions,
+                                          units, units.dimensions)
+            if any(conv_data):
+                new_units, (conversion_factor, offset) = _em_conversion(
+                    self.units, conv_data, units)
+                offset = 0
+            else:
+                new_units = _unit_repr_check_same(self.units, units)
+                (conversion_factor, offset) = self.units.get_conversion_factor(
+                    new_units)
 
             ret = np.asarray(self.ndview * conversion_factor)
             if offset:
@@ -904,14 +890,20 @@ class unyt_array(np.ndarray):
         >>> print(E.in_base("mks"))
         2.5e-07 W
         """
-        em_units, em_us = _check_em_conversion(self.units)
-        if em_units is not None and em_us == unit_system:
-            to_units = em_units
-            equivalence = em_us
+        try:
+            conv_data = _check_em_conversion(self.units)
+        except MKSCGSConversionError:
+            raise UnitsNotReducible(self.units, unit_system)
+        if any(conv_data):
+            to_units, (conv, offset) = _em_conversion(
+                self.units, conv_data, unit_system=unit_system)
         else:
             to_units = self.units.get_base_equivalent(unit_system)
-            equivalence = None
-        return self.in_units(to_units, equivalence=equivalence)
+            conv, offset = self.units.get_conversion_factor(to_units)
+        ret = self.v*conv
+        if offset:
+            ret = ret - offset
+        return type(self)(ret, to_units)
 
     def in_cgs(self):
         """
@@ -974,21 +966,11 @@ class unyt_array(np.ndarray):
         conv_unit = Unit(unit, registry=self.units.registry)
         if self.units.same_dimensions_as(conv_unit):
             self.convert_to_units(conv_unit)
+            return
         this_equiv = equivalence_registry[equivalence](in_place=True)
-        oneway_or_equivalent = (conv_unit.has_equivalent(equivalence)
-                                or this_equiv.one_way)
-        if self.has_equivalent(equivalence) and oneway_or_equivalent:
-            ret = this_equiv._convert(self, conv_unit.dimensions, **kwargs)
-            if isinstance(ret, tuple):
-                # need to subtract off the current value because index
-                # into an array scalar isn't allowed
-                view = self.d
-                view += (ret[0] - view)
-                self.units = Unit(ret[1], registry=self.units.registry)
-            try:
-                self.convert_to_units(conv_unit)
-            except UnitConversionError:
-                InvalidUnitEquivalence(equivalence, self.units, unit)
+        if self.has_equivalent(equivalence):
+            this_equiv.convert(self, conv_unit.dimensions, **kwargs)
+            self.convert_to_units(conv_unit)
         else:
             raise InvalidUnitEquivalence(equivalence, self.units, unit)
 
@@ -1019,19 +1001,10 @@ class unyt_array(np.ndarray):
         if self.units.same_dimensions_as(conv_unit):
             return self.in_units(conv_unit)
         this_equiv = equivalence_registry[equivalence]()
-        oneway_or_equivalent = (
-            conv_unit.has_equivalent(equivalence) or this_equiv.one_way)
-        if self.has_equivalent(equivalence) and oneway_or_equivalent:
-            new_arr = this_equiv._convert(
+        if self.has_equivalent(equivalence):
+            new_arr = this_equiv.convert(
                 self, conv_unit.dimensions, **kwargs)
-            if isinstance(new_arr, tuple):
-                try:
-                    return type(self)(new_arr[0], new_arr[1]).in_units(
-                        conv_unit)
-                except UnitConversionError:
-                    raise InvalidUnitEquivalence(equivalence, self.units, unit)
-            else:
-                return new_arr.in_units(conv_unit)
+            return new_arr.in_units(conv_unit)
         else:
             raise InvalidUnitEquivalence(equivalence, self.units, unit)
 
@@ -1178,10 +1151,6 @@ class unyt_array(np.ndarray):
         >>> data.to_astropy()
         <Quantity [3., 4., 5.] g / cm3>
         """
-        if _astropy.units is None:
-            raise ImportError(
-                "You don't have AstroPy installed, so you can't convert to " +
-                "an AstroPy quantity.")
         return self.value*_astropy.units.Unit(str(self.units), **kwargs)
 
     @classmethod
@@ -1244,9 +1213,8 @@ class unyt_array(np.ndarray):
         >>> a.to_pint()
         <Quantity(4.0, 'centimeter ** 2 / second')>
         """
-        from pint import UnitRegistry
         if unit_registry is None:
-            unit_registry = UnitRegistry()
+            unit_registry = _pint.UnitRegistry()
         powers_dict = self.units.expr.as_powers_dict()
         units = []
         for unit, pow in powers_dict.items():
@@ -1359,10 +1327,7 @@ class unyt_array(np.ndarray):
         dataset = g[dataset_name]
         data = dataset[:]
         units = dataset.attrs.get('units', '')
-        if 'unit_registry' in dataset.attrs.keys():
-            unit_lut = pickle.loads(dataset.attrs['unit_registry'].tostring())
-        else:
-            unit_lut = None
+        unit_lut = pickle.loads(dataset.attrs['unit_registry'].tostring())
         f.close()
         registry = UnitRegistry(lut=unit_lut, add_default_symbols=False)
         return cls(data, units, registry=registry)
@@ -1577,8 +1542,6 @@ class unyt_array(np.ndarray):
             # Unary ufuncs
             inp = inputs[0]
             u = getattr(inp, 'units', None)
-            if u is None:
-                u = NULL_UNIT
             if u.dimensions is angle and ufunc in trigonometric_operators:
                 # ensure np.sin(90*degrees) works as expected
                 inp = inp.in_units('radian').v
@@ -1661,15 +1624,8 @@ class unyt_array(np.ndarray):
                         raise InvalidUnitOperation(
                             "Quantities with units of Fahrenheit or Celsius "
                             "cannot by multiplied, divided, subtracted or "
-                            "added.")
+                            "added with data that has different units.")
                     inp1 = np.asarray(inp1)*conv
-                else:
-                    if ((u0.base_offset and u0.dimensions is temperature or
-                         u1.base_offset and u1.dimensions is temperature)):
-                        raise InvalidUnitOperation(
-                            "Quantities with units of Fahrenheit or Celsius "
-                            "cannot by multiplied, divide, subtracted or "
-                            "added.")
             # get the unit of the result
             unit = unit_operator(u0, u1)
             # actually evaluate the ufunc
@@ -1756,8 +1712,6 @@ class unyt_array(np.ndarray):
         return type(self)(np.copy(np.asarray(self)), self.units)
 
     def __array_finalize__(self, obj):
-        if obj is None and hasattr(self, 'units'):
-            return
         self.units = getattr(obj, 'units', NULL_UNIT)
 
     def __pos__(self):
@@ -1817,19 +1771,7 @@ class unyt_array(np.ndarray):
         metadata extracted in __reduce__ and then serialized by pickle.
         """
         super(unyt_array, self).__setstate__(state[1:])
-        try:
-            unit, lut = state[0]
-        except TypeError:
-            # this case happens when we try to load an old pickle file
-            # created before we serialized the unit symbol lookup table
-            # into the pickle file
-            unit, lut = str(state[0]), default_unit_symbol_lut.copy()
-        # need to fix up the lut if the pickle was saved prior to PR #1728
-        # when the pickle format changed
-        if len(lut['m']) == 2:
-            lut.update(default_unit_symbol_lut)
-            for k, v in [(k, v) for k, v in lut.items() if len(v) == 2]:
-                lut[k] = v + (0.0, r'\rm{' + k.replace('_', '\ ') + '}')
+        unit, lut = state[0]
         registry = UnitRegistry(lut=lut, add_default_symbols=False)
         self.units = Unit(unit, registry=registry)
 
@@ -1838,8 +1780,6 @@ class unyt_array(np.ndarray):
 
         This is necessary for stdlib deepcopy of arrays and quantities.
         """
-        if memodict is None:
-            memodict = {}
         ret = super(unyt_array, self).__deepcopy__(memodict)
         return type(self)(ret, copy.deepcopy(self.units))
 
@@ -2086,6 +2026,14 @@ def ustack(arrs, axis=0):
     This is a wrapper around np.stack that preserves units. See the
     documentation for np.stack for full details.
 
+    Examples
+    --------
+    >>> from unyt import km
+    >>> a = [1, 2, 3]*km
+    >>> b = [2, 3, 4]*km
+    >>> print(ustack([a, b]))
+    [[1. 2. 3.]
+     [2. 3. 4.]] km
     """
     v = np.stack(arrs)
     v = _validate_numpy_wrapper_units(v, arrs)
@@ -2157,23 +2105,25 @@ def loadtxt(fname, dtype='float', delimiter='\t', usecols=None, comments='#'):
                 next_one = True
         else:
             # Here we catch the first line of numbers
-            try:
-                col_words = line.strip().split(delimiter)
-                for word in col_words:
-                    float(word)
-                num_cols = len(col_words)
-                break
-            except ValueError:
-                pass
+            col_words = line.strip().split(delimiter)
+            for word in col_words:
+                float(word)
+            num_cols = len(col_words)
+            break
     f.close()
     if len(units) != num_cols:
         units = ["dimensionless"]*num_cols
     arrays = np.loadtxt(fname, dtype=dtype, comments=comments,
                         delimiter=delimiter, converters=None,
                         unpack=True, usecols=usecols, ndmin=0)
+    if len(arrays.shape) < 2:
+        arrays = [arrays]
     if usecols is not None:
         units = [units[col] for col in usecols]
-    return tuple([unyt_array(arr, unit) for arr, unit in zip(arrays, units)])
+    ret = tuple([unyt_array(arr, unit) for arr, unit in zip(arrays, units)])
+    if len(ret) == 1:
+        return ret[0]
+    return ret
 
 
 def savetxt(fname, arrays, fmt='%.18e', delimiter='\t', header='',
@@ -2218,7 +2168,7 @@ def savetxt(fname, arrays, fmt='%.18e', delimiter='\t', header='',
             units.append(str(array.units))
         else:
             units.append("dimensionless")
-    if header != '':
+    if header != '' and not header.endswith('\n'):
         header += '\n'
     header += " Units\n " + '\t'.join(units)
     np.savetxt(fname, np.transpose(arrays), header=header,
diff --git a/unyt/equivalencies.py b/unyt/equivalencies.py
index 452cf1f..8e0a26e 100644
--- a/unyt/equivalencies.py
+++ b/unyt/equivalencies.py
@@ -24,24 +24,9 @@ from unyt.dimensions import (
     density,
     number_density,
     flux,
-    current_cgs,
-    current_mks,
-    charge_cgs,
-    charge_mks,
-    magnetic_field_cgs,
-    magnetic_field_mks,
-    magnetic_flux_cgs,
-    magnetic_flux_mks,
-    electric_potential_cgs,
-    electric_potential_mks,
-    electric_field_cgs,
-    electric_field_mks,
-    resistance_cgs,
-    resistance_mks
 )
+from unyt.exceptions import InvalidUnitEquivalence
 
-from unyt._physical_ratios import speed_of_light_cm_per_s
-from unyt._unit_lookup_table import unit_prefixes
 from six import add_metaclass
 import numpy as np
 
@@ -53,18 +38,19 @@ class _RegisteredEquivalence(type):
         type.__init__(cls, name, b, d)
         if hasattr(cls, "type_name"):
             equivalence_registry[cls.type_name] = cls
-        if hasattr(cls, "alternate_names"):
-            for name in cls.alternate_names:
-                equivalence_registry[name] = cls
 
 
 @add_metaclass(_RegisteredEquivalence)
 class Equivalence(object):
-    one_way = False
-
     def __init__(self, in_place=False):
         self.in_place = in_place
 
+    def convert(self, x, new_dims, **kwargs):
+        if x.units.dimensions in self._dims and new_dims in self._dims:
+            return self._convert(x, new_dims, **kwargs)
+        else:
+            raise InvalidUnitEquivalence(self, x, new_dims)
+
     def _get_out(self, x):
         if self.in_place:
             return x
@@ -72,7 +58,8 @@ class Equivalence(object):
 
 
 class NumberDensityEquivalence(Equivalence):
-    """Equivalence between mass and number density, given a mean molecular weight.
+    """Equivalence between mass and number density, given a mean molecular
+    weight.
 
     Given a number density :math:`n`, the mass density :math:`\\rho` is:
 
@@ -444,7 +431,7 @@ class ComptonEquivalence(Equivalence):
         return "compton: mass <-> length"
 
 
-class EffectiveTemperature(Equivalence):
+class EffectiveTemperatureEquivalence(Equivalence):
     """Equivalence between the emmitted flux accross all wavelengths and
     temperature of a blackbody
 
@@ -459,6 +446,8 @@ class EffectiveTemperature(Equivalence):
 
     Example
     -------
+    >>> print(EffectiveTemperatureEquivalence())
+    effective_temperature: flux <-> temperature
     >>> from unyt import K, W, m
     >>> (5000*K).to_equivalent('W/m**2', 'effective_temperature')
     unyt_quantity(35439831.25, 'W/m**2')
@@ -481,119 +470,3 @@ class EffectiveTemperature(Equivalence):
 
     def __str__(self):
         return "effective_temperature: flux <-> temperature"
-
-
-em_conversions = {
-    "C": ("esu", 0.1*speed_of_light_cm_per_s, "cgs"),
-    "T": ("gauss", 1.0e4, "cgs"),
-    "Wb": ("Mx", 1.0e8, "cgs"),
-    "A": ("statA", 0.1*speed_of_light_cm_per_s, "cgs"),
-    "V": ("statV", 1.0e-8*speed_of_light_cm_per_s, "cgs"),
-    "ohm": ("statohm", 1.0e9/(speed_of_light_cm_per_s**2), "cgs"),
-    "esu": ("C", 10.0/speed_of_light_cm_per_s, "mks"),
-    "Fr": ("C", 10.0/speed_of_light_cm_per_s, "mks"),
-    "statC": ("C", 10.0/speed_of_light_cm_per_s, "mks"),
-    "gauss": ("T", 1.0e-4, "mks"),
-    "G": ("T", 1.0e-4, "mks"),
-    "Mx": ("Wb", 1.0e-8, "mks"),
-    "statA": ("A", 10.0/speed_of_light_cm_per_s, "mks"),
-    "statV": ("V", 1.0e8/speed_of_light_cm_per_s, "mks"),
-    "statohm": ("ohm", speed_of_light_cm_per_s**2*1.0e-9, "mks"),
-}
-
-
-def _get_em_base_unit(units):
-    unit_str = str(units)
-    if len(unit_str) == 1:
-        return unit_str
-    possible_prefix = unit_str[0]
-    prefix_len = 1
-    if unit_str[:2] == 'da':
-        possible_prefix = 'da'
-        prefix_len += 1
-    if possible_prefix in unit_prefixes:
-        base_unit = unit_str[prefix_len:]
-    else:
-        base_unit = unit_str
-    return base_unit
-
-
-def _check_em_conversion(units):
-    base_unit = _get_em_base_unit(units)
-    em_info = em_conversions.get(base_unit, (None,)*3)
-    return em_info[0], em_info[2]
-
-
-class ElectromagneticSI(Equivalence):
-    """An equivalence between CGS and SI electromagnetic units
-
-    Given data in CGS electromagnetic units (one of esu, gauss, statA, statV,
-    or statohm) this equivelency will convert the data to the appropriate
-    SI electromagnetic unit, using the following mapping:
-
-    * Fr -> C
-    * statC -> C
-    * esu -> C
-    * G -> T
-    * statA -> A
-    * statV -> V
-    * statohm -> ohm
-    * Mx -> Wb
-
-    Example
-    -------
-    >>> from unyt import gauss
-    >>> (10*gauss).to_equivalent('T', 'SI')
-    unyt_quantity(0.001, 'T')
-    """
-    type_name = "SI"
-    alternate_names = ["si", "MKS", "mks"]
-    one_way = True
-    _dims = (current_cgs, charge_cgs, magnetic_field_cgs,
-             electric_field_cgs, electric_potential_cgs,
-             resistance_cgs, magnetic_flux_cgs)
-
-    def _convert(self, x, new_dims):
-        base_unit = _get_em_base_unit(x.units)
-        new_units, convert_factor, _ = em_conversions[base_unit]
-        return x.in_cgs().v*convert_factor, new_units
-
-    def __str__(self):
-        return "SI: EM CGS unit -> EM SI unit"
-
-
-class ElectromagneticCGS(Equivalence):
-    """An equivalence between SI and CGS electromagnetic units
-
-    Given data in SI electromagnetic units (one of C, T, A, V, ohm, or Wb,
-    or a prefixed version of these), this equivalency will convert the
-    data to the appropriate CGS electromagnetic unit, using the following
-    mapping:
-
-    * C -> esu
-    * T -> G
-    * A -> statA
-    * V -> statV
-    * ohm -> statohm
-    * Wb -> Mx
-
-    Example
-    -------
-    >>> from unyt import Tesla
-    >>> (10*Tesla).to_equivalent('G', 'CGS')
-    unyt_quantity(100000., 'G')
-    """
-    type_name = "CGS"
-    alternate_names = ["cgs"]
-    one_way = True
-    _dims = (current_mks, charge_mks, magnetic_field_mks,
-             electric_field_mks, electric_potential_mks,
-             resistance_mks, magnetic_flux_mks)
-
-    def _convert(self, x, new_dims):
-        base_unit = _get_em_base_unit(x.units)
-        new_units, convert_factor, _ = em_conversions[base_unit]
-        return x.in_mks().v*convert_factor, new_units
-
-    def __str__(self):
-        return "CGS: EM SI unit -> EM CGS unit"
diff --git a/unyt/exceptions.py b/unyt/exceptions.py
index e7c50a8..6574ab2 100644
--- a/unyt/exceptions.py
+++ b/unyt/exceptions.py
@@ -41,8 +41,7 @@ class UnitConversionError(Exception):
         Exception.__init__(self)
 
     def __str__(self):
-        err = ("Unit dimensionalities do not match. Tried to convert between "
-               "%s (dim %s) and %s (dim %s)." %
+        err = ("Cannot convert between %s (dim %s) and %s (dim %s)." %
                (self.unit1, self.dimension1, self.unit2, self.dimension2))
         return err
 
@@ -57,6 +56,15 @@ class MissingMKSCurrent(Exception):
         return err
 
 
+class MKSCGSConversionError(Exception):
+    def __init__(self, unit):
+        self.unit = unit
+
+    def __str__(self):
+        err = ("The %s unit cannot be safely converted." % self.unit)
+        return err
+
+
 class UnitsNotReducible(Exception):
     def __init__(self, unit, units_base):
         self.unit = unit
@@ -100,8 +108,14 @@ class InvalidUnitEquivalence(Exception):
         self.unit2 = unit2
 
     def __str__(self):
-        return ("The unit equivalence '%s' does not exist for the units '%s' "
-                "and '%s'." % (self.equiv, self.unit1, self.unit2))
+        from unyt.unit_object import Unit
+        if isinstance(self.unit2, Unit):
+            msg = ("The unit equivalence '%s' does not exist for the units "
+                   "'%s' and '%s'.")
+        else:
+            msg = ("The unit equivalence '%s' does not exist for units '%s' "
+                   "to convert to a new unit with dimensions '%s'.")
+        return msg % (self.equiv, self.unit1, self.unit2)
 
 
 class InvalidUnitOperation(Exception):
diff --git a/unyt/unit_object.py b/unyt/unit_object.py
index dd3a28c..94c54ef 100644
--- a/unyt/unit_object.py
+++ b/unyt/unit_object.py
@@ -13,7 +13,13 @@ A class that represents a unit symbol.
 # -----------------------------------------------------------------------------
 
 
+import copy
+from keyword import iskeyword as _iskeyword
+import numpy as np
+from numbers import Number as numeric_type
 from six import text_type
+import token
+
 from sympy import (
     Expr,
     Mul,
@@ -37,13 +43,21 @@ from sympy.parsing.sympy_parser import (
     auto_number,
     rationalize
 )
-from keyword import iskeyword as _iskeyword
+
 from unyt.dimensions import (
+    angle,
     base_dimensions,
-    temperature,
     dimensionless,
-    angle
+    temperature,
+)
+from unyt.equivalencies import equivalence_registry
+from unyt.exceptions import (
+    InvalidUnitOperation,
+    MKSCGSConversionError,
+    UnitConversionError,
+    UnitsNotReducible
 )
+from unyt._physical_ratios import speed_of_light_cm_per_s
 from unyt._unit_lookup_table import (
     unit_prefixes,
     prefixable_units,
@@ -54,19 +68,6 @@ from unyt.unit_registry import (
     UnitRegistry,
     UnitParseError
 )
-from unyt.exceptions import (
-    InvalidUnitOperation,
-    MissingMKSCurrent,
-    UnitsNotReducible,
-)
-from unyt.equivalencies import (
-    equivalence_registry,
-    _check_em_conversion,
-)
-import copy
-import token
-import numpy as np
-from numbers import Number as numeric_type
 
 sympy_one = sympify(1)
 
@@ -152,14 +153,15 @@ def _get_latex_representation(expr, registry):
         for i in range(1, len(symbols)):
             expr = expr.subs(symbols[i], symbols[0])
     prefix = None
+    l_expr = expr
     if isinstance(expr, Mul):
         coeffs = expr.as_coeff_Mul()
         if coeffs[0] == 1 or not isinstance(coeffs[0], Float):
-            pass
+            l_expr = coeffs[1]
         else:
-            expr = coeffs[1]
+            l_expr = coeffs[1]
             prefix = Float(coeffs[0], 2)
-    latex_repr = latex(expr, symbol_names=symbol_table, mul_symbol="dot",
+    latex_repr = latex(l_expr, symbol_names=symbol_table, mul_symbol="dot",
                        fold_frac_powers=True, fold_short_frac=True)
 
     if prefix is not None:
@@ -174,7 +176,7 @@ def _get_latex_representation(expr, registry):
 unit_text_transform = (_auto_positive_symbol, rationalize, auto_number)
 
 
-class Unit(Expr):
+class Unit(object):
     """
     A symbolic unit, using sympy functionality. We only add "dimensions" so
     that sympy understands relations between different units.
@@ -193,8 +195,7 @@ class Unit(Expr):
     __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):
+                dimensions=None, registry=None, latex_repr=None):
         """
         Create a new unit. May be an atomic unit (like a gram) or combinations
         of atomic units (like g / cm**3).
@@ -217,31 +218,23 @@ class Unit(Expr):
         latex_repr : string
             A string to render the unit as LaTeX
 
-        Additional keyword arguments are passed as assumptions to the Sympy
-        Expr initializer
-
         """
         # Simplest case. If user passes a Unit object, just use the expr.
-        unit_key = None
         if isinstance(unit_expr, (str, bytes, text_type)):
             if isinstance(unit_expr, bytes):
                 unit_expr = unit_expr.decode("utf-8")
 
-            if registry and unit_expr in registry.unit_objs:
-                return registry.unit_objs[unit_expr]
-            else:
-                unit_key = unit_expr
-                if not unit_expr:
-                    # Bug catch...
-                    # if unit_expr is an empty string, parse_expr fails hard...
-                    unit_expr = "1"
-                try:
-                    unit_expr = parse_expr(unit_expr, global_dict=global_dict,
-                                           transformations=unit_text_transform)
-                except SyntaxError as e:
-                    msg = ("Unit expression %s raised an error "
-                           "during parsing:\n%s" % (unit_expr, repr(e)))
-                    raise UnitParseError(msg)
+            if not unit_expr:
+                # Bug catch...
+                # if unit_expr is an empty string, parse_expr fails hard...
+                unit_expr = "1"
+            try:
+                unit_expr = parse_expr(unit_expr, global_dict=global_dict,
+                                       transformations=unit_text_transform)
+            except SyntaxError as e:
+                msg = ("Unit expression %s raised an error "
+                       "during parsing:\n%s" % (unit_expr, repr(e)))
+                raise UnitParseError(msg)
         elif isinstance(unit_expr, Unit):
             # grab the unit object's sympy expression.
             unit_expr = unit_expr.expr
@@ -307,7 +300,7 @@ class Unit(Expr):
                 base_offset = 0.0
 
         # Create obj with superclass construct.
-        obj = Expr.__new__(cls, **assumptions)
+        obj = super(Unit, cls).__new__(cls)
 
         # Attach attributes to obj.
         obj.expr = unit_expr
@@ -318,9 +311,6 @@ class Unit(Expr):
         obj._latex_repr = latex_repr
         obj.registry = registry
 
-        if unit_key is not None:
-            registry.unit_objs[unit_key] = obj
-
         # Return `obj` so __init__ can handle it.
 
         return obj
@@ -348,18 +338,9 @@ class Unit(Expr):
     def units(self):
         return self
 
-    # Some sympy conventions
-    def __getnewargs__(self):
-        return (self.expr, self.is_atomic, self.base_value, self.dimensions,
-                self.registry)
-
     def __hash__(self):
         return super(Unit, self).__hash__()
 
-    def _hashable_content(self):
-        return (self.expr, self.is_atomic, self.base_value, self.dimensions,
-                self.registry)
-
     # end sympy conventions
 
     def __repr__(self):
@@ -374,10 +355,6 @@ class Unit(Expr):
         # @todo: don't use dunder method?
         return self.expr.__str__()
 
-    # for sympy.printing
-    def _sympystr(self, *args):
-        return str(self.expr)
-
     #
     # Start unit operations
     #
@@ -476,15 +453,12 @@ class Unit(Expr):
 
         base_offset = 0.0
         if self.base_offset or u.base_offset:
-            if u.dimensions in (temperature, angle) and self.is_dimensionless:
-                base_offset = u.base_offset
-            elif (self.dimensions in (temperature, angle) and
-                  u.is_dimensionless):
+            if self.dimensions in (temperature, angle) and u.is_dimensionless:
                 base_offset = self.base_offset
             else:
                 raise InvalidUnitOperation(
                     "Quantities with units of Farhenheit "
-                    "and Celsius cannot be multiplied.")
+                    "and Celsius cannot be divided.")
 
         return Unit(self.expr / u.expr,
                     base_value=(self.base_value / u.base_value),
@@ -536,8 +510,6 @@ class Unit(Expr):
         return copy.deepcopy(self)
 
     def __deepcopy__(self, memodict=None):
-        if memodict is None:
-            memodict = {}
         expr = str(self.expr)
         base_value = copy.deepcopy(self.base_value)
         base_offset = copy.deepcopy(self.base_offset)
@@ -590,13 +562,12 @@ class Unit(Expr):
 
         Returns
         -------
-        True if the unit name begins with "code" False otherwise
+        True if the unit consists of atom units that being with "code".
+        False otherwise
 
         """
         for atom in self.expr.atoms():
-            if str(atom).startswith("code") or atom.is_Number:
-                pass
-            else:
+            if not (str(atom).startswith("code") or atom.is_Number):
                 return False
         return True
 
@@ -623,8 +594,10 @@ class Unit(Expr):
         Example
         -------
         >>> from unyt import km
-        >>> km.units.has_equivalent('spectral')
+        >>> km.has_equivalent('spectral')
         True
+        >>> km.has_equivalent('mass_energy')
+        False
         """
         try:
             this_equiv = equivalence_registry[equiv]()
@@ -637,9 +610,9 @@ class Unit(Expr):
         """Create and return dimensionally-equivalent units in a specified base.
 
         >>> from unyt import g, cm
-        >>> (g/cm**3).units.get_base_equivalent('mks')
+        >>> (g/cm**3).get_base_equivalent('mks')
         kg/m**3
-        >>> (g/cm**3).units.get_base_equivalent('solar')
+        >>> (g/cm**3).get_base_equivalent('solar')
         Mearth/AU**3
         """
         from unyt.unit_systems import unit_system_registry
@@ -649,14 +622,15 @@ class Unit(Expr):
             unit_system = self.registry.unit_system_id
         unit_system = unit_system_registry[str(unit_system)]
         try:
-            em_units, em_us = _check_em_conversion(self.units)
-            if em_units is not None and em_us == str(unit_system):
-                units_string = em_units
-            else:
-                units_string = unit_system[self.dimensions]
-        except MissingMKSCurrent:
+            conv_data = _check_em_conversion(self.units)
+        except MKSCGSConversionError:
             raise UnitsNotReducible(self.units, unit_system)
-        return Unit(units_string, registry=self.registry)
+        if any(conv_data):
+            new_units, _ = _em_conversion(
+                self, conv_data, unit_system=unit_system)
+        else:
+            new_units = unit_system[self.dimensions]
+        return Unit(new_units, registry=self.registry)
 
     def get_cgs_equivalent(self):
         """Create and return dimensionally-equivalent cgs units.
@@ -664,7 +638,7 @@ class Unit(Expr):
         Example
         -------
         >>> from unyt import kg, m
-        >>> (kg/m**3).units.get_cgs_equivalent()
+        >>> (kg/m**3).get_cgs_equivalent()
         g/cm**3
         """
         return self.get_base_equivalent(unit_system="cgs")
@@ -675,7 +649,7 @@ class Unit(Expr):
         Example
         -------
         >>> from unyt import g, cm
-        >>> (g/cm**3).units.get_mks_equivalent()
+        >>> (g/cm**3).get_mks_equivalent()
         kg/m**3
         """
         return self.get_base_equivalent(unit_system="mks")
@@ -699,9 +673,9 @@ class Unit(Expr):
         Examples
         --------
         >>> from unyt import km, cm, degree_fahrenheit, degree_celsius
-        >>> km.units.get_conversion_factor(cm.units)
+        >>> km.get_conversion_factor(cm)
         (100000.0, None)
-        >>> degree_celsius.units.get_conversion_factor(degree_fahrenheit.units)
+        >>> degree_celsius.get_conversion_factor(degree_fahrenheit)
         (1.7999999999999998, -31.999999999999886)
         """
         return _get_conversion_factor(self, other_units)
@@ -712,7 +686,7 @@ class Unit(Expr):
         Examples
         --------
         >>> from unyt import g, cm
-        >>> (g/cm**3).units.latex_representation()
+        >>> (g/cm**3).latex_representation()
         '\\\\frac{\\\\rm{g}}{\\\\rm{cm}^{3}}'
         """
         return self.latex_repr
@@ -722,6 +696,98 @@ class Unit(Expr):
 #
 
 
+em_conversions = {
+    "C": ("esu", 0.1*speed_of_light_cm_per_s, "cgs"),
+    "T": ("gauss", 1.0e4, "cgs"),
+    "Wb": ("Mx", 1.0e8, "cgs"),
+    "A": ("statA", 0.1*speed_of_light_cm_per_s, "cgs"),
+    "V": ("statV", 1.0e-8*speed_of_light_cm_per_s, "cgs"),
+    "ohm": ("statohm", 1.0e9/(speed_of_light_cm_per_s**2), "cgs"),
+    "esu": ("C", 10.0/speed_of_light_cm_per_s, "mks"),
+    "Fr": ("C", 10.0/speed_of_light_cm_per_s, "mks"),
+    "statC": ("C", 10.0/speed_of_light_cm_per_s, "mks"),
+    "gauss": ("T", 1.0e-4, "mks"),
+    "G": ("T", 1.0e-4, "mks"),
+    "Mx": ("Wb", 1.0e-8, "mks"),
+    "statA": ("A", 10.0/speed_of_light_cm_per_s, "mks"),
+    "statV": ("V", 1.0e8/speed_of_light_cm_per_s, "mks"),
+    "statohm": ("ohm", speed_of_light_cm_per_s**2*1.0e-9, "mks"),
+}
+
+
+def _em_conversion(orig_units, conv_data, to_units=None, unit_system=None):
+    """Convert between E&M & MKS base units.
+
+    If orig_units is a CGS (or MKS) E&M unit, conv_data contains the
+    corresponding MKS (or CGS) unit and scale factor converting between them.
+    This must be done by replacing the expression of the original unit
+    with the new one in the unit expression and multiplying by the scale
+    factor.
+    """
+    conv_unit, scale = conv_data
+    new_expr = orig_units.copy().expr.replace(
+        orig_units.expr, scale*conv_unit.expr)
+    if unit_system is not None:
+        # we don't know the to_units, so we get it directly from the
+        # conv_data
+        inter_expr = orig_units.copy().expr.replace(
+            orig_units.expr, conv_unit.expr)
+        to_units = Unit(inter_expr, registry=orig_units.registry)
+    new_units = Unit(new_expr, registry=orig_units.registry)
+    conv = new_units.get_conversion_factor(to_units)
+    return to_units, conv
+
+
+def _get_em_base_unit(units):
+    unit_str = str(units)
+    if len(unit_str) == 1:
+        return unit_str
+    possible_prefix = unit_str[0]
+    prefix_len = 1
+    if unit_str[:2] == 'da':
+        possible_prefix = 'da'
+        prefix_len += 1
+    if possible_prefix in unit_prefixes:
+        base_unit = unit_str[prefix_len:]
+    else:
+        base_unit = unit_str
+    return base_unit
+
+
+def _check_em_conversion(units, to_units=None):
+    """Check to see if the units contain E&M units
+
+    This function supports unyt's ability to convert data to and from E&M
+    electromagnetic units. However, this support is limited and only very
+    simple unit expressions can be readily converted. This function
+    to see if the unit is an atomic base unit that is present in the
+    em_conversions dict. If it does not contain E&M units, the function
+    returns an empty tuple. If it does contain an atomic E&M unit in
+    the em_conversions dict, it returns a tuple containing the unit to convert
+    to and scale factor. If it contains a more complicated E&M unit and we are
+    trying to convert between CGS & MKS E&M units, it raises an error.
+    """
+    em_map = ()
+    base_unit = _get_em_base_unit(str(units))
+    base_to_unit = _get_em_base_unit(str(to_units))
+    if base_unit == base_to_unit:
+        return em_map
+    if base_unit in em_conversions:
+        em_info = em_conversions.get(base_unit, (None,)*2)
+        to_unit = Unit(em_info[0], registry=units.registry)
+        if to_units is None or to_units.dimensions == to_unit.dimensions:
+            em_map = (to_unit, em_info[1])
+    if em_map:
+        return em_map
+    for unit in units.expr.atoms():
+        if unit.is_Number:
+            pass
+        base_unit = _get_em_base_unit(str(unit))
+        if base_unit in em_conversions:
+            raise MKSCGSConversionError(units)
+    return em_map
+
+
 def _get_conversion_factor(old_units, new_units):
     """
     Get the conversion factor between two units of equivalent dimensions. This
@@ -743,16 +809,16 @@ def _get_conversion_factor(old_units, new_units):
         Offset between the old unit and new unit.
 
     """
+    if old_units.dimensions != new_units.dimensions:
+        raise UnitConversionError(old_units, old_units.dimensions,
+                                  new_units, new_units.dimensions)
     ratio = old_units.base_value / new_units.base_value
     if old_units.base_offset == 0 and new_units.base_offset == 0:
         return (ratio, None)
     else:
-        if old_units.dimensions in (temperature, angle):
-            return ratio, ratio*old_units.base_offset - new_units.base_offset
-        else:
-            raise InvalidUnitOperation(
-                "Fahrenheit and Celsius are not absolute temperature scales "
-                "and cannot be used in compound unit symbols.")
+        # the dimensions are the same, so both are temperatures, where
+        # it's legal to convert units so no need to do error checking
+        return ratio, ratio*old_units.base_offset - new_units.base_offset
 
 #
 # Helper functions
@@ -771,10 +837,6 @@ def _get_unit_data_from_expr(unit_expr, unit_symbol_lut):
         Provides the unit data for each valid unit symbol.
 
     """
-    # The simplest case first
-    if isinstance(unit_expr, Unit):
-        return (unit_expr.base_value, unit_expr.dimensions)
-
     # Now for the sympy possibilities
     if isinstance(unit_expr, Number):
         if unit_expr is sympy_one:
diff --git a/unyt/unit_registry.py b/unyt/unit_registry.py
index e1d19db..b231d9c 100644
--- a/unyt/unit_registry.py
+++ b/unyt/unit_registry.py
@@ -14,9 +14,7 @@ A registry for units that can be added to and modified.
 
 
 import json
-import re
 
-from distutils.version import LooseVersion
 from unyt.exceptions import (
     SymbolNotFoundError,
     UnitParseError,
@@ -27,15 +25,8 @@ import six
 from sympy import (
     sympify,
     srepr,
-    __version__ as sympy_version
 )
 
-SYMPY_VERSION = LooseVersion(sympy_version)
-
-
-def _positive_symbol_replacer(match):
-    return match.group().replace(')\')', ')\', positive=True)')
-
 
 class UnitRegistry:
     """A registry for unit symbols"""
@@ -46,7 +37,6 @@ class UnitRegistry:
             self.lut = lut
         else:
             self.lut = {}
-        self.unit_objs = {}
 
         if add_default_symbols:
             self.lut.update(default_unit_symbol_lut)
@@ -58,8 +48,6 @@ class UnitRegistry:
         return item in self.lut
 
     def pop(self, item):
-        if item in self.unit_objs:
-            del self.unit_objs[item]
         if item in self.lut:
             del self.lut[item]
 
@@ -124,8 +112,6 @@ class UnitRegistry:
                 "in this registry." % symbol)
 
         del self.lut[symbol]
-        if symbol in self.unit_objs:
-            del self.unit_objs[symbol]
 
     def modify(self, symbol, base_value):
         """
@@ -150,9 +136,6 @@ class UnitRegistry:
         self.lut[symbol] = ((float(base_value), new_dimensions) +
                             self.lut[symbol][2:])
 
-        if symbol in self.unit_objs:
-            del self.unit_objs[symbol]
-
     def keys(self):
         """
         Print out the units contained in the lookup table.
@@ -168,10 +151,6 @@ class UnitRegistry:
         for k, v in six.iteritems(self.lut):
             san_v = list(v)
             repr_dims = srepr(v[1])
-            if SYMPY_VERSION < LooseVersion("1.0.0"):
-                # see https://github.com/sympy/sympy/issues/6131
-                repr_dims = re.sub("Symbol\('\([a-z_]*\)'\)",
-                                   _positive_symbol_replacer, repr_dims)
             san_v[1] = repr_dims
             sanitized_lut[k] = tuple(san_v)
 
@@ -198,8 +177,6 @@ class UnitRegistry:
         """
         equiv = [k for k, v in self.lut.items()
                  if v[1] is unit_object.dimensions]
-        equiv += [n for n, u in self.unit_objs.items()
-                  if u.dimensions is unit_object.dimensions]
         equiv = list(sorted(set(equiv)))
         return equiv
 
diff --git a/unyt/unit_symbols.py b/unyt/unit_symbols.py
index 7a96bc4..fc2b024 100644
--- a/unyt/unit_symbols.py
+++ b/unyt/unit_symbols.py
@@ -311,6 +311,8 @@ HA = hourangle = Unit("hourangle")
 electrostatic_unit = esu = Unit("esu")
 #: Gauss
 gauss = G = Unit("gauss")
+#: Statcoulomb
+statcoulomb = statC = Unit("statC")
 #: Statampere
 statampere = statA = Unit("statA")
 #: Statvolt
diff --git a/unyt/unit_systems.py b/unyt/unit_systems.py
index 9da9a37..a0f4a53 100644
--- a/unyt/unit_systems.py
+++ b/unyt/unit_systems.py
@@ -22,31 +22,12 @@ from unyt.unit_object import (
     Unit,
     _get_system_unit_string
 )
-from unyt import physical_constants as pc
 
 unit_system_registry = {}
 
 cmks = dimensions.current_mks
 
 
-class _UnitSystemConstants(object):
-    """
-    A class to facilitate conversions of physical constants into a given unit
-    system specified by *name*.
-    """
-    def __init__(self, name):
-        self.name = name
-
-    def __repr__(self):
-        return "Physical constants in %s units." % self.name
-
-    def __str__(self):
-        return self.name
-
-    def __getattr__(self, item):
-        return getattr(pc, item).in_base(self.name)
-
-
 class UnitSystem(object):
     """
     Create a UnitSystem for facilitating conversions to a default set of units.
@@ -108,13 +89,12 @@ class UnitSystem(object):
         self.base_units = self.units_map.copy()
         unit_system_registry[name] = self
         self.name = name
-        self.constants = _UnitSystemConstants(self.name)
 
     def __getitem__(self, key):
         if isinstance(key, string_types):
             key = getattr(dimensions, key)
         um = self.units_map
-        if key not in um or um[key].dimensions is not key:
+        if key not in um or um[key] is None or um[key].dimensions is not key:
             if cmks in key.free_symbols and self.units_map[cmks] is None:
                 raise MissingMKSCurrent(self.name)
             units = _get_system_unit_string(key, self.units_map)


Test patch:
diff --git a/unyt/tests/test_unit_registry.py b/unyt/tests/test_unit_registry.py
new file mode 100644
index 0000000..c3d4c60
--- /dev/null
+++ b/unyt/tests/test_unit_registry.py
@@ -0,0 +1,54 @@
+"""
+Test unit lookup tables and registry
+
+
+
+
+"""
+
+# -----------------------------------------------------------------------------
+# Copyright (c) 2018, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the LICENSE file, distributed with this software.
+# -----------------------------------------------------------------------------
+
+import pytest
+
+from unyt.dimensions import length
+from unyt.exceptions import (
+    SymbolNotFoundError,
+    UnitParseError,
+)
+from unyt.unit_registry import UnitRegistry
+
+
+def test_add_modify_error():
+    from unyt import m
+    ureg = UnitRegistry()
+
+    with pytest.raises(UnitParseError):
+        ureg.add('tayne', 1, length)
+    with pytest.raises(UnitParseError):
+        ureg.add('tayne', 'blah', length)
+    with pytest.raises(UnitParseError):
+        ureg.add('tayne', 1.0, length, offset=1)
+    with pytest.raises(UnitParseError):
+        ureg.add('tayne', 1.0, length, offset='blah')
+
+    ureg.add('tayne', 1.1, length)
+
+    with pytest.raises(SymbolNotFoundError):
+        ureg.remove('tayn')
+    with pytest.raises(SymbolNotFoundError):
+        ureg.modify('tayn', 1.2)
+
+    ureg.modify('tayne', 1.0*m)
+
+    assert ureg['tayne'][:3] == ureg['m'][:3]
+
+
+def test_keys():
+    ureg = UnitRegistry()
+    assert sorted(ureg.keys()) == sorted(ureg.lut.keys())
diff --git a/unyt/tests/test_unit_systems.py b/unyt/tests/test_unit_systems.py
index ddf4d01..5f2e30c 100644
--- a/unyt/tests/test_unit_systems.py
+++ b/unyt/tests/test_unit_systems.py
@@ -13,18 +13,18 @@ Test unit systems.
 
 import pytest
 
-from unyt.exceptions import IllDefinedUnitSystem
+from unyt.exceptions import (
+    IllDefinedUnitSystem,
+    MissingMKSCurrent,
+)
 from unyt.unit_object import Unit
 from unyt.unit_systems import (
     UnitSystem,
+    cgs_unit_system,
     unit_system_registry,
 )
 from unyt.unit_registry import UnitRegistry
 from unyt import dimensions
-from unyt import (
-    gauss,
-    Tesla,
-)
 
 
 def test_unit_systems():
@@ -57,17 +57,17 @@ def test_unit_system_id():
     assert reg1.unit_system_id == reg2.unit_system_id
 
 
-def test_cgs_mks_unit_conversions():
-    t = 1*Tesla
-    g = 1*gauss
-    assert t.to_equivalent('G', 'cgs') == 1e4*gauss
-    assert t.to_equivalent('G', 'CGS') == 1e4*gauss
-    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():
     with pytest.raises(IllDefinedUnitSystem):
         UnitSystem('atomic', 'nm', 'fs', 'nK', 'rad')
+
+
+def test_mks_current():
+    with pytest.raises(MissingMKSCurrent):
+        cgs_unit_system[dimensions.current_mks]
+    with pytest.raises(MissingMKSCurrent):
+        cgs_unit_system[dimensions.magnetic_field]
+    with pytest.raises(MissingMKSCurrent):
+        cgs_unit_system[dimensions.current_mks] = 'foo'
+    with pytest.raises(MissingMKSCurrent):
+        cgs_unit_system[dimensions.magnetic_field] = 'bar'
diff --git a/unyt/tests/test_units.py b/unyt/tests/test_units.py
index 5bef117..eea75cc 100644
--- a/unyt/tests/test_units.py
+++ b/unyt/tests/test_units.py
@@ -24,6 +24,7 @@ from numpy.testing import (
 )
 import operator
 import pytest
+from six.moves import cPickle as pickle
 from sympy import Symbol
 from unyt._testing import assert_allclose_units
 from unyt.unit_registry import UnitRegistry
@@ -38,13 +39,17 @@ from unyt.dimensions import (
     power,
     rate
 )
+from unyt.exceptions import (
+    InvalidUnitOperation,
+    UnitsNotReducible,
+    UnitConversionError,
+)
 from unyt.unit_object import (
     default_unit_registry,
-    _get_conversion_factor,
     Unit,
     UnitParseError,
-    InvalidUnitOperation
 )
+from unyt.unit_systems import UnitSystem
 from unyt._unit_lookup_table import (
     default_unit_symbol_lut,
     unit_prefixes,
@@ -77,9 +82,8 @@ def test_no_conflicting_symbols():
             new_symbol = "%s%s" % (prefix, symbol)
 
             # test if we have seen this symbol
-            if new_symbol in full_set:
-                print("Duplicate symbol: %s" % new_symbol)
-                raise RuntimeError
+            assert new_symbol not in full_set, ("Duplicate symbol: %s" %
+                                                new_symbol)
 
             full_set.add(new_symbol)
 
@@ -95,6 +99,8 @@ def test_dimensionless():
     assert u1.expr == 1
     assert u1.base_value == 1
     assert u1.dimensions == 1
+    assert u1 != 'hello!'
+    assert (u1 == 'hello') is False
 
     u2 = Unit("")
 
@@ -151,6 +157,13 @@ def test_create_from_string():
         Unit('m+g')
     with pytest.raises(UnitParseError):
         Unit('m-g')
+    with pytest.raises(UnitParseError):
+        Unit('hello!')
+
+    cm = Unit('cm')
+    data = 1*cm
+
+    assert Unit(data) == cm
 
 
 def test_create_from_expr():
@@ -246,12 +259,8 @@ def test_create_fail_on_unknown_symbol():
     Fail to create unit with unknown symbol, without base_value and dimensions.
 
     """
-    try:
+    with pytest.raises(UnitParseError):
         Unit(Symbol("jigawatts"))
-    except UnitParseError:
-        assert True
-    else:
-        assert False
 
 
 def test_create_fail_on_bad_symbol_type():
@@ -259,12 +268,8 @@ def test_create_fail_on_bad_symbol_type():
     Fail to create unit with bad symbol type.
 
     """
-    try:
+    with pytest.raises(UnitParseError):
         Unit([1])  # something other than Expr and str
-    except UnitParseError:
-        assert True
-    else:
-        assert False
 
 
 def test_create_fail_on_bad_dimensions_type():
@@ -272,12 +277,8 @@ def test_create_fail_on_bad_dimensions_type():
     Fail to create unit with bad dimensions type.
 
     """
-    try:
+    with pytest.raises(UnitParseError):
         Unit("a", base_value=1, dimensions="(mass)")
-    except UnitParseError:
-        assert True
-    else:
-        assert False
 
 
 def test_create_fail_on_dimensions_content():
@@ -286,13 +287,8 @@ def test_create_fail_on_dimensions_content():
 
     """
     a = Symbol("a")
-
-    try:
+    with pytest.raises(UnitParseError):
         Unit("a", base_value=1, dimensions=a)
-    except UnitParseError:
-        pass
-    else:
-        assert False
 
 
 def test_create_fail_on_base_value_type():
@@ -300,12 +296,8 @@ def test_create_fail_on_base_value_type():
     Fail to create unit with bad base_value type.
 
     """
-    try:
+    with pytest.raises(UnitParseError):
         Unit("a", base_value="a", dimensions=(mass/time))
-    except UnitParseError:
-        assert True
-    else:
-        assert False
 
 
 def test_string_representation():
@@ -323,6 +315,7 @@ def test_string_representation():
     assert str(speed) == "pc/Myr"
     assert repr(speed) == "pc/Myr"
     assert str(dimensionless) == "dimensionless"
+    assert repr(dimensionless) == "(dimensionless)"
 
 
 def test_multiplication():
@@ -420,6 +413,36 @@ def test_equality():
     assert u1 == u2
 
 
+def test_invalid_operations():
+    u1 = Unit('cm')
+    u2 = Unit('m')
+
+    with pytest.raises(InvalidUnitOperation):
+        u1+u2
+    with pytest.raises(InvalidUnitOperation):
+        u1 += u2
+    with pytest.raises(InvalidUnitOperation):
+        1+u1
+    with pytest.raises(InvalidUnitOperation):
+        u1+1
+    with pytest.raises(InvalidUnitOperation):
+        u1-u2
+    with pytest.raises(InvalidUnitOperation):
+        u1 -= u2
+    with pytest.raises(InvalidUnitOperation):
+        1-u1
+    with pytest.raises(InvalidUnitOperation):
+        u1-1
+    with pytest.raises(InvalidUnitOperation):
+        u1 *= u2
+    with pytest.raises(InvalidUnitOperation):
+        u1*'hello!'
+    with pytest.raises(InvalidUnitOperation):
+        u1 /= u2
+    with pytest.raises(InvalidUnitOperation):
+        u1/'hello!'
+
+
 def test_base_equivalent():
     """
     Check base equivalent of a unit.
@@ -445,9 +468,15 @@ def test_base_equivalent():
     assert u2.dimensions == mass_density
     assert u3.dimensions == mass_density
 
-    assert_allclose_units(_get_conversion_factor(u1, u3)[0],
+    assert_allclose_units(u1.get_conversion_factor(u3)[0],
                           Msun_mks / Mpc_mks**3, 1e-12)
 
+    with pytest.raises(UnitConversionError):
+        u1.get_conversion_factor(Unit('m'))
+
+    with pytest.raises(UnitConversionError):
+        u1.get_conversion_factor(Unit('degF'))
+
 
 def test_temperature_offsets():
     u1 = Unit('degC')
@@ -479,6 +508,9 @@ def test_latex_repr():
     test_unit = Unit('1e9*cm')
     assert_equal(test_unit.latex_repr, '1.0 \\times 10^{9}\\ \\rm{cm}')
 
+    test_unit = Unit('1.0*cm')
+    assert_equal(test_unit.latex_repr, '\\rm{cm}')
+
 
 def test_latitude_longitude():
     lat = unit_symbols.lat
@@ -529,19 +561,98 @@ def test_creation_from_ytarray():
 
 
 def test_list_same_dimensions():
+    from unyt import m
     reg = default_unit_registry
-    for name1, u1 in reg.unit_objs.items():
-        for name2 in reg.list_same_dimensions(u1):
-            if name2 == name1:
-                continue
-            if name2 in reg.unit_objs:
-                dim2 = reg.unit_objs[name2].dimensions
-            else:
-                _, dim2, _, _ = reg.lut[name2]
-            assert u1.dimensions is dim2
+    for equiv in reg.list_same_dimensions(m):
+        assert Unit(equiv).dimensions is length
 
 
 def test_decagram():
     dag = Unit('dag')
     g = Unit('g')
     assert dag.get_conversion_factor(g) == (10.0, None)
+
+
+def test_pickle():
+    cm = Unit('cm')
+    assert cm == pickle.loads(pickle.dumps(cm))
+
+
+def test_preserve_offset():
+    from unyt import degF, dimensionless
+
+    new_unit = degF*dimensionless
+
+    assert new_unit is not degF
+    assert new_unit == degF
+    assert new_unit.base_offset == degF.base_offset
+
+    new_unit = degF/dimensionless
+
+    assert new_unit is not degF
+    assert new_unit == degF
+    assert new_unit.base_offset == degF.base_offset
+
+    with pytest.raises(InvalidUnitOperation):
+        dimensionless/degF
+
+
+def test_code_unit():
+    from unyt import UnitRegistry
+
+    ureg = UnitRegistry()
+    ureg.add('code_length', 10., length)
+    u = Unit('code_length', registry=ureg)
+    assert u.is_code_unit is True
+    assert u.get_base_equivalent() == Unit('m')
+    u = Unit('cm')
+    assert u.is_code_unit is False
+
+    UnitSystem(ureg.unit_system_id, 'code_length', 'kg', 's', registry=ureg)
+
+    u = Unit('cm', registry=ureg)
+    ue = u.get_base_equivalent('code')
+
+    assert str(ue) == 'code_length'
+    assert ue.base_value == 10
+    assert ue.dimensions is length
+
+    class FakeDataset(object):
+        unit_registry = ureg
+
+    ds = FakeDataset()
+
+    UnitSystem(ds, 'code_length', 'kg', 's', registry=ureg)
+
+    u = Unit('cm', registry=ureg)
+    ue = u.get_base_equivalent(ds)
+
+    assert str(ue) == 'code_length'
+    assert ue.base_value == 10
+    assert ue.dimensions is length
+
+    with pytest.raises(UnitParseError):
+        Unit('code_length')
+
+
+def test_bad_equivalence():
+    from unyt import cm
+
+    with pytest.raises(KeyError):
+        cm.has_equivalent('dne')
+
+
+def test_em_unit_base_equivalent():
+    from unyt import A, cm
+
+    with pytest.raises(UnitsNotReducible):
+        (A/cm).get_base_equivalent('cgs')
+
+
+def test_define_unit_error():
+    from unyt import define_unit
+
+    with pytest.raises(RuntimeError):
+        define_unit('foobar', 'baz')
+    with pytest.raises(RuntimeError):
+        define_unit('foobar', 12)
diff --git a/unyt/tests/test_unyt_array.py b/unyt/tests/test_unyt_array.py
index 8cecb89..a7de894 100644
--- a/unyt/tests/test_unyt_array.py
+++ b/unyt/tests/test_unyt_array.py
@@ -38,15 +38,22 @@ from unyt.array import (
     unary_operators,
     binary_operators,
     uconcatenate,
+    ucross,
+    udot,
     uintersect1d,
+    unorm,
     uunion1d,
     loadtxt,
     savetxt
 )
 from unyt.exceptions import (
+    InvalidUnitEquivalence,
     InvalidUnitOperation,
+    IterableUnitCoercionError,
+    UnitConversionError,
     UnitOperationError,
     UnitParseError,
+    UnitsNotReducible,
 )
 from unyt._testing import assert_allclose_units
 from unyt.unit_symbols import (
@@ -58,7 +65,9 @@ from unyt.unit_symbols import (
 from unyt.unit_registry import UnitRegistry
 from unyt._on_demand_imports import (
     _astropy,
-    _pint
+    _h5py,
+    _pint,
+    NotAModule,
 )
 from unyt._physical_ratios import (
     metallicity_sun,
@@ -684,11 +693,41 @@ def test_temperature_conversions():
     with pytest.raises(InvalidUnitOperation):
         np.multiply(balmy, balmy)
     with pytest.raises(InvalidUnitOperation):
-        np.add(balmy_F, balmy_F)
+        np.multiply(balmy_F, balmy_F)
     with pytest.raises(InvalidUnitOperation):
-        np.subtract(balmy_C, balmy_C)
+        np.multiply(balmy_F, balmy_C)
+    with pytest.raises(InvalidUnitOperation):
+        np.divide(balmy, balmy)
+    with pytest.raises(InvalidUnitOperation):
+        np.divide(balmy_F, balmy_F)
     with pytest.raises(InvalidUnitOperation):
         np.divide(balmy_F, balmy_C)
+    with pytest.raises(InvalidUnitOperation):
+        balmy*km
+    with pytest.raises(InvalidUnitOperation):
+        balmy*balmy
+    with pytest.raises(InvalidUnitOperation):
+        balmy_F*balmy_F
+    with pytest.raises(InvalidUnitOperation):
+        balmy_F*balmy_C
+    with pytest.raises(InvalidUnitOperation):
+        2*balmy_F
+    with pytest.raises(InvalidUnitOperation):
+        balmy/balmy
+    with pytest.raises(InvalidUnitOperation):
+        balmy_F/balmy_F
+    with pytest.raises(InvalidUnitOperation):
+        balmy_F/balmy_C
+    assert_equal(np.add(balmy_F, balmy_F), unyt_quantity(80.33*2, 'degF'))
+    with pytest.raises(InvalidUnitOperation):
+        np.add(balmy_F, balmy_C)
+    with pytest.raises(InvalidUnitOperation):
+        balmy_F + balmy_C
+    assert_equal(np.subtract(balmy_C, balmy_C), unyt_quantity(0, 'degC'))
+    with pytest.raises(InvalidUnitOperation):
+        np.subtract(balmy_F, balmy_C)
+    with pytest.raises(InvalidUnitOperation):
+        balmy_F - balmy_C
 
     # Does CGS conversion from F to K work?
     assert_array_almost_equal(balmy.in_cgs(), unyt_quantity(300, 'K'))
@@ -785,6 +824,10 @@ def test_copy():
     assert_equal(copy.deepcopy(quan), quan)
     assert_array_equal(copy.deepcopy(arr), arr)
 
+    memo = {}
+    assert_equal(copy.deepcopy(quan, memo), quan)
+    assert_array_equal(copy.deepcopy(arr), arr)
+
     assert_equal(quan.copy(), quan)
     assert_array_equal(arr.copy(), arr)
 
@@ -817,11 +860,7 @@ def unary_ufunc_comparison(ufunc, a):
             'isinf', 'isnan', 'signbit', 'sign', 'rint', 'logical_not']):
         # These operations should return identical results compared to numpy.
         with np.errstate(invalid='ignore'):
-            try:
-                ret = ufunc(a, out=out)
-            except UnitOperationError:
-                assert ufunc in (np.deg2rad, np.rad2deg)
-                ret = ufunc(unyt_array(a, '1'))
+            ret = ufunc(a, out=out)
 
             assert_array_equal(ret, out)
             assert_array_equal(ret, ufunc(a_array))
@@ -871,14 +910,21 @@ def unary_ufunc_comparison(ufunc, a):
         assert_array_equal(ret2, npret2)
     elif ufunc is np.invert:
         with pytest.raises(TypeError):
-            ufunc(a)
+            ufunc(a.astype('int64'))
     elif hasattr(np, 'isnat') and ufunc is np.isnat:
         # numpy 1.13 raises ValueError, numpy 1.14 and newer raise TypeError
         with pytest.raises((TypeError, ValueError)):
             ufunc(a)
-    else:
-        # There shouldn't be any untested ufuncs.
-        assert False
+    # no untested ufuncs
+    assert ufunc in yield_np_ufuncs([
+        'isreal', 'iscomplex', 'exp', 'exp2', 'log', 'log2', 'log10', 'expm1',
+        'log1p', 'sin', 'cos', 'tan', 'arcsin', 'arccos', 'arctan', 'sinh',
+        'cosh', 'tanh', 'arccosh', 'arcsinh', 'arctanh', 'deg2rad', 'rad2deg',
+        'isfinite', 'isinf', 'isnan', 'signbit', 'sign', 'rint', 'logical_not',
+        'absolute', 'fabs', 'conjugate', 'floor', 'ceil', 'trunc', 'negative',
+        'spacing', 'positive', 'ones_like', 'square', 'sqrt', 'reciprocal',
+        'invert', 'isnat', 'modf', 'frexp',
+    ])
 
 
 def binary_ufunc_comparison(ufunc, a, b):
@@ -922,15 +968,11 @@ def binary_ufunc_comparison(ufunc, a, b):
 
 def test_ufuncs():
     for ufunc in unary_operators:
-        if ufunc is None:
-            continue
         unary_ufunc_comparison(ufunc, unyt_array([.3, .4, .5], 'cm'))
         unary_ufunc_comparison(ufunc, unyt_array([12, 23, 47], 'g'))
         unary_ufunc_comparison(ufunc, unyt_array([2, 4, -6], 'erg/m**3'))
 
     for ufunc in binary_operators:
-        if ufunc is None:
-            continue
         # arr**arr is undefined for arrays with units because
         # each element of the result would have different units.
         if ufunc is np.power:
@@ -965,6 +1007,12 @@ def test_reductions():
     assert_equal(ev_result.units, res.units)
     assert_isinstance(ev_result, unyt_array)
 
+    qv = unyt_array([1, 2, 3], 'cm').dot(unyt_array([1, 2, 3], 'cm'))
+    qa = unyt_quantity(14, 'cm**2')
+    assert qv == qa
+    assert qv.units == qa.units
+    assert_isinstance(qv, unyt_quantity)
+
     answers = {
         'prod': (unyt_quantity(720, 'cm**6'),
                  unyt_array([4, 10, 18], 'cm**2'),
@@ -1077,6 +1125,8 @@ def test_to_value():
 
 
 def test_astropy():
+    if isinstance(_astropy.__version__, NotAModule):
+        return
     ap_arr = np.arange(10)*_astropy.units.km/_astropy.units.hr
     yt_arr = unyt_array(np.arange(10), "km/hr")
     yt_arr2 = unyt_array.from_astropy(ap_arr)
@@ -1098,10 +1148,12 @@ def test_astropy():
 
 
 def test_pint():
+    if isinstance(_pint.UnitRegistry, NotAModule):
+        return
     ureg = _pint.UnitRegistry()
 
-    p_arr = np.arange(10)*ureg.km/ureg.hr
-    yt_arr = unyt_array(np.arange(10), "km/hr")
+    p_arr = np.arange(10)*ureg.km/ureg.year
+    yt_arr = unyt_array(np.arange(10), "km/yr")
     yt_arr2 = unyt_array.from_pint(p_arr)
 
     p_quan = 10.*ureg.g**0.5/(ureg.mm**3)
@@ -1162,8 +1214,14 @@ def test_subclass():
     assert_isinstance(a[:2], unyt_a_subclass)
     assert_isinstance(unyt_a_subclass(yta), unyt_a_subclass)
 
+    with pytest.raises(RuntimeError):
+        a + 'hello'
+
 
 def test_h5_io():
+    if isinstance(_h5py.__version__, NotAModule):
+        return
+
     tmpdir = tempfile.mkdtemp()
     curdir = os.getcwd()
     os.chdir(tmpdir)
@@ -1183,6 +1241,45 @@ def test_h5_io():
     assert_equal(warr.units.registry['code_length'],
                  iarr.units.registry['code_length'])
 
+    # test code to overwrite existing dataset
+
+    warr.write_hdf5('test.h5')
+
+    giarr = unyt_array.from_hdf5('test.h5')
+
+    assert_equal(warr, giarr)
+
+    # test code to overwrite existing dataset with data that has a different
+    # shape
+
+    warr = unyt_array(np.random.random((255, 255)), 'code_length',
+                      registry=reg)
+
+    warr.write_hdf5('test.h5')
+
+    giarr = unyt_array.from_hdf5('test.h5')
+
+    assert_equal(warr, giarr)
+
+    os.remove('test.h5')
+
+    # write to a group that doesn't exist
+
+    warr.write_hdf5('test.h5', dataset_name="test_dset",
+                    group_name='/arrays/test_group')
+
+    giarr = unyt_array.from_hdf5('test.h5', dataset_name="test_dset",
+                                 group_name='/arrays/test_group')
+
+    assert_equal(warr, giarr)
+
+    os.remove('test.h5')
+
+    # write to a group that does exist
+
+    with _h5py.File('test.h5') as f:
+        f.create_group('/arrays/test_group')
+
     warr.write_hdf5('test.h5', dataset_name="test_dset",
                     group_name='/arrays/test_group')
 
@@ -1191,6 +1288,8 @@ def test_h5_io():
 
     assert_equal(warr, giarr)
 
+    os.remove('test.h5')
+
     os.chdir(curdir)
     shutil.rmtree(tmpdir)
 
@@ -1198,6 +1297,28 @@ def test_h5_io():
 def test_equivalencies():
     import unyt as u
 
+    # equivalence is ignored if the conversion doesn't need one
+    data = 12*u.g
+    data.convert_to_equivalent('kg', None)
+    assert data.value == .012
+    assert data.units == u.kg
+
+    data = 12*u.g
+    data = data.to_equivalent('kg', None)
+    assert data.value == .012
+    assert data.units == u.kg
+
+    # incorrect usage of an equivalence raises errors
+
+    with pytest.raises(InvalidUnitEquivalence):
+        data.convert_to_equivalent('erg', 'thermal')
+    with pytest.raises(InvalidUnitEquivalence):
+        data.convert_to_equivalent('m', 'mass_energy')
+    with pytest.raises(InvalidUnitEquivalence):
+        data.to_equivalent('erg', 'thermal')
+    with pytest.raises(InvalidUnitEquivalence):
+        data.to_equivalent('m', 'mass_energy')
+
     # Mass-energy
 
     mp = u.mp.copy()
@@ -1433,34 +1554,47 @@ def test_electromagnetic():
 
     # Various tests of SI and CGS electromagnetic units
 
-    qp_cgs = u.qp_cgs.in_units("C", "SI")
+    t = 1*u.Tesla
+    g = 1*u.gauss
+    assert t.to('gauss') == 1e4*u.gauss
+    assert g.to('T') == 1e-4*u.Tesla
+
+    qp_cgs = u.qp_cgs.in_units("C")
     assert_equal(qp_cgs.units.dimensions, dimensions.charge_mks)
     assert_almost_equal(qp_cgs.v, 10.0*u.qp.v/speed_of_light_cm_per_s)
     qp = 1*u.qp_cgs
-    qp.convert_to_units("C", "SI")
+    qp.convert_to_units("C")
     assert_equal(qp.units.dimensions, dimensions.charge_mks)
     assert_almost_equal(qp.v, 10*u.qp.v/u.clight.v)
 
-    qp_cgs = u.qp.in_units("esu", "CGS")
+    qp_cgs = u.qp.in_units("esu")
     assert_array_almost_equal(qp_cgs, u.qp_cgs)
     assert_equal(qp_cgs.units.dimensions, u.qp_cgs.units.dimensions)
     qp = u.qp.copy()
-    qp.convert_to_units('esu', 'CGS')
+    qp.convert_to_units('esu')
     assert_almost_equal(qp_cgs, qp_cgs)
     assert qp.units == u.esu.units
-    qp.convert_to_units('C', 'MKS')
+    qp.convert_to_units('C')
     assert_almost_equal(u.qp, qp)
     assert qp.units == u.C.units
 
-    qp_mks_k = u.qp_cgs.in_units("kC", "SI")
+    qp_mks_k = u.qp_cgs.in_units("kC")
     assert_array_almost_equal(
         qp_mks_k.v, 1.0e-2*u.qp_cgs.v/speed_of_light_cm_per_s)
     qp = 1*u.qp_cgs
-    qp.convert_to_units('kC', 'SI')
+    qp.convert_to_units('kC')
     assert_almost_equal(qp, qp_mks_k)
 
     B = 1.0*u.T
-    B_cgs = B.in_units("gauss", "CGS")
+    B_cgs = B.in_units("gauss")
+    assert_equal(B.units.dimensions, dimensions.magnetic_field_mks)
+    assert_equal(B_cgs.units.dimensions, dimensions.magnetic_field_cgs)
+    assert_array_almost_equal(B_cgs, unyt_quantity(1.0e4, "gauss"))
+    B_cgs = B.in_cgs()
+    assert_equal(B.units.dimensions, dimensions.magnetic_field_mks)
+    assert_equal(B_cgs.units.dimensions, dimensions.magnetic_field_cgs)
+    assert_array_almost_equal(B_cgs, unyt_quantity(1.0e4, "gauss"))
+    B_cgs = B.in_base('cgs')
     assert_equal(B.units.dimensions, dimensions.magnetic_field_mks)
     assert_equal(B_cgs.units.dimensions, dimensions.magnetic_field_cgs)
     assert_array_almost_equal(B_cgs, unyt_quantity(1.0e4, "gauss"))
@@ -1477,18 +1611,24 @@ def test_electromagnetic():
     assert_equal(u_mks.units.dimensions, dimensions.pressure)
     u_cgs = B_cgs*B_cgs/(8*np.pi)
     assert_equal(u_cgs.units.dimensions, dimensions.pressure)
-    assert_almost_equal(u_mks.in_cgs(), u_cgs)
+    with pytest.raises(UnitConversionError):
+        u_cgs.to(u_mks.units)
+    with pytest.raises(UnitConversionError):
+        u_mks.to(u_cgs.units)
+    with pytest.raises(UnitsNotReducible):
+        u_mks.in_cgs()
+    with pytest.raises(UnitsNotReducible):
+        u_cgs.in_mks()
 
     current = 1.0*u.A
-    I_cgs = current.in_units("statA", equivalence="CGS")
+    I_cgs = current.in_units("statA")
     assert_array_almost_equal(
         I_cgs, unyt_quantity(0.1*speed_of_light_cm_per_s, "statA"))
-    assert_array_almost_equal(
-        I_cgs.in_units("mA", equivalence="SI"), current.in_units("mA"))
+    assert_array_almost_equal(I_cgs.in_units("mA"), current.in_units("mA"))
     assert_equal(I_cgs.units.dimensions, dimensions.current_cgs)
-    current.convert_to_units('statA', equivalence='CGS')
+    current.convert_to_units('statA')
     assert current.units == u.statA.units
-    current.convert_to_units('A', equivalence='MKS')
+    current.convert_to_units('A')
     assert current.units == u.A.units
     I_cgs2 = current.to('statA')
     assert I_cgs2.units == u.statA.units
@@ -1497,18 +1637,36 @@ def test_electromagnetic():
 
     current = 1.0*u.A
     R = unyt_quantity(1.0, "ohm")
-    R_cgs = R.in_units("statohm", "CGS")
+    R_cgs = R.in_units("statohm")
     P_mks = current*current*R
     P_cgs = I_cgs*I_cgs*R_cgs
     assert_equal(P_mks.units.dimensions, dimensions.power)
     assert_equal(P_cgs.units.dimensions, dimensions.power)
-    assert_array_almost_equal(P_cgs.in_cgs(), P_mks.in_cgs())
-    assert_array_almost_equal(P_cgs.in_mks(), unyt_quantity(1.0, "W"))
+    with pytest.raises(UnitsNotReducible):
+        P_cgs.in_cgs()
+    with pytest.raises(UnitsNotReducible):
+        P_mks.in_cgs()
 
     V = unyt_quantity(1.0, "statV")
-    V_mks = V.in_units("V", "SI")
+    V_mks = V.in_units("V")
     assert_array_almost_equal(V_mks.v, 1.0e8*V.v/speed_of_light_cm_per_s)
 
+    data = 1.0*u.C*u.T*u.V
+    with pytest.raises(UnitConversionError):
+        data.to('statC*G*statV')
+    with pytest.raises(UnitConversionError):
+        data.convert_to_units('statC*G*statV')
+    with pytest.raises(UnitsNotReducible):
+        data.in_cgs()
+
+    data = 1.0*u.statC*u.G*u.statV
+    with pytest.raises(UnitConversionError):
+        data.to('C*T*V')
+    with pytest.raises(UnitConversionError):
+        data.convert_to_units('C*T*V')
+    with pytest.raises(UnitsNotReducible):
+        data.in_mks()
+
 
 def test_ytarray_coercion():
     a = unyt_array([1, 2, 3], 'cm')
@@ -1525,9 +1683,14 @@ def test_ytarray_coercion():
 def test_numpy_wrappers():
     a1 = unyt_array([1, 2, 3], 'cm')
     a2 = unyt_array([2, 3, 4, 5, 6], 'cm')
+    a3 = unyt_array([[1, 2, 3], [4, 5, 6]], 'cm')
     catenate_answer = [1, 2, 3, 2, 3, 4, 5, 6]
     intersect_answer = [2, 3]
     union_answer = [1, 2, 3, 4, 5, 6]
+    cross_answer = [-2, 4, -2]
+    norm_answer = np.sqrt(1**2 + 2**2 + 3**2)
+    arr_norm_answer = [norm_answer, np.sqrt(4**2 + 5**2 + 6**2)]
+    dot_answer = 14
 
     assert_array_equal(unyt_array(catenate_answer, 'cm'),
                        uconcatenate((a1, a2)))
@@ -1540,6 +1703,23 @@ def test_numpy_wrappers():
     assert_array_equal(unyt_array(union_answer, 'cm'), uunion1d(a1, a2))
     assert_array_equal(union_answer, np.union1d(a1, a2))
 
+    assert_array_equal(unyt_array(cross_answer, 'cm**2'),
+                       ucross(a1, a1+(2*a1.units)))
+    assert_array_equal(cross_answer, np.cross(a1.v, a1.v+2))
+
+    assert_array_equal(unorm(a1), unyt_quantity(norm_answer, 'cm'))
+    assert_array_equal(np.linalg.norm(a1), norm_answer)
+    assert_array_equal(unorm(a3, axis=1), unyt_array(arr_norm_answer, 'cm'))
+    assert_array_equal(np.linalg.norm(a3, axis=1), arr_norm_answer)
+
+    assert_array_equal(udot(a1, a1), unyt_quantity(dot_answer, 'cm**2'))
+
+    assert_array_equal(np.array(catenate_answer), uconcatenate((a1.v, a2.v)))
+    with pytest.raises(RuntimeError):
+        uconcatenate((a1, a2.v))
+    with pytest.raises(RuntimeError):
+        uconcatenate((a1.to('m'), a2))
+
 
 def test_dimensionless_conversion():
     a = unyt_quantity(1, 'Zsun')
@@ -1564,7 +1744,7 @@ def test_modified_unit_division():
     assert ret.units.base_value == 1.0
 
 
-def test_load_and_save():
+def test_loadtxt_and_savetxt():
     tmpdir = tempfile.mkdtemp()
     curdir = os.getcwd()
     os.chdir(tmpdir)
@@ -1580,6 +1760,47 @@ def test_load_and_save():
     assert_array_equal(b, d)
     assert_array_equal(c, e)
 
+    # adding newlines to the file doesn't matter
+
+    savetxt("arrays.dat", [a, b, c], delimiter=",")
+
+    with open('arrays.dat', 'r+') as f:
+        content = f.read()
+        f.seek(0, 0)
+        f.write('\n' + content)
+
+    d, e = loadtxt("arrays.dat", usecols=(1, 2), delimiter=",")
+
+    assert_array_equal(b, d)
+    assert_array_equal(c, e)
+
+    # data saved by numpy savetxt are loaded without units
+
+    np.savetxt("arrays.dat", np.squeeze(np.transpose([a.v, b.v, c.v])),
+               delimiter=",")
+
+    d, e = loadtxt("arrays.dat", usecols=(1, 2), delimiter=",")
+
+    assert_array_equal(b.v, d)
+    assert_array_equal(c.v, e)
+
+    # save a single array
+
+    savetxt("arrays.dat", a)
+
+    d = loadtxt("arrays.dat")
+
+    assert_array_equal(a, d)
+
+    # save an array with no units and an array with units with a header
+
+    savetxt("arrays.dat", [a.v, b], header='this is a header!')
+
+    d, e = loadtxt("arrays.dat")
+
+    assert_array_equal(a.v, d)
+    assert_array_equal(b, e)
+
     os.chdir(curdir)
     shutil.rmtree(tmpdir)
 
@@ -1622,3 +1843,62 @@ def test_ones_and_zeros_like():
     assert_equal(zd.units, data.units)
     assert_equal(od, unyt_array([1, 1, 1], 'cm'))
     assert_equal(od.units, data.units)
+
+
+def test_coerce_iterable():
+    from unyt import cm, km
+
+    a = unyt_array([1, 2, 3], 'cm')
+    b = [1*cm, 2*km, 3*cm]
+
+    with pytest.raises(IterableUnitCoercionError):
+        a + b
+    with pytest.raises(IterableUnitCoercionError):
+        b + a
+
+
+def test_bypass_validation():
+    from unyt import unyt_array, cm, UnitRegistry
+
+    obj = unyt_array(np.array([1., 2., 3.]), cm, bypass_validation=True)
+    assert obj.units is cm
+
+    reg = UnitRegistry()
+    obj = unyt_array(np.array([1., 2., 3.]), cm, registry=reg,
+                     bypass_validation=True)
+    assert obj.units == cm
+    assert obj.units.registry is reg
+
+
+def test_creation():
+    from unyt import cm, UnitRegistry
+
+    data = [1, 2, 3]*cm
+
+    new_data = unyt_array(data)
+
+    assert new_data.units is cm
+    assert_array_equal(new_data.v, np.array([1, 2, 3], dtype='float64'))
+
+    reg = UnitRegistry()
+
+    new_data = unyt_array(data, registry=reg)
+    assert_array_equal(new_data.v, np.array([1, 2, 3], dtype='float64'))
+    assert new_data.units is not cm
+    assert new_data.units == cm
+    assert new_data.units.registry is reg
+
+    new_data = unyt_array([1, 2, 3], cm)
+    assert_array_equal(new_data.v, np.array([1, 2, 3], dtype='float64'))
+    assert new_data.units is cm
+
+    new_data = unyt_array([1, 2, 3], cm, registry=reg)
+    assert_array_equal(new_data.v, np.array([1, 2, 3], dtype='float64'))
+    assert new_data.units is not cm
+    assert new_data.units == cm
+    assert new_data.units.registry is reg
+
+    with pytest.raises(RuntimeError):
+        unyt_quantity('hello', 'cm')
+    with pytest.raises(RuntimeError):
+        unyt_quantity(np.array([1, 2, 3]), 'cm')