Repository: yt-project/unyt
Pull request #2 (2018-04-15 15:22:37)

Change the default unit system to MKS
This PR changes the default unit system of `unyt` to SI/MKS. This requires touching nearly all the physical ratios, most of the physical constants, and the unit lookup table, as well as a bunch of other code.

It also:

* Adds candela (dimensions of luminous intensity). 
* Adds some derived units related to candela.
* Adds units for magnetic flux in MKS and CGS.

I think it's probably appropriate to also add to the docs a bit of a discussion about the SI system, haven't done that yet. Should probably also add some more testing coverage which uses these new units and dimensions.

ngoldbaum: Also sorry it looks like the tests are broken. That might be my fault. I pushed a bunch of changes to the equivalency system last night that I guess you needed to rebase on?

On the plus side there are in-place equivalencies now :) And lots more tests. There are a couple bugfixes that I want to backport for yt.units that the new tests caught.

Patch:
diff --git a/README.rst b/README.rst
index 80a185e..42d50ba 100644
--- a/README.rst
+++ b/README.rst
@@ -41,14 +41,6 @@ The ``unyt`` package provides a subclass of NumPy's ``ndarray`` class that knows
 
 And a whole lot more! See `the documentation <http://unyt.readthedocs.io>`_ for more examples as well as full API docs.
 
-Code of Conduct
----------------
-
-The ``unyt`` package is part of `The yt Project
-<https://yt-project.org>`. Participating in ``unyt`` development therefore
-happens under the auspices of the `yt community code of conduct
-<http://yt-project.org/doc/developing/developing.html#yt-community-code-of-conduct>`. If for any reason you feel that the code of conduct has been violated, please send an e-mail to confidential@yt-project.org with details describing the incident. All emails sent to this address will be treated with the strictest confidence.
-
 License
 -------
 
diff --git a/docs/index.rst b/docs/index.rst
index 9a05e41..4fc21c7 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -19,7 +19,7 @@ that can be directly imported from from the ``unyt`` namespace:
   >>> from unyt import G, Mearth, Rearth
   >>> v_esc = (2*G*Mearth/Rearth)**(1./2)
   >>> print(v_esc.to('km/s'))
-  11.254342299356157 km/s
+  11.254342299356159 km/s
 
 
 .. toctree::
diff --git a/docs/usage.rst b/docs/usage.rst
index 823ceae..80a0767 100644
--- a/docs/usage.rst
+++ b/docs/usage.rst
@@ -83,7 +83,7 @@ end, days::
 
   >>> period = 2*pi*(semimajor_axis**3/(G*Mjup))**0.5
   >>> period
-  unyt_array([0.00483381, 0.00970288, 0.01954837, 0.04559936], 'km**(3/2)*s/cm**(3/2)')
+  unyt_array([ 4.83380797,  9.70288268, 19.54836529, 45.5993645 ], 'km**(3/2)*s/m**(3/2)')
   >>> period.to('d')
   unyt_array([ 1.76919479,  3.55129736,  7.1547869 , 16.68956617], 'd')
 
@@ -154,10 +154,10 @@ where :math:`\rho` is the average density of the Earth.
   ...
   >>> rho = Mearth / (4./3 * pi* Rearth**3)
   >>> print(rho.to('g/cm**3'))
-  5.581225129861077 g/cm**3
+  5.581225129861083 g/cm**3
   >>> tff = np.sqrt(3*pi/(32*G*rho))
   >>> print(tff.to('min'))
-  14.8202885145703 min
+  14.820288514570295 min
 
 If you make a mistake by adding two things that have different dimensions,
 :mod:`unyt` will raise an error to let you know that you have a bug in your
@@ -393,6 +393,8 @@ atoms and molecules:
      time: fs
      temperature: nK
      angle: rad
+     current_mks: A
+     luminous_intensity: cd
     Other Units:
      energy: eV
 
@@ -458,7 +460,7 @@ equivalence to calculate the mass of a black hole with a radius of one AU:
 
   >>> from unyt import AU
   >>> (1.0*AU).to('Msun', equivalence='schwarzschild')
-  unyt_quantity(50658673.46804737, 'Msun')
+  unyt_quantity(50658673.46804734, 'Msun')
 
 Both the methods that convert data in-place and the ones that return a copy
 support optionally specifying equivalence. In addition to the methods described
@@ -484,8 +486,8 @@ For full API documnentation and an autogenerated listing of the built-in
 equivalencies in :mod:`unyt` as well as a short usage example for each, see the
 :mod:`unyt.equivalencies` API listing.
 
-Dealing with code that doesn't use :mod:`unyt`
-++++++++++++++++++++++++++++++++++++++++++++++
+Dealing with code that doesn't understand units
++++++++++++++++++++++++++++++++++++++++++++++++
 
 Optimally, a function will work the same irrespective of whether the data passed in has units attached or not:
 
@@ -566,101 +568,11 @@ If you would like to create a view rather than a copy, you can apply units like
   >>> data_with_units.base is data
   True
 
-Any set of units can be used for either of these operations. For example, if
-you already have an existing array, you could do this to create a new array
-with the same units:
-
-  >>> more_data = [4, 5, 6]*data_with_units.units
-  >>> more_data
-  unyt_array([4., 5., 6.], 'g')
-
-Working with code that uses ``astropy.units``
----------------------------------------------
-
-The :mod:`unyt` library can convert data contained inside of an Astropy
-``Quantity`` instance. It can also produce a ``Quantity`` from an existing
-:class:`unyt_array <unyt.array.unyt_array>` instance. To convert data from
-``astropy.units`` to :mod:`unyt` use the :func:`unyt_array.from_astropy
-<unyt.array.unyt_array.from_astropy>` function:
-
-  >>> from astropy.units import km
-  >>> from unyt import unyt_quantity
-  >>> unyt_quantity.from_astropy(km)
-  unyt_quantity(1., 'km')
-  >>> a = [1, 2, 3]*km
-  >>> a
-  <Quantity [1., 2., 3.] km>
-  >>> unyt_array.from_astropy(a)
-  unyt_array([1., 2., 3.], 'km')
-
-To convert data *to* ``astropy.units`` use the :meth:`unyt_array.to_astropy <unyt.array.unyt_array.to_astropy>` method:
-
-  >>> from unyt import g, cm
-  >>> data = [3, 4, 5]*g/cm**3
-  >>> data.to_astropy()
-  <Quantity [3., 4., 5.] g / cm3>
-  >>> (4*cm).to_astropy()
-  <Quantity 4. cm>
-
-
-Working with code that uses ``Pint``
-------------------------------------
-
-The :mod:`unyt` library can also convert data contained in ``Pint`` ``Quantity``
-instances. To convert data from ``Pint`` to :mod:`unyt`, use the :func:`unyt_array.from_pint <unyt.array.unyt_array.from_pint>` function:
-
-  >>> from pint import UnitRegistry
-  >>> import numpy as np
-  >>> ureg = UnitRegistry()
-  >>> a = np.arange(4)
-  >>> b = ureg.Quantity(a, "erg/cm**3")
-  >>> b
-  <Quantity([0 1 2 3], 'erg / centimeter ** 3')>
-  >>> c = unyt_array.from_pint(b)
-  >>> c
-  unyt_array([0., 1., 2., 3.], 'erg/cm**3')
-
-And to convert data contained in a :class:`unyt_array <unyt.array.unyt_array>`
-instance, use the :meth:`unyt_array.to_pint <unyt.array.unyt_array.to_pint>`
-method:
-
-  >>> from unyt import cm, s
-  >>> a = 4*cm**2/s
-  >>> print(a)
-  4.0 cm**2/s
-  >>> a.to_pint()
-  <Quantity(4.0, 'centimeter ** 2 / second')>
-  >>> b = [1, 2, 3]*cm
-  >>> b.to_pint()
-  <Quantity([1. 2. 3.], 'centimeter')>
-
-
 Integrating :mod:`unyt` Into a Python Library
 +++++++++++++++++++++++++++++++++++++++++++++
 
-The :mod:`unyt` library began life as the unit system for the ``yt`` data
-analysis and visualization package, in the form of ``yt.units``. In this role,
-:mod:`unyt` was deeply integrated into a larger python library. Due to these
-origins, it is straightforward to build applications that ensure unit
-consistency by making use of :mod:`unyt`. Below we discuss a few topics that
-most often come up when integrating :mod:`unyt` into a new or existing Python library.
-
 Unit Registries
 ---------------
 
-Often it is convenient to define new custom units.
-
-Writing Data with Units to Disk
--------------------------------
-
-Pickles
-*******
-
-HDF5 Files
-**********
-
-Text Files
-**********
-
 Performance Considerations
 --------------------------
diff --git a/unyt/_physical_ratios.py b/unyt/_physical_ratios.py
index d13917c..3e7c88e 100644
--- a/unyt/_physical_ratios.py
+++ b/unyt/_physical_ratios.py
@@ -10,15 +10,16 @@ import numpy as np
 # http://nssdc.gsfc.nasa.gov/planetary/factsheet/jupiterfact.html
 
 # Elementary masses
-mass_electron_grams = 9.10938291e-28
-amu_grams = 1.660538921e-24
-mass_hydrogen_grams = 1.007947*amu_grams
+mass_electron_kg = 9.10938291e-31
+amu_kg = 1.660538921e-27
+amu_grams = amu_kg * 1.0e3
+mass_hydrogen_kg = 1.007947*amu_kg
 
 # Solar values (see Mamajek 2012)
 # https://sites.google.com/site/mamajeksstarnotes/bc-scale
-mass_sun_grams = 1.98841586e33
+mass_sun_kg = 1.98841586e30
 temp_sun_kelvin = 5870.0
-luminosity_sun_ergs_per_sec = 3.8270e33
+luminosity_sun_watts = 3.8270e26
 
 # Consistent with solar abundances used in Cloudy
 metallicity_sun = 0.01295
@@ -34,19 +35,21 @@ mpc_per_rearth = 2.06470307893e-16
 mpc_per_rjup = 2.26566120943e-15
 mpc_per_miles = 5.21552871e-20
 mpc_per_km = 3.24077929e-20
-mpc_per_cm = 3.24077929e-25
-kpc_per_cm = mpc_per_cm / mpc_per_kpc
-pc_per_cm = mpc_per_cm / mpc_per_pc
+mpc_per_m = 3.24077929e-23
+kpc_per_m = mpc_per_m / mpc_per_kpc
+pc_per_m = mpc_per_m / mpc_per_pc
 km_per_pc = 3.08567758e13
+cm_per_pc = 3.08567758e18
+cm_per_mpc = 3.08567758e21
 km_per_m = 1e-3
 km_per_cm = 1e-5
 m_per_cm = 1e-2
-ly_per_cm = 1.05702341e-18
-rsun_per_cm = 1.4378145e-11
-rearth_per_cm = 1.56961033e-9  # Mean (volumetric) radius
-rjup_per_cm = 1.43039006737e-10  # Mean (volumetric) radius
-au_per_cm = 6.68458712e-14
-ang_per_cm = 1.0e8
+ly_per_m = 1.05702341e-20
+rsun_per_m = 1.4378145e-9
+rearth_per_m = 1.56961033e-7  # Mean (volumetric) radius
+rjup_per_m = 1.43039006737e-8  # Mean (volumetric) radius
+au_per_m = 6.68458712e-12
+ang_per_m = 1.0e10
 
 m_per_fpc = 0.0324077929
 
@@ -58,18 +61,19 @@ rearth_per_mpc = 1.0 / mpc_per_rearth
 rjup_per_mpc = 1.0 / mpc_per_rjup
 miles_per_mpc = 1.0 / mpc_per_miles
 km_per_mpc = 1.0 / mpc_per_km
-cm_per_mpc = 1.0 / mpc_per_cm
-cm_per_kpc = 1.0 / kpc_per_cm
+m_per_mpc = 1.0 / mpc_per_m
+m_per_kpc = 1.0 / kpc_per_m
+m_per_km = 1.0 / km_per_m
 cm_per_km = 1.0 / km_per_cm
 cm_per_m = 1.0 / m_per_cm
 pc_per_km = 1.0 / km_per_pc
-cm_per_pc = 1.0 / pc_per_cm
-cm_per_ly = 1.0 / ly_per_cm
-cm_per_rsun = 1.0 / rsun_per_cm
-cm_per_rearth = 1.0 / rearth_per_cm
-cm_per_rjup = 1.0 / rjup_per_cm
-cm_per_au = 1.0 / au_per_cm
-cm_per_ang = 1.0 / ang_per_cm
+m_per_pc = 1.0 / pc_per_m
+m_per_ly = 1.0 / ly_per_m
+m_per_rsun = 1.0 / rsun_per_m
+m_per_rearth = 1.0 / rearth_per_m
+m_per_rjup = 1.0 / rjup_per_m
+m_per_au = 1.0 / au_per_m
+m_per_ang = 1.0 / ang_per_m
 
 # time
 # "IAU Style Manual" by G.A. Wilkins, Comm. 5, in IAU Transactions XXB (1989)
@@ -83,18 +87,25 @@ sec_per_min = 60.0
 day_per_year = 365.25
 
 # velocities, accelerations
-speed_of_light_cm_per_s = 2.99792458e10
-standard_gravity_cm_per_s2 = 9.80665e2
+speed_of_light_m_per_s = 2.99792458e8
+speed_of_light_cm_per_s = speed_of_light_m_per_s * 100.0
+standard_gravity_m_per_s2 = 9.80665
 
 # some constants
-newton_cgs = 6.67384e-8
-planck_cgs = 6.62606957e-27
+newton_mks = 6.67384e-11
+planck_mks = 6.62606957e-34
+# permeability of Free Space
+mu_0 = 4.0e-7*np.pi
+# permittivity of Free Space
+eps_0 = 1.0/(speed_of_light_m_per_s**2*mu_0)
 
 # temperature / energy
-boltzmann_constant_erg_per_K = 1.3806488e-16
+boltzmann_constant_J_per_K = 1.3806488e-23
 erg_per_eV = 1.602176562e-12
+J_per_eV = erg_per_eV * 1.0e-7
 erg_per_keV = erg_per_eV * 1.0e3
-K_per_keV = erg_per_keV / boltzmann_constant_erg_per_K
+J_per_keV = J_per_eV * 1.0e3
+K_per_keV = J_per_keV / boltzmann_constant_J_per_K
 keV_per_K = 1.0 / K_per_keV
 keV_per_erg = 1.0 / erg_per_keV
 eV_per_erg = 1.0 / erg_per_eV
@@ -107,17 +118,17 @@ erg_per_s_per_watt = 1e7
 # Standards", in Highlights of Astronomy (I. Appenzeller, ed.), Table 1,
 # Kluwer Academic Publishers, Dordrecht.
 # REMARK: following masses include whole systems (planet + moons)
-mass_jupiter_grams = mass_sun_grams / 1047.3486
-mass_mercury_grams = mass_sun_grams / 6023600.0
-mass_venus_grams = mass_sun_grams / 408523.71
-mass_earth_grams = mass_sun_grams / 328900.56
-mass_mars_grams = mass_sun_grams / 3098708.0
-mass_saturn_grams = mass_sun_grams / 3497.898
-mass_uranus_grams = mass_sun_grams / 22902.98
-mass_neptune_grams = mass_sun_grams / 19412.24
+mass_jupiter_kg = mass_sun_kg / 1047.3486
+mass_mercury_kg = mass_sun_kg / 6023600.0
+mass_venus_kg = mass_sun_kg / 408523.71
+mass_earth_kg = mass_sun_kg / 328900.56
+mass_mars_kg = mass_sun_kg / 3098708.0
+mass_saturn_kg = mass_sun_kg / 3497.898
+mass_uranus_kg = mass_sun_kg / 22902.98
+mass_neptune_kg = mass_sun_kg / 19412.24
 
 # flux
-jansky_cgs = 1.0e-23
+jansky_mks = 1.0e-26
 # Cosmological constants
 # Calculated with H = 100 km/s/Mpc, value given in units of h^2 g cm^-3
 # Multiply by h^2 to get the critical density in units of g cm^-3
@@ -133,15 +144,15 @@ HUGE = 1.0e90
 TINY = 1.0e-40
 
 # Planck units
-hbar_cgs = 0.5*planck_cgs/np.pi
-planck_mass_grams = np.sqrt(hbar_cgs*speed_of_light_cm_per_s/newton_cgs)
-planck_length_cm = np.sqrt(hbar_cgs*newton_cgs/speed_of_light_cm_per_s**3)
-planck_time_s = planck_length_cm / speed_of_light_cm_per_s
-planck_energy_erg = (planck_mass_grams * speed_of_light_cm_per_s *
-                     speed_of_light_cm_per_s)
-planck_temperature_K = planck_energy_erg / boltzmann_constant_erg_per_K
-planck_charge_esu = np.sqrt(hbar_cgs*speed_of_light_cm_per_s)
+hbar_mks = 0.5*planck_mks/np.pi
+planck_mass_kg = np.sqrt(hbar_mks*speed_of_light_m_per_s/newton_mks)
+planck_length_m = np.sqrt(hbar_mks*newton_mks/speed_of_light_m_per_s**3)
+planck_time_s = planck_length_m / speed_of_light_m_per_s
+planck_energy_J = (planck_mass_kg * speed_of_light_m_per_s *
+                   speed_of_light_m_per_s)
+planck_temperature_K = planck_energy_J / boltzmann_constant_J_per_K
+planck_charge_C = np.sqrt(4.0*np.pi*eps_0*hbar_mks*speed_of_light_m_per_s)
 
 # Imperial and other non-metric units
-grams_per_pound = 453.59237
+kg_per_pound = 0.45359237
 pascal_per_atm = 101325.0
diff --git a/unyt/_unit_lookup_table.py b/unyt/_unit_lookup_table.py
index 963276e..acbf600 100644
--- a/unyt/_unit_lookup_table.py
+++ b/unyt/_unit_lookup_table.py
@@ -14,43 +14,42 @@ The default unit symbol lookup table.
 
 from unyt import dimensions
 from unyt._physical_ratios import (
-    cm_per_pc,
-    cm_per_ly,
-    cm_per_au,
-    cm_per_rsun,
-    cm_per_m,
-    erg_per_s_per_watt,
+    m_per_pc,
+    m_per_ly,
+    m_per_au,
+    m_per_rsun,
     watt_per_horsepower,
-    mass_sun_grams,
+    mass_sun_kg,
     sec_per_year,
     sec_per_day,
     sec_per_hr,
     sec_per_min,
     temp_sun_kelvin,
-    luminosity_sun_ergs_per_sec,
+    luminosity_sun_watts,
     metallicity_sun,
-    erg_per_eV,
+    J_per_eV,
+    amu_kg,
     amu_grams,
-    mass_electron_grams,
-    cm_per_ang,
-    jansky_cgs,
-    mass_jupiter_grams,
-    mass_earth_grams,
+    mass_electron_kg,
+    m_per_ang,
+    jansky_mks,
+    mass_jupiter_kg,
+    mass_earth_kg,
     kelvin_per_rankine,
-    speed_of_light_cm_per_s,
-    planck_length_cm,
-    planck_charge_esu,
-    planck_energy_erg,
-    planck_mass_grams,
+    speed_of_light_m_per_s,
+    planck_length_m,
+    planck_charge_C,
+    planck_energy_J,
+    planck_mass_kg,
     planck_temperature_K,
     planck_time_s,
-    mass_hydrogen_grams,
-    grams_per_pound,
-    standard_gravity_cm_per_s2,
+    mass_hydrogen_kg,
+    kg_per_pound,
+    standard_gravity_m_per_s2,
     pascal_per_atm,
-    newton_cgs,
-    cm_per_rearth,
-    cm_per_rjup
+    newton_mks,
+    m_per_rearth,
+    m_per_rjup
 )
 import numpy as np
 
@@ -59,49 +58,56 @@ import numpy as np
 
 default_unit_symbol_lut = {
     # base
-    "g":  (1.0, dimensions.mass, 0.0, r"\rm{g}"),
+    "m": (1.0, dimensions.length, 0.0, r"\rm{m}"),
+    "g":  (1.0e-3, dimensions.mass, 0.0, r"\rm{g}"),
     "s":  (1.0, dimensions.time, 0.0, r"\rm{s}"),
     "K":  (1.0, dimensions.temperature, 0.0, r"\rm{K}"),
     "radian": (1.0, dimensions.angle, 0.0, r"\rm{radian}"),
+    "A": (1.0, dimensions.current_mks, 0.0, r"\rm{A}"),
+    "cd": (1.0, dimensions.luminous_intensity, 0.0, r"\rm{cd}"),
+    "mol": (1.0 / amu_grams, dimensions.dimensionless, 0.0, r"\rm{mol}"),
 
-    # other cgs
-    "dyne": (1.0, dimensions.force, 0.0, r"\rm{dyn}"),
-    "erg":  (1.0, dimensions.energy, 0.0, r"\rm{erg}"),
-    "Ba": (1.0, dimensions.pressure, 0.0, r"\rm{Ba}"),
-    "esu":  (1.0, dimensions.charge_cgs, 0.0, r"\rm{esu}"),
-    "gauss": (1.0, dimensions.magnetic_field_cgs, 0.0, r"\rm{G}"),
+    # some cgs
+    "dyne": (1.0e-5, dimensions.force, 0.0, r"\rm{dyn}"),
+    "erg":  (1.0e-7, dimensions.energy, 0.0, r"\rm{erg}"),
+    "Ba": (0.1, dimensions.pressure, 0.0, r"\rm{Ba}"),
+    "esu":  (1.0e-3**1.5, dimensions.charge_cgs, 0.0, r"\rm{esu}"),
+    "gauss": (0.1**0.5, dimensions.magnetic_field_cgs, 0.0, r"\rm{G}"),
     "degC": (1.0, dimensions.temperature, -273.15, r"^\circ\rm{C}"),
-    "statA": (1.0, dimensions.current_cgs, 0.0, r"\rm{statA}"),
-    "statV": (1.0, dimensions.electric_potential_cgs, 0.0, r"\rm{statV}"),
-    "statohm": (1.0, dimensions.resistance_cgs, 0.0, r"\rm{statohm}"),
+    "statA": (1.0e-3**1.5, dimensions.current_cgs, 0.0, r"\rm{statA}"),
+    "statV": (0.1*1.0e-3**0.5, dimensions.electric_potential_cgs,
+              0.0, r"\rm{statV}"),
+    "statohm": (100.0, dimensions.resistance_cgs, 0.0, r"\rm{statohm}"),
+    "Mx": (1.0e-3**1.5, dimensions.magnetic_flux_cgs, 0.0, r"\rm{Mx}"),
 
     # some SI
-    "m": (1.0e2, dimensions.length, 0.0, r"\rm{m}"),
-    "J": (1.0e7, dimensions.energy, 0.0, r"\rm{J}"),
-    "W": (1.0e7, dimensions.power, 0.0, r"\rm{W}"),
+    "J": (1.0, dimensions.energy, 0.0, r"\rm{J}"),
+    "W": (1.0, dimensions.power, 0.0, r"\rm{W}"),
     "Hz": (1.0, dimensions.rate, 0.0, r"\rm{Hz}"),
-    "N": (1.0e5, dimensions.force, 0.0, r"\rm{N}"),
+    "N": (1.0, dimensions.force, 0.0, r"\rm{N}"),
     "C": (1.0, dimensions.charge_mks, 0.0, r"\rm{C}"),
-    "A": (1.0, dimensions.current_mks, 0.0, r"\rm{A}"),
-    "T": (1000.0, dimensions.magnetic_field_mks, 0.0, r"\rm{T}"),
-    "Pa": (10.0, dimensions.pressure, 0.0, r"\rm{Pa}"),
-    "V": (1.0e7, dimensions.electric_potential_mks, 0.0, r"\rm{V}"),
-    "ohm": (1.0e7, dimensions.resistance_mks, 0.0, r"\Omega"),
+    "T": (1.0, dimensions.magnetic_field_mks, 0.0, r"\rm{T}"),
+    "Pa": (1.0, dimensions.pressure, 0.0, r"\rm{Pa}"),
+    "V": (1.0, dimensions.electric_potential_mks, 0.0, r"\rm{V}"),
+    "ohm": (1.0, dimensions.resistance_mks, 0.0, r"\Omega"),
+    "Wb": (1.0, dimensions.magnetic_flux_mks, 0.0, r"\rm{Wb}"),
+    "lm": (1.0, dimensions.luminous_flux, 0.0, r"\rm{lm}"),
+    "lx": (1.0, dimensions.luminous_flux/dimensions.area, 0.0, r"\rm{lx}"),
 
     # Imperial and other non-metric units
-    "ft": (30.48, dimensions.length, 0.0, r"\rm{ft}"),
-    "yd": (91.44, dimensions.length, 0.0, r"\rm{yd}"),
-    "mile": (160934.4, dimensions.length, 0.0, r"\rm{mile}"),
+    "in": (0.0254, dimensions.length, 0.0, r"\rm{in}"),
+    "ft": (0.3048, dimensions.length, 0.0, r"\rm{ft}"),
+    "yd": (0.9144, dimensions.length, 0.0, r"\rm{yd}"),
+    "mile": (1609.344, dimensions.length, 0.0, r"\rm{mile}"),
     "degF": (kelvin_per_rankine, dimensions.temperature, -459.67,
              "^\circ\rm{F}"),
     "R": (kelvin_per_rankine, dimensions.temperature, 0.0, r"^\circ\rm{R}"),
-    "lbf": (grams_per_pound*standard_gravity_cm_per_s2, dimensions.force, 0.0,
-            r"\rm{lbf}"),
-    "lb": (grams_per_pound, dimensions.mass, 0.0, r"\rm{lb}"),
-    "lbm": (grams_per_pound, dimensions.mass, 0.0, r"\rm{lbm}"),
-    "atm": (pascal_per_atm*10., dimensions.pressure, 0.0, r"\rm{atm}"),
-    "hp": (watt_per_horsepower*erg_per_s_per_watt, dimensions.power, 0.0,
-           r"\rm{hp}"),
+    "lbf": (kg_per_pound*standard_gravity_m_per_s2,
+            dimensions.force, 0.0, r"\rm{lbf}"),
+    "lb": (kg_per_pound, dimensions.mass, 0.0, r"\rm{lb}"),
+    "lbm": (kg_per_pound, dimensions.mass, 0.0, r"\rm{lbm}"),
+    "atm": (pascal_per_atm, dimensions.pressure, 0.0, r"\rm{atm}"),
+    "hp": (watt_per_horsepower, dimensions.power, 0.0, r"\rm{hp}"),
 
     # dimensionless stuff
     "h": (1.0, dimensions.dimensionless, 0.0, r"h"),
@@ -115,26 +121,26 @@ default_unit_symbol_lut = {
     "yr":  (sec_per_year, dimensions.time, 0.0, r"\rm{yr}"),
 
     # Velocities
-    "c": (speed_of_light_cm_per_s, dimensions.velocity, 0.0, r"\rm{c}"),
+    "c": (speed_of_light_m_per_s, dimensions.velocity, 0.0, r"\rm{c}"),
 
     # Solar units
-    "Msun": (mass_sun_grams, dimensions.mass, 0.0, r"M_\odot"),
-    "msun": (mass_sun_grams, dimensions.mass, 0.0, r"M_\odot"),
-    "Rsun": (cm_per_rsun, dimensions.length, 0.0, r"R_\odot"),
-    "rsun": (cm_per_rsun, dimensions.length, 0.0, r"R_\odot"),
-    "R_sun": (cm_per_rsun, dimensions.length, 0.0, r"R_\odot"),
-    "r_sun": (cm_per_rsun, dimensions.length, 0.0, r"R_\odot"),
-    "Lsun": (luminosity_sun_ergs_per_sec, dimensions.power, 0.0, r"L_\odot"),
+    "Msun": (mass_sun_kg, dimensions.mass, 0.0, r"M_\odot"),
+    "msun": (mass_sun_kg, dimensions.mass, 0.0, r"M_\odot"),
+    "Rsun": (m_per_rsun, dimensions.length, 0.0, r"R_\odot"),
+    "rsun": (m_per_rsun, dimensions.length, 0.0, r"R_\odot"),
+    "R_sun": (m_per_rsun, dimensions.length, 0.0, r"R_\odot"),
+    "r_sun": (m_per_rsun, dimensions.length, 0.0, r"R_\odot"),
+    "Lsun": (luminosity_sun_watts, dimensions.power, 0.0, r"L_\odot"),
     "Tsun": (temp_sun_kelvin, dimensions.temperature, 0.0, r"T_\odot"),
     "Zsun": (metallicity_sun, dimensions.dimensionless, 0.0, r"Z_\odot"),
-    "Mjup": (mass_jupiter_grams, dimensions.mass, 0.0, r"M_{\rm{Jup}}"),
-    "Mearth": (mass_earth_grams, dimensions.mass, 0.0, r"M_\oplus"),
+    "Mjup": (mass_jupiter_kg, dimensions.mass, 0.0, r"M_{\rm{Jup}}"),
+    "Mearth": (mass_earth_kg, dimensions.mass, 0.0, r"M_\oplus"),
 
     # astro distances
-    "AU": (cm_per_au, dimensions.length, 0.0, r"\rm{AU}"),
-    "au": (cm_per_au, dimensions.length, 0.0, r"\rm{AU}"),
-    "ly": (cm_per_ly, dimensions.length, 0.0, r"\rm{ly}"),
-    "pc": (cm_per_pc, dimensions.length, 0.0, r"\rm{pc}"),
+    "AU": (m_per_au, dimensions.length, 0.0, r"\rm{AU}"),
+    "au": (m_per_au, dimensions.length, 0.0, r"\rm{AU}"),
+    "ly": (m_per_ly, dimensions.length, 0.0, r"\rm{ly}"),
+    "pc": (m_per_pc, dimensions.length, 0.0, r"\rm{pc}"),
 
     # angles
     "degree": (np.pi/180., dimensions.angle, 0.0, r"\rm{deg}"),  # degrees
@@ -150,51 +156,52 @@ default_unit_symbol_lut = {
     "lon": (np.pi/180.0, dimensions.angle, -180.0, r"\rm{Longitude}"),
 
     # misc
-    "eV": (erg_per_eV, dimensions.energy, 0.0, r"\rm{eV}"),
-    "amu": (amu_grams, dimensions.mass, 0.0, r"\rm{amu}"),
-    "angstrom": (cm_per_ang, dimensions.length, 0.0, r"\AA"),
-    "Jy": (jansky_cgs, dimensions.specific_flux, 0.0, r"\rm{Jy}"),
+    "eV": (J_per_eV, dimensions.energy, 0.0, r"\rm{eV}"),
+    "amu": (amu_kg, dimensions.mass, 0.0, r"\rm{amu}"),
+    "angstrom": (m_per_ang, dimensions.length, 0.0, r"\AA"),
+    "Jy": (jansky_mks, dimensions.specific_flux, 0.0, r"\rm{Jy}"),
     "counts": (1.0, dimensions.dimensionless, 0.0, r"\rm{counts}"),
     "photons": (1.0, dimensions.dimensionless, 0.0, r"\rm{photons}"),
-    "me": (mass_electron_grams, dimensions.mass, 0.0, r"m_e"),
-    "mp": (mass_hydrogen_grams, dimensions.mass, 0.0, r"m_p"),
-    "mol": (1.0/amu_grams, dimensions.dimensionless, 0.0, r"\rm{mol}"),
-    'Sv': (cm_per_m**2, dimensions.specific_energy, 0.0, r"\rm{Sv}"),
-    "rayleigh": (0.25e6/np.pi, dimensions.count_intensity, 0.0, r"\rm{R}"),
+    "me": (mass_electron_kg, dimensions.mass, 0.0, r"m_e"),
+    "mp": (mass_hydrogen_kg, dimensions.mass, 0.0, r"m_p"),
+    'Sv': (1.0, dimensions.specific_energy, 0.0, r"\rm{Sv}"),
+    "rayleigh": (2.5e9/np.pi, dimensions.count_intensity, 0.0, r"\rm{R}"),
+    "lambert": (1.0e4/np.pi, dimensions.luminance, 0.0, r"\rm{L}"),
+    "nt": (1.0, dimensions.luminance, 0.0, r"\rm{nt}"),
 
     # for AstroPy compatibility
-    "solMass": (mass_sun_grams, dimensions.mass, 0.0, r"M_\odot"),
-    "solRad": (cm_per_rsun, dimensions.length, 0.0, r"R_\odot"),
-    "solLum": (luminosity_sun_ergs_per_sec, dimensions.power, 0.0, r"L_\odot"),
-    "dyn": (1.0, dimensions.force, 0.0, r"\rm{dyn}"),
+    "solMass": (mass_sun_kg, dimensions.mass, 0.0, r"M_\odot"),
+    "solRad": (m_per_rsun, dimensions.length, 0.0, r"R_\odot"),
+    "solLum": (luminosity_sun_watts, dimensions.power, 0.0, r"L_\odot"),
+    "dyn": (1.0e-5, dimensions.force, 0.0, r"\rm{dyn}"),
     "sr": (1.0, dimensions.solid_angle, 0.0, r"\rm{sr}"),
     "rad": (1.0, dimensions.angle, 0.0, r"\rm{rad}"),
     "deg": (np.pi/180., dimensions.angle, 0.0, r"\rm{deg}"),
-    "Fr":  (1.0, dimensions.charge_cgs, 0.0, r"\rm{Fr}"),
-    "G": (1.0, dimensions.magnetic_field_cgs, 0.0, r"\rm{G}"),
-    "Angstrom": (cm_per_ang, dimensions.length, 0.0, r"\AA"),
-    "statC": (1.0, dimensions.charge_cgs, 0.0, r"\rm{statC}"),
+    "Fr":  (1.0e-3**1.5, dimensions.charge_cgs, 0.0, r"\rm{Fr}"),
+    "G": (0.1**0.5, dimensions.magnetic_field_cgs, 0.0, r"\rm{G}"),
+    "Angstrom": (m_per_ang, dimensions.length, 0.0, r"\AA"),
+    "statC": (1.0e-3**1.5, dimensions.charge_cgs, 0.0, r"\rm{statC}"),
 
     # Planck units
-    "m_pl": (planck_mass_grams, dimensions.mass, 0.0, r"m_{\rm{P}}"),
-    "l_pl": (planck_length_cm, dimensions.length, 0.0, r"\ell_\rm{P}"),
+    "m_pl": (planck_mass_kg, dimensions.mass, 0.0, r"m_{\rm{P}}"),
+    "l_pl": (planck_length_m, dimensions.length, 0.0, r"\ell_\rm{P}"),
     "t_pl": (planck_time_s, dimensions.time, 0.0, r"t_{\rm{P}}"),
     "T_pl": (planck_temperature_K, dimensions.temperature, 0.0, r"T_{\rm{P}}"),
-    "q_pl": (planck_charge_esu, dimensions.charge_cgs, 0.0, r"q_{\rm{P}}"),
-    "E_pl": (planck_energy_erg, dimensions.energy, 0.0, r"E_{\rm{P}}"),
+    "q_pl": (planck_charge_C, dimensions.charge_mks, 0.0, r"q_{\rm{P}}"),
+    "E_pl": (planck_energy_J, dimensions.energy, 0.0, r"E_{\rm{P}}"),
 
     # Geometrized units
-    "m_geom": (mass_sun_grams, dimensions.mass, 0.0, r"M_\odot"),
-    "l_geom": (newton_cgs*mass_sun_grams/speed_of_light_cm_per_s**2,
+    "m_geom": (mass_sun_kg, dimensions.mass, 0.0, r"M_\odot"),
+    "l_geom": (newton_mks*mass_sun_kg/speed_of_light_m_per_s**2,
                dimensions.length, 0.0, r"M_\odot"),
-    "t_geom": (newton_cgs*mass_sun_grams/speed_of_light_cm_per_s**3,
+    "t_geom": (newton_mks*mass_sun_kg/speed_of_light_m_per_s**3,
                dimensions.time, 0.0, r"M_\odot"),
 
     # Some Solar System units
-    "R_earth": (cm_per_rearth, dimensions.length, 0.0, r"R_\oplus"),
-    "r_earth": (cm_per_rearth, dimensions.length, 0.0, r"R_\oplus"),
-    "R_jup": (cm_per_rjup, dimensions.length, 0.0, r"R_\mathrm{Jup}"),
-    "r_jup": (cm_per_rjup, dimensions.length, 0.0, r"R_\mathrm{Jup}"),
+    "R_earth": (m_per_rearth, dimensions.length, 0.0, r"R_\oplus"),
+    "r_earth": (m_per_rearth, dimensions.length, 0.0, r"R_\oplus"),
+    "R_jup": (m_per_rjup, dimensions.length, 0.0, r"R_\mathrm{Jup}"),
+    "r_jup": (m_per_rjup, dimensions.length, 0.0, r"R_\mathrm{Jup}"),
 }
 
 # This dictionary formatting from magnitude package, credit to Juan Reyero.
@@ -255,13 +262,20 @@ prefixable_units = [
     "ohm",
     "statohm",
     "Sv",
+    "mol",
+    "cd",
+    "lm",
+    "lx",
+    "Wb",
+    "Mx",
 ]
 
 default_base_units = {
-    dimensions.mass: 'g',
-    dimensions.length: 'cm',
+    dimensions.mass: 'kg',
+    dimensions.length: 'm',
     dimensions.time: 's',
     dimensions.temperature: 'K',
     dimensions.angle: 'radian',
     dimensions.current_mks: 'A',
+    dimensions.luminous_intensity: 'cd',
 }
diff --git a/unyt/array.py b/unyt/array.py
index 66e4f15..9923d44 100644
--- a/unyt/array.py
+++ b/unyt/array.py
@@ -122,11 +122,10 @@ from unyt.exceptions import (
     InvalidUnitOperation,
     UnitConversionError,
     UnitOperationError,
-    UnitsNotReducible,
 )
 from unyt.equivalencies import (
     equivalence_registry,
-    em_conversions
+    _check_em_conversion,
 )
 try:
     from functools import lru_cache
@@ -621,6 +620,12 @@ class unyt_array(np.ndarray):
         """
         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)
@@ -634,7 +639,7 @@ class unyt_array(np.ndarray):
         else:
             self.convert_to_equivalent(units, equivalence, **kwargs)
 
-    def convert_to_base(self, unit_system="cgs", equivalence=None, **kwargs):
+    def convert_to_base(self, unit_system="mks", equivalence=None, **kwargs):
         """
         Convert the array in-place to the equivalent base units in
         the specified unit system.
@@ -646,7 +651,7 @@ class unyt_array(np.ndarray):
         ----------
         unit_system : string, optional
             The unit system to be used in the conversion. If not specified,
-            the default base units of cgs are used.
+            the default base units of mks are used.
         equivalence : string, optional
             The equivalence you wish to use. To see which equivalencies
             are supported for this object, try the ``list_equivalencies``
@@ -774,6 +779,12 @@ 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)
@@ -875,7 +886,7 @@ class unyt_array(np.ndarray):
         else:
             return v
 
-    def in_base(self, unit_system="cgs", equivalence=None, **kwargs):
+    def in_base(self, unit_system="mks"):
         """
         Creates a copy of this array with the data in the specified unit
         system, and returns it in that system's base units.
@@ -884,7 +895,7 @@ class unyt_array(np.ndarray):
         ----------
         unit_system : string, optional
             The unit system to be used in the conversion. If not specified,
-            the default base units of cgs are used.
+            the default base units of mks are used.
 
         Examples
         --------
@@ -893,14 +904,16 @@ class unyt_array(np.ndarray):
         >>> print(E.in_base("mks"))
         2.5e-07 W
         """
-        try:
+        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
+        else:
             to_units = self.units.get_base_equivalent(unit_system)
-        except UnitsNotReducible:
-            to_units, _ = em_conversions[self.units.dimensions]
-            raise EquivalentDimsError(self.units, to_units, unit_system)
-        return self.in_units(to_units, equivalence=equivalence, **kwargs)
+            equivalence = None
+        return self.in_units(to_units, equivalence=equivalence)
 
-    def in_cgs(self, equivalence=None, **kwargs):
+    def in_cgs(self):
         """
         Creates a copy of this array with the data in the equivalent cgs units,
         and returns it.
@@ -915,15 +928,9 @@ class unyt_array(np.ndarray):
         >>> print((10*Newton/km).in_cgs())
         10.0 g/s**2
         """
-        try:
-            to_units = self.units.get_cgs_equivalent()
-        except UnitsNotReducible:
-            to_units, _ = em_conversions[self.units.dimensions]
-            if equivalence is not None:
-                raise EquivalentDimsError(self.units, to_units, 'cgs')
-        return self.in_units(to_units, equivalence=equivalence, **kwargs)
+        return self.in_base("cgs")
 
-    def in_mks(self, equivalence=None, **kwargs):
+    def in_mks(self):
         """
         Creates a copy of this array with the data in the equivalent mks units,
         and returns it.
@@ -938,13 +945,7 @@ class unyt_array(np.ndarray):
         >>> print((1*mile).in_mks())
         1609.344 m
         """
-        try:
-            to_units = self.units.get_mks_equivalent()
-        except UnitsNotReducible:
-            to_units, _ = em_conversions[self.units.dimensions]
-            if equivalence is not None:
-                raise EquivalentDimsError(self.units, to_units, 'mks')
-        return self.in_units(to_units, equivalence=equivalence, **kwargs)
+        return self.in_base("mks")
 
     def convert_to_equivalent(self, unit, equivalence, **kwargs):
         """
@@ -1012,7 +1013,7 @@ class unyt_array(np.ndarray):
         >>> from unyt import K
         >>> a = 1.0e7*K
         >>> print(a.to_equivalent("keV", "thermal"))
-        0.8617332401096502 keV
+        0.8617332401096504 keV
         """
         conv_unit = Unit(unit, registry=self.units.registry)
         if self.units.same_dimensions_as(conv_unit):
@@ -1568,40 +1569,28 @@ class unyt_array(np.ndarray):
             out = None
             out_func = None
         else:
-            # we need to get both the actual "out" object and a view onto it
-            # in case we need to do in-place operations
             out = kwargs.pop('out')[0]
             out_func = out.view(np.ndarray)
         if len(inputs) == 1:
-            # 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
-            # evaluate the ufunc
             out_arr = func(np.asarray(inp), out=out_func, **kwargs)
             if ufunc in (multiply, divide) and method == 'reduce':
-                # a reduction of a multiply or divide corresponds to
-                # a repeated product which we implement as an exponent
                 power_sign = POWER_SIGN_MAPPING[ufunc]
                 if 'axis' in kwargs and kwargs['axis'] is not None:
                     unit = u**(power_sign*inp.shape[kwargs['axis']])
                 else:
                     unit = u**(power_sign*inp.size)
             else:
-                # get unit of result
                 unit = self._ufunc_registry[ufunc](u)
-            # use type(self) here so we can support user-defined
-            # subclasses of unyt_array
             ret_class = type(self)
         elif len(inputs) == 2:
-            # binary ufuncs
             i0 = inputs[0]
             i1 = inputs[1]
-            # coerce inputs to be ndarrays if they aren't already
             inp0 = _coerce_iterable_units(i0)
             inp1 = _coerce_iterable_units(i1)
             u0 = getattr(i0, 'units', None) or getattr(inp0, 'units', None)
@@ -1628,13 +1617,6 @@ class unyt_array(np.ndarray):
                                  _arctan2_unit):
                 # check "is" equality first for speed
                 if u0 is not u1 and u0 != u1:
-                    # we allow adding, multiplying, comparisons with
-                    # zero-filled arrays, lists, etc or scalar zero. We
-                    # do not allow zero-filled unyt_array instances for
-                    # performance reasons. If we did allow it, every
-                    # binary operation would need to scan over all the
-                    # elements of both arrays to check for arrays filled
-                    # with zeros
                     if ((not isinstance(i0, unyt_array) or
                          not isinstance(i1, unyt_array))):
                         any_nonzero = [np.count_nonzero(i0),
@@ -1645,8 +1627,6 @@ class unyt_array(np.ndarray):
                             u1 = u0
                     if not u0.same_dimensions_as(u1):
                         if unit_operator is _comparison_unit:
-                            # we allow comparisons between data with
-                            # units and dimensionless data
                             if u0.is_dimensionless:
                                 u0 = u1
                             elif u1.is_dimensionless:
@@ -1669,9 +1649,7 @@ class unyt_array(np.ndarray):
                             "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
             out_arr = func(inp0.view(np.ndarray), inp1.view(np.ndarray),
                            out=out_func, **kwargs)
             if unit_operator in (_multiply_units, _divide_units):
diff --git a/unyt/dimensions.py b/unyt/dimensions.py
index cf29e17..a1d3665 100644
--- a/unyt/dimensions.py
+++ b/unyt/dimensions.py
@@ -27,12 +27,14 @@ temperature = Symbol("(temperature)", positive=True)
 angle = Symbol("(angle)", positive=True)
 #: current_mks
 current_mks = Symbol("(current_mks)", positive=True)
+#: luminous_intensity
+luminous_intensity = Symbol("(luminous_intensity)", positive=True)
 #: dimensionless
 dimensionless = sympify(1)
 
 #: A list of all of the base dimensions
 base_dimensions = [mass, length, time, temperature, angle, current_mks,
-                   dimensionless]
+                   dimensionless, luminous_intensity]
 
 #
 # Derived dimensions
@@ -90,41 +92,50 @@ specific_energy = energy / mass
 count_flux = 1 / (area*time)
 #: count_intensity
 count_intensity = count_flux / solid_angle
+#: luminous_flux
+luminous_flux = luminous_intensity * solid_angle
+#: luminance
+luminance = luminous_intensity / area
 
 # Gaussian electromagnetic units
 #: charge_cgs
-charge = charge_cgs = (energy * length)**Rational(1, 2)  # proper 1/2 power
+charge_cgs = (energy * length)**Rational(1, 2)  # proper 1/2 power
 #: current_cgs
-current = current_cgs = charge_cgs / time
+current_cgs = charge_cgs / time
 #: electric_field_cgs
-electric_field = electric_field_cgs = charge_cgs / length**2
+electric_field_cgs = charge_cgs / length**2
 #: magnetic_field_cgs
-magnetic_field = magnetic_field_cgs = electric_field_cgs
+magnetic_field_cgs = electric_field_cgs
 #: electric_potential_cgs
-electric_potential = electric_potential_cgs = energy / charge_cgs
+electric_potential_cgs = energy / charge_cgs
 #: resistance_cgs
-resistance = resistance_cgs = electric_potential_cgs / current_cgs
+resistance_cgs = electric_potential_cgs / current_cgs
+#: magnetic_flux_cgs
+magnetic_flux_cgs = magnetic_field_cgs * area
 
 # SI electromagnetic units
 #: charge_mks
-charge_mks = current_mks * time
+charge = charge_mks = current_mks * time
 #: electric_field_mks
-electric_field_mks = force / charge_mks
+electric_field = electric_field_mks = force / charge_mks
 #: magnetic_field_mks
-magnetic_field_mks = electric_field_mks / velocity
+magnetic_field = magnetic_field_mks = electric_field_mks / velocity
 #: electric_potential_mks
-electric_potential_mks = energy / charge_mks
+electric_potential = electric_potential_mks = energy / charge_mks
 #: resistance_mks
-resistance_mks = electric_potential_mks / current_mks
+resistance = resistance_mks = electric_potential_mks / current_mks
+#: magnetic_flux_mks
+magnetic_flux = magnetic_flux_mks = magnetic_field_mks * area
 
 #: a list containing all derived_dimensions
 derived_dimensions = [
     rate, velocity, acceleration, jerk, snap, crackle, pop,
     momentum, force, energy, power, charge_cgs, electric_field_cgs,
     magnetic_field_cgs, solid_angle, flux, specific_flux, volume,
-    area, current_cgs, charge_mks, electric_field_mks,
+    luminous_flux, area, current_cgs, charge_mks, electric_field_mks,
     magnetic_field_mks, electric_potential_cgs, electric_potential_mks,
-    resistance_cgs, resistance_mks]
+    resistance_cgs, resistance_mks, magnetic_flux_mks, magnetic_flux_cgs,
+    luminance]
 
 
 #: a list containing all dimensions
@@ -133,6 +144,7 @@ dimensions = base_dimensions + derived_dimensions
 #: a dict containing a bidirectional mapping from
 #: mks dimension to cgs dimension
 em_dimensions = {magnetic_field_mks: magnetic_field_cgs,
+                 magnetic_flux_mks: magnetic_flux_cgs,
                  charge_mks: charge_cgs,
                  current_mks: current_cgs,
                  electric_potential_mks: electric_potential_cgs,
diff --git a/unyt/equivalencies.py b/unyt/equivalencies.py
index a653f12..ec1b915 100644
--- a/unyt/equivalencies.py
+++ b/unyt/equivalencies.py
@@ -30,6 +30,8 @@ from unyt.dimensions import (
     charge_mks,
     magnetic_field_cgs,
     magnetic_field_mks,
+    magnetic_flux_cgs,
+    magnetic_flux_mks,
     electric_potential_cgs,
     electric_potential_mks,
     electric_field_cgs,
@@ -39,6 +41,7 @@ from unyt.dimensions import (
 )
 
 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
 
@@ -219,19 +222,19 @@ class SpectralEquivalence(Equivalence):
         from unyt import physical_constants as pc
         if new_dims == energy:
             if x.units.dimensions == length:
-                return np.divide(pc.clight*pc.hcgs, x, out=self._get_out(x))
+                return np.divide(pc.clight*pc.hmks, x, out=self._get_out(x))
             elif x.units.dimensions == rate:
-                return np.multiply(x, pc.hcgs, out=self._get_out(x))
+                return np.multiply(x, pc.hmks, out=self._get_out(x))
         elif new_dims == length:
             if x.units.dimensions == rate:
                 return np.divide(pc.clight, x, out=self._get_out(x))
             elif x.units.dimensions == energy:
-                return np.divide(pc.hcgs*pc.clight, x, out=self._get_out(x))
+                return np.divide(pc.hmks*pc.clight, x, out=self._get_out(x))
         elif new_dims == rate:
             if x.units.dimensions == length:
                 return np.divide(pc.clight, x, out=self._get_out(x))
             elif x.units.dimensions == energy:
-                return np.divide(x, pc.hcgs, out=self._get_out(x))
+                return np.divide(x, pc.hmks, out=self._get_out(x))
 
     def __str__(self):
         return "spectral: length <-> frequency <-> energy"
@@ -338,7 +341,7 @@ class LorentzEquivalence(Equivalence):
     >>> from unyt import c, dimensionless
     >>> v = 0.99*c
     >>> print(v.to_equivalent('', 'lorentz'))
-    7.088812050083354 dimensionless
+    7.088812050083393 dimensionless
     >>> fast = 99.9*dimensionless
     >>> fast.to_equivalent('c', 'lorentz')
     unyt_quantity(0.9999499, 'c')
@@ -394,7 +397,7 @@ class SchwarzschildEquivalence(Equivalence):
     >>> (10*Msun).to_equivalent('km', 'schwarzschild')
     unyt_quantity(29.5305543, 'km')
     >>> (1*AU).to_equivalent('Msun', 'schwarzschild')
-    unyt_quantity(50658673.46804737, 'Msun')
+    unyt_quantity(50658673.46804734, 'Msun')
     """
     type_name = "schwarzschild"
     _dims = (mass, length,)
@@ -413,8 +416,12 @@ class SchwarzschildEquivalence(Equivalence):
 
 
 class ComptonEquivalence(Equivalence):
-    """Equivalence between the wavelength change of a compton scattered photon
-    and the mass of the particle it scatters off of.
+    """Equivalence between the Compton wavelength
+    of a particle and its mass.
+
+    .. math::
+
+      \\lambda_c = h/mc
 
     Example
     -------
@@ -431,7 +438,7 @@ class ComptonEquivalence(Equivalence):
 
     def _convert(self, x, new_dims):
         from unyt import physical_constants as pc
-        return np.divide(pc.hcgs/pc.clight, x, out=self._get_out(x))
+        return np.divide(pc.hmks/pc.clight, x, out=self._get_out(x))
 
     def __str__(self):
         return "compton: mass <-> length"
@@ -466,29 +473,52 @@ class EffectiveTemperature(Equivalence):
         if new_dims == flux:
             x4 = np.power(x, 4, out=self._get_out(x))
             return np.multiply(
-                pc.stefan_boltzmann_constant_cgs, x4, out=self._get_out(x))
+                pc.stefan_boltzmann_constant_mks, x4, out=self._get_out(x))
         elif new_dims == temperature:
             T4 = np.divide(
-                x, pc.stefan_boltzmann_constant_cgs, out=self._get_out(x))
+                x, pc.stefan_boltzmann_constant_mks, out=self._get_out(x))
             return np.power(T4, 0.25, out=self._get_out(x))
 
     def __str__(self):
         return "effective_temperature: flux <-> temperature"
 
-
 em_conversions = {
-    charge_mks: ("esu", 0.1*speed_of_light_cm_per_s),
-    magnetic_field_mks: ("gauss", 1.0e4),
-    current_mks: ("statA", 0.1*speed_of_light_cm_per_s),
-    electric_potential_mks: ("statV", 1.0e-8*speed_of_light_cm_per_s),
-    resistance_mks: ("statohm", 1.0e9/(speed_of_light_cm_per_s**2)),
-    charge_cgs: ("C", 10.0/speed_of_light_cm_per_s),
-    magnetic_field_cgs: ("T", 1.0e-4),
-    current_cgs: ("A", 10.0/speed_of_light_cm_per_s),
-    electric_potential_cgs: ("V", 1.0e8/speed_of_light_cm_per_s),
-    resistance_cgs: ("ohm", speed_of_light_cm_per_s**2*1.0e-9),
+    "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
@@ -497,11 +527,14 @@ class ElectromagneticSI(Equivalence):
     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
     -------
@@ -514,11 +547,11 @@ class ElectromagneticSI(Equivalence):
     one_way = True
     _dims = (current_cgs, charge_cgs, magnetic_field_cgs,
              electric_field_cgs, electric_potential_cgs,
-             resistance_cgs)
+             resistance_cgs, magnetic_flux_cgs)
 
     def _convert(self, x, new_dims):
-        old_dims = x.units.dimensions
-        new_units, convert_factor = em_conversions[old_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):
@@ -528,15 +561,17 @@ class ElectromagneticSI(Equivalence):
 class ElectromagneticCGS(Equivalence):
     """An equivalence between SI and CGS electromagnetic units
 
-    Given data in SI electromagnetic units (one of C, T, A, V, or ohm), this
-    equivalency will convert the data to the appropriate CGS electromagnetic
-    unit, using the following mapping:
+    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
     -------
@@ -549,11 +584,11 @@ class ElectromagneticCGS(Equivalence):
     one_way = True
     _dims = (current_mks, charge_mks, magnetic_field_mks,
              electric_field_mks, electric_potential_mks,
-             resistance_mks)
+             resistance_mks, magnetic_flux_mks)
 
     def _convert(self, x, new_dims):
-        old_dims = x.units.dimensions
-        new_units, convert_factor = em_conversions[old_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):
diff --git a/unyt/physical_constants.py b/unyt/physical_constants.py
index 0567340..0535fd2 100644
--- a/unyt/physical_constants.py
+++ b/unyt/physical_constants.py
@@ -20,143 +20,149 @@ For example::
 from math import pi as _pi
 
 from unyt._physical_ratios import (
-    mass_electron_grams as _mass_electron_grams,
-    amu_grams as _amu_grams,
-    mass_hydrogen_grams as _mass_hydrogen_grams,
-    speed_of_light_cm_per_s as _speed_of_light_cm_per_s,
-    boltzmann_constant_erg_per_K as _boltzmann_constant_erg_per_K,
-    mass_sun_grams as _mass_sun_grams,
-    mass_jupiter_grams as _mass_jupiter_grams,
-    mass_mercury_grams as _mass_mercury_grams,
-    mass_venus_grams as _mass_venus_grams,
-    mass_earth_grams as _mass_earth_grams,
-    mass_mars_grams as _mass_mars_grams,
-    mass_saturn_grams as _mass_saturn_grams,
-    mass_uranus_grams as _mass_uranus_grams,
-    mass_neptune_grams as _mass_neptune_grams,
-    planck_mass_grams as _planck_mass_grams,
-    planck_length_cm as _planck_length_cm,
+    mass_electron_kg as _mass_electron_kg,
+    amu_kg as _amu_kg,
+    mass_hydrogen_kg as _mass_hydrogen_kg,
+    speed_of_light_m_per_s as _speed_of_light_m_per_s,
+    boltzmann_constant_J_per_K as _boltzmann_constant_J_per_K,
+    mass_sun_kg as _mass_sun_kg,
+    mass_jupiter_kg as _mass_jupiter_kg,
+    mass_mercury_kg as _mass_mercury_kg,
+    mass_venus_kg as _mass_venus_kg,
+    mass_earth_kg as _mass_earth_kg,
+    mass_mars_kg as _mass_mars_kg,
+    mass_saturn_kg as _mass_saturn_kg,
+    mass_uranus_kg as _mass_uranus_kg,
+    mass_neptune_kg as _mass_neptune_kg,
+    planck_mass_kg as _planck_mass_kg,
+    planck_length_m as _planck_length_m,
     planck_time_s as _planck_time_s,
-    planck_energy_erg as _planck_energy_erg,
-    planck_charge_esu as _planck_charge_esu,
+    planck_energy_J as _planck_energy_J,
+    planck_charge_C as _planck_charge_C,
     planck_temperature_K as _planck_temperature_K,
-    standard_gravity_cm_per_s2 as _standard_gravity_cm_per_s2,
-    newton_cgs as _newton_cgs,
-    planck_cgs as _planck_cgs
+    standard_gravity_m_per_s2 as _standard_gravity_m_per_s2,
+    newton_mks as _newton_mks,
+    planck_mks as _planck_mks,
+    eps_0 as _eps_0,
+    mu_0 as _mu_0
 )
 from unyt.array import unyt_quantity
 
 #: mass of the electron
-mass_electron_cgs = mass_electron = me = \
-    unyt_quantity(_mass_electron_grams, 'g')
+mass_electron_mks = mass_electron = me = \
+    unyt_quantity(_mass_electron_kg, 'kg')
 
 #: atomic mass unit
-amu_cgs = amu = unyt_quantity(_amu_grams, 'g')
+amu_mks = amu = unyt_quantity(_amu_kg, 'kg')
 
 #: Avogadro's number
-Na = avogadros_number = unyt_quantity(6.02214085774*10**23)
+Na = avogadros_number = unyt_quantity(6.02214085774*10**23, 'mol**-1')
 
 #: Mass of hydrogen
-mp = mh = mass_hydrogen = mass_hydrogen_cgs = \
-    unyt_quantity(_mass_hydrogen_grams, 'g')
+mp = mh = mass_hydrogen = mass_hydrogen_mks = \
+    unyt_quantity(_mass_hydrogen_kg, 'kg')
 
 # Velocities
 #: speed of light
-c = clight = speed_of_light = speed_of_light_cgs = \
-    unyt_quantity(_speed_of_light_cm_per_s, 'cm/s')
+c = clight = speed_of_light = speed_of_light_mks = \
+    unyt_quantity(_speed_of_light_m_per_s, 'm/s')
 
 # Cross Sections
 #: Thompson cross section (8*pi/3 (alpha*hbar*c/(2*pi))**2)
 sigma_thompson = thompson_cross_section = cross_section_thompson = \
-    cross_section_thompson_cgs = unyt_quantity(6.65245854533e-25, 'cm**2')
+    cross_section_thompson_mks = unyt_quantity(6.65245854533e-29, 'm**2')
 
 # Charge
 #: Charge of the proton
-qp = elementary_charge = proton_charge = charge_prothon = charge_proton_cgs = \
-    unyt_quantity(4.8032056e-10, 'esu')
+qp = elementary_charge = proton_charge = charge_proton = charge_proton_mks = \
+    unyt_quantity(1.6021766208e-19, 'C')
 #: Charge of the electron
-electron_charge = charge_electron = charge_electron_cgs = -qp
+electron_charge = charge_electron = charge_electron_mks = -qp
+
+qp_cgs = elementary_charge_cgs = proton_charge_cgs = charge_proton_cgs = \
+    unyt_quantity(4.8032056e-10, 'esu')
+electron_charge_cgs = charge_electron_cgs = -charge_proton_cgs
 
 # Physical Constants
 #: Boltzmann constant
-kb = kboltz = boltzmann_constant = boltzmann_constant_cgs = unyt_quantity(
-    _boltzmann_constant_erg_per_K, 'erg/K')
+kb = kboltz = boltzmann_constant = boltzmann_constant_mks = unyt_quantity(
+    _boltzmann_constant_J_per_K, 'J/K')
 
 #: Gravitational constant
-G = newtons_constant = gravitational_constant = gravitational_constant_cgs = \
-    unyt_quantity(_newton_cgs, 'cm**3/g/s**2')
+G = newtons_constant = gravitational_constant = gravitational_constant_mks = \
+    unyt_quantity(_newton_mks, 'm**3/kg/s**2')
 
 #: Planck constant
-hcgs = planck_constant = planck_constant_cgs = \
-    unyt_quantity(_planck_cgs, 'erg*s')
+hmks = planck_constant = planck_constant_mks = \
+    unyt_quantity(_planck_mks, 'J*s')
 
 #: Reduced Planck constant
-hbar = reduced_planck_constant = 0.5*hcgs/_pi
+hbar = reduced_planck_constant = 0.5*hmks/_pi
 
 #: Stefan-Boltzmann constant
-stefan_boltzmann_constant = stefan_boltzmann_constant_cgs = unyt_quantity(
-    5.670373e-5, 'erg/cm**2/s**1/K**4')
+stefan_boltzmann_constant = stefan_boltzmann_constant_mks = unyt_quantity(
+    5.670373e-8, 'W/m**2/K**4')
 
 #: Current CMB temperature
 CMD_temperature = Tcmb = unyt_quantity(2.726, 'K')
 
 # Solar System
 #: Mass of the sun
-msun = solar_mass = mass_sun = mass_sun_cgs = unyt_quantity(
-    _mass_sun_grams, 'g')
+msun = solar_mass = mass_sun = mass_sun_mks = unyt_quantity(
+    _mass_sun_kg, 'kg')
 
 #: Mass of Jupiter
-Mjup = mjup = jupiter_mass = mass_jupiter = mass_jupiter_cgs = unyt_quantity(
-    _mass_jupiter_grams, 'g')
+Mjup = mjup = jupiter_mass = mass_jupiter = mass_jupiter_mks = unyt_quantity(
+    _mass_jupiter_kg, 'kg')
 
 #: Mass of Mercury
-mercury_mass = mass_mercury = mass_mercury_cgs = unyt_quantity(
-    _mass_mercury_grams, 'g')
+mercury_mass = mass_mercury = mass_mercury_mks = unyt_quantity(
+    _mass_mercury_kg, 'kg')
 
 #: Mass of Venus
-venus_mass = mass_venus = mass_venus_cgs = unyt_quantity(
-    _mass_venus_grams, 'g')
+venus_mass = mass_venus = mass_venus_mks = unyt_quantity(
+    _mass_venus_kg, 'kg')
 
 #: Mass of Earth
-Mearth = mearth = earth_mass = mass_earth = mass_earth_cgs = unyt_quantity(
-    _mass_earth_grams, 'g')
+Mearth = mearth = earth_mass = mass_earth = mass_earth_mks = unyt_quantity(
+    _mass_earth_kg, 'kg')
 
 #: Mass of Mars
-mars_mass = mass_margs = mass_mars_cgs = unyt_quantity(_mass_mars_grams, 'g')
+mars_mass = mass_mars = mass_mars_mks = unyt_quantity(_mass_mars_kg, 'kg')
 
 #: Mass of Saturn
-saturn_mass = mass_saturn = mass_saturn_cgs = unyt_quantity(
-    _mass_saturn_grams, 'g')
+saturn_mass = mass_saturn = mass_saturn_mks = unyt_quantity(
+    _mass_saturn_kg, 'kg')
 
 #: Mass of Uranus
-uranus_mass = mass_uranus = mass_uranus_cgs = unyt_quantity(
-    _mass_uranus_grams, 'g')
+uranus_mass = mass_uranus = mass_uranus_mks = unyt_quantity(
+    _mass_uranus_kg, 'kg')
 
 #: Mass of Neptune
-neptune_mass = mass_neptune = mass_neptune_cgs = unyt_quantity(
-    _mass_neptune_grams, 'g')
+neptune_mass = mass_neptune = mass_neptune_mks = unyt_quantity(
+    _mass_neptune_kg, 'kg')
 
 # Planck units
 #: Planck mass
-m_pl = planck_mass = unyt_quantity(_planck_mass_grams, "g")
+m_pl = planck_mass = unyt_quantity(_planck_mass_kg, "g")
 #: Planck length
-l_pl = planck_length = unyt_quantity(_planck_length_cm, "cm")
+l_pl = planck_length = unyt_quantity(_planck_length_m, "cm")
 #: Planck time
 t_pl = planck_time = unyt_quantity(_planck_time_s, "s")
 #: Planck energy
-E_pl = planck_energy = unyt_quantity(_planck_energy_erg, "erg")
+E_pl = planck_energy = unyt_quantity(_planck_energy_J, "erg")
 #: Planck charge
-q_pl = planck_charge = unyt_quantity(_planck_charge_esu, "esu")
+q_pl = planck_charge = unyt_quantity(_planck_charge_C, "C")
 #: Planck tempearture
 T_pl = planck_temperature = unyt_quantity(_planck_temperature_K, "K")
 
-# MKS E&M units
+# MKS E&M constants
 #: Permeability of Free Space
-mu_0 = unyt_quantity(4.0e-7*_pi, "N/A**2")
+mu_0 = unyt_quantity(_mu_0, "N/A**2")
 #: Permittivity of Free Space
-eps_0 = (1.0/(clight**2*mu_0)).in_units("C**2/N/m**2")
+eps_0 = unyt_quantity(_eps_0, "C**2/N/m**2")
 
 # Misc
 #: Standard gravitational acceleration
-standard_gravity = standard_gravity_cgs = unyt_quantity(
-    _standard_gravity_cm_per_s2, "cm/s**2")
+standard_gravity = standard_gravity_mks = unyt_quantity(
+    _standard_gravity_m_per_s2, "m/s**2")
diff --git a/unyt/unit_object.py b/unyt/unit_object.py
index 3f49545..167d83c 100644
--- a/unyt/unit_object.py
+++ b/unyt/unit_object.py
@@ -57,10 +57,11 @@ from unyt.unit_registry import (
 from unyt.exceptions import (
     InvalidUnitOperation,
     MissingMKSCurrent,
+    UnitsNotReducible,
 )
 from unyt.equivalencies import (
     equivalence_registry,
-    em_conversions,
+    _check_em_conversion,
 )
 import copy
 import token
@@ -632,7 +633,7 @@ class Unit(Expr):
         old_dims = self.dimensions
         return old_dims in this_equiv._dims
 
-    def get_base_equivalent(self, unit_system="cgs"):
+    def get_base_equivalent(self, unit_system="mks"):
         """Create and return dimensionally-equivalent units in a specified base.
 
         >>> from unyt import g, cm
@@ -648,9 +649,13 @@ class Unit(Expr):
             unit_system = self.registry.unit_system_id
         unit_system = unit_system_registry[str(unit_system)]
         try:
-            units_string = unit_system[self.dimensions]
+            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:
-            units_string = em_conversions[self.dimensions][0]
+            raise UnitsNotReducible(self.units, unit_system)
         return Unit(units_string, registry=self.registry)
 
     def get_cgs_equivalent(self):
@@ -933,7 +938,7 @@ def _define_unit(registry, symbol, value, tex_repr=None, offset=None,
             value = unyt_quantity(value[0], value[1])
         else:
             raise RuntimeError("\"value\" needs to be a (value, unit) tuple!")
-    base_value = float(value.in_base(unit_system='cgs-ampere'))
+    base_value = float(value.in_base(unit_system='mks'))
     dimensions = value.units.dimensions
     registry.add(symbol, base_value, dimensions, tex_repr=tex_repr,
                  offset=offset)
diff --git a/unyt/unit_registry.py b/unyt/unit_registry.py
index 69681e6..2dc5fdf 100644
--- a/unyt/unit_registry.py
+++ b/unyt/unit_registry.py
@@ -136,7 +136,7 @@ class UnitRegistry:
 
         if hasattr(base_value, "in_base"):
             new_dimensions = base_value.units.dimensions
-            base_value = base_value.in_base('cgs-ampere')
+            base_value = base_value.in_base('mks')
             base_value = base_value.value
         else:
             new_dimensions = self.lut[symbol][1]
diff --git a/unyt/unit_symbols.py b/unyt/unit_symbols.py
index d8b45ca..3fdfdc3 100644
--- a/unyt/unit_symbols.py
+++ b/unyt/unit_symbols.py
@@ -180,6 +180,16 @@ THz = thz = terahertz = Unit("THz")
 Pa = pascal = Unit("Pa")
 #: sievert
 Sv = sievert = Unit("Sv")
+#: candela
+cd = Candela = candela = Unit("cd")
+#: lumen
+lm = Lumen = lumen = Unit("lm")
+#: lux
+lx = Lux = lux = Unit('lx')
+#: lambert
+lambert = Lambert = Unit('lambert')
+#: nit
+nt = nit = Nit = Unit("nt")
 
 #
 # Imperial units
@@ -296,6 +306,8 @@ statampere = statA = Unit("statA")
 statvolt = statV = Unit("statV")
 #: Statohm
 statohm = Unit("statohm")
+#: Maxwell
+Mx = maxwell = Maxwell = Unit("Mx")
 
 #
 # SI electromagnetic units
@@ -311,6 +323,8 @@ A = ampere = Ampere = Unit("A")
 V = volt = Volt = Unit("V")
 #: Ohm
 ohm = Ohm = Unit("ohm")
+#: Weber
+Wb = weber = Weber = Unit("Wb")
 
 #
 # Geographic units
diff --git a/unyt/unit_systems.py b/unyt/unit_systems.py
index 843ace2..11cbb7d 100644
--- a/unyt/unit_systems.py
+++ b/unyt/unit_systems.py
@@ -26,6 +26,7 @@ from unyt import physical_constants as pc
 
 unit_system_registry = {}
 
+cmks = dimensions.current_mks
 
 class _UnitSystemConstants(object):
     """
@@ -64,32 +65,44 @@ class UnitSystem(object):
         The base temperature unit of this unit system. Defaults to "K".
     angle_unit : string, optional
         The base angle unit of this unit system. Defaults to "rad".
+    mks_system: boolean, optional
+        Whether or not this unit system has SI-specific units.
+        Default: False
     current_mks_unit : string, optional
-        The base current unit of this unit system. Only used in MKS
-        or MKS-based unit systems.
+        The base current unit of this unit system. Defaults to "A".
+    luminous_intensity_unit : string, optional
+        The base luminous intensity unit of this unit system.
+        Defaults to "cd".
     registry : :class:`yt.units.unit_registry.UnitRegistry` object
         The unit registry associated with this unit system. Only
         useful for defining unit systems based on code units.
     """
     def __init__(self, name, length_unit, mass_unit, time_unit,
-                 temperature_unit="K", angle_unit="rad", current_mks_unit=None,
+                 temperature_unit="K", angle_unit="rad",
+                 current_mks_unit="A", luminous_intensity_unit="cd",
                  registry=None):
         self.registry = registry
+        if current_mks_unit is not None:
+            current_mks_unit = Unit(current_mks_unit,
+                                    registry=self.registry)
         self.units_map = OrderedDict([
             (dimensions.length, Unit(length_unit, registry=self.registry)),
             (dimensions.mass, Unit(mass_unit, registry=self.registry)),
             (dimensions.time, Unit(time_unit, registry=self.registry)),
             (dimensions.temperature, Unit(
                 temperature_unit, registry=self.registry)),
-            (dimensions.angle, Unit(angle_unit, registry=self.registry))])
+            (dimensions.angle, Unit(angle_unit, registry=self.registry)),
+            (dimensions.current_mks, current_mks_unit),
+            (dimensions.luminous_intensity, Unit(
+                luminous_intensity_unit, registry=self.registry))])
         for dimension, unit in self.units_map.items():
+            # CGS sets the current_mks unit to none, so catch it here
+            if unit is None:
+                continue
             if unit.dimensions is not dimension:
                 raise IllDefinedUnitSystem(self.units_map)
-        self._dims = ["length", "mass", "time", "temperature", "angle"]
-        if current_mks_unit is not None:
-            self.units_map[dimensions.current_mks] = Unit(
-                current_mks_unit, registry=self.registry)
-            self._dims.append("current_mks")
+        self._dims = ["length", "mass", "time", "temperature", "angle",
+                      "current_mks", "luminous_intensity"]
         self.registry = registry
         self.base_units = self.units_map.copy()
         unit_system_registry[name] = self
@@ -101,8 +114,7 @@ class UnitSystem(object):
             key = getattr(dimensions, key)
         um = self.units_map
         if key not in um or um[key].dimensions is not key:
-            cmks = dimensions.current_mks
-            if cmks in key.free_symbols and cmks not in self.units_map:
+            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)
             self.units_map[key] = Unit(units, registry=self.registry)
@@ -113,6 +125,8 @@ class UnitSystem(object):
             if key not in self._dims:
                 self._dims.append(key)
             key = getattr(dimensions, key)
+        if self.units_map[cmks] is None and cmks in key.free_symbols:
+            raise MissingMKSCurrent(self.name)
         self.units_map[key] = Unit(value, registry=self.registry)
 
     def __str__(self):
@@ -122,7 +136,8 @@ class UnitSystem(object):
         repr = "%s Unit System\n" % self.name
         repr += " Base Units:\n"
         for dim in self.base_units:
-            repr += "  %s: %s\n" % (str(dim).strip("()"), self.base_units[dim])
+            if self.base_units[dim] is not None:
+                repr += "  %s: %s\n" % (str(dim).strip("()"), self.base_units[dim])
         repr += " Other Units:\n"
         for key in self._dims:
             dim = getattr(dimensions, key)
@@ -130,9 +145,8 @@ class UnitSystem(object):
                 repr += "  %s: %s\n" % (key, self.units_map[dim])
         return repr[:-1]
 
-
 #: The CGS unit system
-cgs_unit_system = UnitSystem("cgs", "cm", "g", "s")
+cgs_unit_system = UnitSystem("cgs", "cm", "g", "s", current_mks_unit=None)
 cgs_unit_system["energy"] = "erg"
 cgs_unit_system["specific_energy"] = "erg/g"
 cgs_unit_system["pressure"] = "dyne/cm**2"
@@ -143,7 +157,7 @@ cgs_unit_system["current_cgs"] = "statA"
 cgs_unit_system["power"] = "erg/s"
 
 #: The MKS unit system
-mks_unit_system = UnitSystem("mks", "m", "kg", "s", current_mks_unit="A")
+mks_unit_system = UnitSystem("mks", "m", "kg", "s")
 mks_unit_system["energy"] = "J"
 mks_unit_system["specific_energy"] = "J/kg"
 mks_unit_system["pressure"] = "Pa"
@@ -176,15 +190,4 @@ geometrized_unit_system = UnitSystem("geometrized", "l_geom",
 planck_unit_system = UnitSystem("planck", "l_pl", "m_pl", "t_pl",
                                 temperature_unit="T_pl")
 planck_unit_system["energy"] = "E_pl"
-planck_unit_system["charge_cgs"] = "q_pl"
-
-#: A CGS unit system with Ampere (SI-like E&M units)
-cgs_ampere_unit_system = UnitSystem('cgs-ampere', 'cm', 'g', 's',
-                                    current_mks_unit='A')
-cgs_ampere_unit_system["energy"] = "erg"
-cgs_ampere_unit_system["specific_energy"] = "erg/g"
-cgs_ampere_unit_system["pressure"] = "dyne/cm**2"
-cgs_ampere_unit_system["force"] = "dyne"
-cgs_ampere_unit_system["magnetic_field_cgs"] = "gauss"
-cgs_ampere_unit_system["charge_cgs"] = "esu"
-cgs_ampere_unit_system["current_cgs"] = "statA"
+planck_unit_system["charge_mks"] = "q_pl"


Test patch:
diff --git a/unyt/tests/test_unit_systems.py b/unyt/tests/test_unit_systems.py
index bf778ba..29f658f 100644
--- a/unyt/tests/test_unit_systems.py
+++ b/unyt/tests/test_unit_systems.py
@@ -33,7 +33,8 @@ from unyt import (
 def test_unit_systems():
     goofy_unit_system = UnitSystem(
         "goofy", "ly", "lbm", "hr", temperature_unit="R",
-        angle_unit="arcsec", current_mks_unit="mA")
+        angle_unit="arcsec", current_mks_unit="mA",
+        luminous_intensity_unit="cd")
     assert goofy_unit_system["temperature"] == Unit("R")
     assert goofy_unit_system[dimensions.solid_angle] == Unit("arcsec**2")
     assert goofy_unit_system["energy"] == Unit("lbm*ly**2/hr**2")
@@ -47,27 +48,23 @@ def test_unit_system_id():
     reg1 = UnitRegistry()
     reg2 = UnitRegistry()
     assert reg1.unit_system_id == reg2.unit_system_id
-    reg1.modify('m', 2.0)
+    reg1.modify('g', 2.0)
     assert reg1.unit_system_id != reg2.unit_system_id
     reg1 = UnitRegistry()
     reg1.add('dinosaurs', 12.0, dimensions.length)
     assert reg1.unit_system_id != reg2.unit_system_id
     reg1 = UnitRegistry()
-    reg1.remove('m')
+    reg1.remove('g')
     assert reg1.unit_system_id != reg2.unit_system_id
-    reg1.add('m', 100., dimensions.length)
+    reg1.add('g', 1.0e-3, dimensions.mass)
     assert reg1.unit_system_id == reg2.unit_system_id
 
 
 def test_cgs_mks_unit_conversions():
     t = 1*Tesla
     g = 1*gauss
-    with pytest.raises(EquivalentDimsError):
-        t.to('G')
     assert t.to_equivalent('G', 'cgs') == 1e4*gauss
     assert t.to_equivalent('G', 'CGS') == 1e4*gauss
-    with pytest.raises(EquivalentDimsError):
-        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
diff --git a/unyt/tests/test_units.py b/unyt/tests/test_units.py
index ed17051..1c5f081 100644
--- a/unyt/tests/test_units.py
+++ b/unyt/tests/test_units.py
@@ -18,6 +18,7 @@ Test symbolic unit handling.
 import numpy as np
 from numpy.testing import (
     assert_almost_equal,
+    assert_allclose,
     assert_array_almost_equal_nulp,
     assert_equal
 )
@@ -26,14 +27,14 @@ import pytest
 from sympy import Symbol
 from unyt._testing import assert_allclose_units
 from unyt.unit_registry import UnitRegistry
-from unyt import electrostatic_unit, elementary_charge
+from unyt import electrostatic_unit, elementary_charge_cgs
 from unyt.dimensions import (
     mass,
     length,
     time,
     temperature,
     energy,
-    magnetic_field,
+    magnetic_field_cgs,
     power,
     rate
 )
@@ -51,11 +52,11 @@ from unyt._unit_lookup_table import (
 )
 import unyt.unit_symbols as unit_symbols
 from unyt._physical_ratios import (
-    cm_per_pc,
+    m_per_pc,
     sec_per_year,
-    cm_per_km,
-    cm_per_mpc,
-    mass_sun_grams
+    m_per_km,
+    m_per_mpc,
+    mass_sun_kg
 )
 
 
@@ -112,35 +113,35 @@ def test_create_from_string():
 
     """
 
-    u1 = Unit("g * cm**2 * s**-2")
+    u1 = Unit("kg * m**2 * s**-2")
     assert u1.dimensions == energy
     assert u1.base_value == 1.0
 
     # make sure order doesn't matter
-    u2 = Unit("cm**2 * s**-2 * g")
+    u2 = Unit("m**2 * s**-2 * kg")
     assert u2.dimensions == energy
     assert u2.base_value == 1.0
 
     # Test rationals
-    u3 = Unit("g**0.5 * cm**-0.5 * s**-1")
-    assert u3.dimensions == magnetic_field
+    u3 = Unit("kg**0.5 * m**-0.5 * s**-1")
+    assert u3.dimensions == magnetic_field_cgs
     assert u3.base_value == 1.0
 
     # sqrt functions
-    u4 = Unit("sqrt(g)/sqrt(cm)/s")
-    assert u4.dimensions == magnetic_field
+    u4 = Unit("sqrt(kg)/sqrt(m)/s")
+    assert u4.dimensions == magnetic_field_cgs
     assert u4.base_value == 1.0
 
     # commutative sqrt function
-    u5 = Unit("sqrt(g/cm)/s")
-    assert u5.dimensions == magnetic_field
+    u5 = Unit("sqrt(kg/m)/s")
+    assert u5.dimensions == magnetic_field_cgs
     assert u5.base_value == 1.0
 
     # nonzero CGS conversion factor
     u6 = Unit("Msun/pc**3")
     assert u6.dimensions == mass/length**3
     assert_array_almost_equal_nulp(np.array([u6.base_value]),
-                                   np.array([mass_sun_grams/cm_per_pc**3]))
+                                   np.array([mass_sun_kg/m_per_pc**3]))
 
     with pytest.raises(UnitParseError):
         Unit('m**m')
@@ -157,8 +158,8 @@ def test_create_from_expr():
     Create units from sympy Exprs and check attributes.
 
     """
-    pc_cgs = cm_per_pc
-    yr_cgs = sec_per_year
+    pc_mks = m_per_pc
+    yr_mks = sec_per_year
 
     # Symbol expr
     s1 = Symbol("pc", positive=True)
@@ -178,10 +179,10 @@ def test_create_from_expr():
     assert u3.expr == s3
     assert u4.expr == s4
 
-    assert_allclose_units(u1.base_value, pc_cgs, 1e-12)
-    assert_allclose_units(u2.base_value, yr_cgs, 1e-12)
-    assert_allclose_units(u3.base_value, pc_cgs * yr_cgs, 1e-12)
-    assert_allclose_units(u4.base_value, pc_cgs**2 / yr_cgs, 1e-12)
+    assert_allclose_units(u1.base_value, pc_mks, 1e-12)
+    assert_allclose_units(u2.base_value, yr_mks, 1e-12)
+    assert_allclose_units(u3.base_value, pc_mks * yr_mks, 1e-12)
+    assert_allclose_units(u4.base_value, pc_mks**2 / yr_mks, 1e-12)
 
     assert u1.dimensions == length
     assert u2.dimensions == time
@@ -195,15 +196,15 @@ def test_create_with_duplicate_dimensions():
 
     """
 
-    u1 = Unit("erg * s**-1")
+    u1 = Unit("J * s**-1")
     u2 = Unit("km/s/Mpc")
-    km_cgs = cm_per_km
-    Mpc_cgs = cm_per_mpc
+    km_mks = m_per_km
+    Mpc_mks = m_per_mpc
 
     assert u1.base_value == 1
     assert u1.dimensions == power
 
-    assert_allclose_units(u2.base_value, km_cgs / Mpc_cgs, 1e-12)
+    assert_allclose_units(u2.base_value, km_mks / Mpc_mks, 1e-12)
     assert u2.dimensions == rate
 
 
@@ -329,8 +330,8 @@ def test_multiplication():
     Multiply two units.
 
     """
-    msun_cgs = mass_sun_grams
-    pc_cgs = cm_per_pc
+    msun_mks = mass_sun_kg
+    pc_mks = m_per_pc
 
     # Create symbols
     msun_sym = Symbol("Msun", positive=True)
@@ -345,7 +346,7 @@ def test_multiplication():
     u3 = u1 * u2
 
     assert u3.expr == msun_sym * pc_sym
-    assert_allclose_units(u3.base_value, msun_cgs * pc_cgs, 1e-12)
+    assert_allclose_units(u3.base_value, msun_mks * pc_mks, 1e-12)
     assert u3.dimensions == mass * length
 
     # Pow and Mul operations
@@ -355,7 +356,7 @@ def test_multiplication():
     u6 = u4 * u5
 
     assert u6.expr == pc_sym**2 * msun_sym * s_sym
-    assert_allclose_units(u6.base_value, pc_cgs**2 * msun_cgs, 1e-12)
+    assert_allclose_units(u6.base_value, pc_mks**2 * msun_mks, 1e-12)
     assert u6.dimensions == length**2 * mass * time
 
 
@@ -364,8 +365,8 @@ def test_division():
     Divide two units.
 
     """
-    pc_cgs = cm_per_pc
-    km_cgs = cm_per_km
+    pc_mks = m_per_pc
+    km_mks = m_per_km
 
     # Create symbols
     pc_sym = Symbol("pc", positive=True)
@@ -379,7 +380,7 @@ def test_division():
     u3 = u1 / u2
 
     assert u3.expr == pc_sym / (km_sym * s_sym)
-    assert_allclose_units(u3.base_value, pc_cgs / km_cgs, 1e-12)
+    assert_allclose_units(u3.base_value, pc_mks / km_mks, 1e-12)
     assert u3.dimensions == 1 / time
 
 
@@ -390,20 +391,20 @@ def test_power():
     """
     from sympy import nsimplify
 
-    pc_cgs = cm_per_pc
-    mK_cgs = 1e-3
+    pc_mks = m_per_pc
+    mK_mks = 1e-3
     u1_dims = mass * length**2 * time**-3 * temperature**4
-    u1 = Unit("g * pc**2 * s**-3 * mK**4")
+    u1 = Unit("kg * pc**2 * s**-3 * mK**4")
 
     u2 = u1**2
 
     assert u2.dimensions == u1_dims**2
-    assert_allclose_units(u2.base_value, (pc_cgs**2 * mK_cgs**4)**2, 1e-12)
+    assert_allclose_units(u2.base_value, (pc_mks**2 * mK_mks**4)**2, 1e-12)
 
     u3 = u1**(-1.0/3)
 
     assert u3.dimensions == nsimplify(u1_dims**(-1.0/3))
-    assert_allclose_units(u3.base_value, (pc_cgs**2 * mK_cgs**4)**(-1.0/3),
+    assert_allclose_units(u3.base_value, (pc_mks**2 * mK_mks**4)**(-1.0/3),
                           1e-12)
 
 
@@ -424,17 +425,17 @@ def test_base_equivalent():
     Check base equivalent of a unit.
 
     """
-    Msun_cgs = mass_sun_grams
-    Mpc_cgs = cm_per_mpc
+    Msun_mks = mass_sun_kg
+    Mpc_mks = m_per_mpc
 
     u1 = Unit("Msun * Mpc**-3")
-    u2 = Unit("g * cm**-3")
+    u2 = Unit("kg * m**-3")
     u3 = u1.get_base_equivalent()
 
     assert u2.expr == u3.expr
     assert u2 == u3
 
-    assert_allclose_units(u1.base_value, Msun_cgs / Mpc_cgs**3, 1e-12)
+    assert_allclose_units(u1.base_value, Msun_mks / Mpc_mks**3, 1e-12)
     assert u2.base_value == 1
     assert u3.base_value == 1
 
@@ -445,7 +446,7 @@ def test_base_equivalent():
     assert u3.dimensions == mass_density
 
     assert_allclose_units(_get_conversion_factor(u1, u3)[0],
-                          Msun_cgs / Mpc_cgs**3, 1e-12)
+                          Msun_mks / Mpc_mks**3, 1e-12)
 
 
 def test_temperature_offsets():
@@ -466,7 +467,7 @@ def test_latex_repr():
                  "\\rm{pc}/(1+z)")
 
     test_unit = Unit('Mpccm', registry=registry)
-    assert_almost_equal(test_unit.base_value, cm_per_mpc/3)
+    assert_almost_equal(test_unit.base_value, m_per_mpc/3)
     assert_equal(test_unit.latex_repr, r'\rm{Mpc}/(1+z)')
 
     test_unit = Unit('cm**-3', base_value=1.0, registry=registry)
@@ -515,15 +516,15 @@ def test_creation_from_ytarray():
     assert_equal(u1, Unit('esu'))
     assert_equal(u1, electrostatic_unit.units)
 
-    u2 = Unit(elementary_charge)
+    u2 = Unit(elementary_charge_cgs)
     assert_equal(str(u2), '4.8032056e-10*esu')
     assert_equal(u2, Unit('4.8032056e-10*esu'))
-    assert_equal(u1, elementary_charge.units)
+    assert_equal(u1, elementary_charge_cgs.units)
 
-    assert_equal((u1/u2).base_value, electrostatic_unit/elementary_charge)
+    assert_allclose((u1/u2).base_value, electrostatic_unit/elementary_charge_cgs)
 
     with pytest.raises(UnitParseError):
-        Unit([1, 2, 3]*elementary_charge)
+        Unit([1, 2, 3]*elementary_charge_cgs)
 
 
 def test_list_same_dimensions():
diff --git a/unyt/tests/test_unyt_array.py b/unyt/tests/test_unyt_array.py
index e25977e..8cecb89 100644
--- a/unyt/tests/test_unyt_array.py
+++ b/unyt/tests/test_unyt_array.py
@@ -1236,7 +1236,7 @@ def test_equivalencies():
 
     lam = 4000*u.angstrom
     hnu = lam.in_units('erg', 'spectral')
-    assert_allclose_units(hnu, u.hcgs*u.clight/lam)
+    assert_allclose_units(hnu, u.hmks*u.clight/lam)
     lam.convert_to_units('eV', 'spectral')
     assert_allclose_units(lam, hnu)
     assert lam.units == u.eV.units
@@ -1256,7 +1256,7 @@ def test_equivalencies():
 
     nu = 1*u.MHz
     E = nu.to('erg', 'spectral')
-    assert_allclose_units(E, u.hcgs*nu)
+    assert_allclose_units(E, u.hmks*nu)
     nu.convert_to_units('J', 'spectral')
     assert_allclose_units(nu, E)
     assert nu.units == u.J.units
@@ -1266,7 +1266,7 @@ def test_equivalencies():
 
     E = 13.6*u.eV
     nu = E.to('Hz', 'spectral')
-    assert_allclose_units(nu, E/u.hcgs)
+    assert_allclose_units(nu, E/u.hmks)
     E.convert_to_units('MHz', 'spectral')
     assert_allclose_units(nu, E)
     assert E.units == u.MHz.units
@@ -1276,7 +1276,7 @@ def test_equivalencies():
 
     E = 13.6*u.eV
     lam = E.to('nm', 'spectral')
-    assert_allclose_units(lam, u.hcgs*u.clight/E)
+    assert_allclose_units(lam, u.hmks*u.clight/E)
     E.convert_to_units('angstrom', 'spectral')
     assert_allclose_units(E, lam)
     assert E.units == u.angstrom.units
@@ -1290,7 +1290,7 @@ def test_equivalencies():
     gg = 5./3.
     T = 1e8*u.K
     c_s = T.in_units("km/s", equivalence="sound_speed")
-    assert_equal(c_s, np.sqrt(gg*u.kboltz*T/(mu*u.mh)))
+    assert_allclose_units(c_s, np.sqrt(gg*u.kboltz*T/(mu*u.mh)))
     assert_allclose_units(T, c_s.in_units("K", "sound_speed"))
     T.convert_to_units('m/s', 'sound_speed')
     assert_allclose_units(c_s, T)
@@ -1301,7 +1301,7 @@ def test_equivalencies():
     gg = 4./3.
     T = 1e8*u.K
     c_s = T.in_units("km/s", "sound_speed", mu=mu, gamma=gg)
-    assert_equal(c_s, np.sqrt(gg*u.kboltz*T/(mu*u.mh)))
+    assert_allclose_units(c_s, np.sqrt(gg*u.kboltz*T/(mu*u.mh)))
     assert_allclose_units(T, c_s.in_units("K", "sound_speed", mu=mu, gamma=gg))
     T.convert_to_units('m/s', 'sound_speed', mu=mu, gamma=gg)
     assert_allclose_units(c_s, T)
@@ -1354,12 +1354,12 @@ def test_equivalencies():
 
     msun = 1*u.Msun
     msun.convert_to_equivalent('km', 'schwarzschild')
-    R = u.mass_sun_cgs.in_units("kpc", "schwarzschild")
+    R = u.mass_sun_mks.in_units("kpc", "schwarzschild")
     assert_allclose_units(msun, R)
-    assert_allclose_units(R.in_cgs(), 2*u.G*u.mass_sun_cgs/(u.clight**2))
-    assert_allclose_units(u.mass_sun_cgs, R.in_units("g", "schwarzschild"))
+    assert_allclose_units(R.in_mks(), 2*u.G*u.mass_sun_mks/(u.clight**2))
+    assert_allclose_units(u.mass_sun_mks, R.in_units("kg", "schwarzschild"))
     R.convert_to_units('Msun', 'schwarzschild')
-    assert_allclose_units(u.mass_sun_cgs, R)
+    assert_allclose_units(u.mass_sun_mks, R)
     assert R.units == u.Msun.units
     assert msun.units == u.km.units
 
@@ -1369,7 +1369,7 @@ def test_equivalencies():
     me.convert_to_units('nm', 'compton')
     length = u.me.in_units("angstrom", "compton")
     assert_allclose_units(length, me)
-    assert_allclose_units(length, u.hcgs/(u.me*u.clight))
+    assert_allclose_units(length, u.hmks/(u.me*u.clight))
     assert_allclose_units(u.me, length.in_units("g", "compton"))
     assert me.units == u.nm.units
     assert length.units == u.angstrom.units
@@ -1378,42 +1378,42 @@ def test_equivalencies():
 
     # Number density
 
-    rho = u.mp/u.cm**3
-    n = rho.in_units("cm**-3", "number_density")
+    rho = u.mp/u.m**3
+    n = rho.in_units("m**-3", "number_density")
     assert_allclose_units(n, rho/(u.mh*0.6))
-    assert_allclose_units(rho, n.in_units("g/cm**3", "number_density"))
-    rho.convert_to_units('m**-3', 'number_density')
-    assert rho.units == (1/u.m**3).units
-    assert n.units == (1/u.cm**3).units
+    assert_allclose_units(rho, n.in_units("kg/m**3", "number_density"))
+    rho.convert_to_units('cm**-3', 'number_density')
+    assert rho.units == (1/u.cm**3).units
+    assert n.units == (1/u.m**3).units
     assert_allclose_units(n, rho)
-    rho.convert_to_units('g/cm**3', 'number_density')
-    assert_allclose_units(u.mp/u.cm**3, rho)
-    assert rho.units == (u.g/u.cm**3).units
+    rho.convert_to_units('kg/m**3', 'number_density')
+    assert_allclose_units(u.mp/u.m**3, rho)
+    assert rho.units == (u.kg/u.m**3).units
 
-    rho = u.mp/u.cm**3
-    n = rho.in_units("cm**-3", equivalence="number_density", mu=0.75)
+    rho = u.mp/u.m**3
+    n = rho.in_units("m**-3", equivalence="number_density", mu=0.75)
     assert_allclose_units(n, rho/(u.mh*0.75))
     assert_allclose_units(
-        rho, n.in_units("g/cm**3", equivalence="number_density", mu=0.75))
-    rho.convert_to_units('m**-3', 'number_density', mu=0.75)
-    assert rho.units == (1/u.m**3).units
-    assert n.units == (1/u.cm**3).units
+        rho, n.in_units("kg/m**3", equivalence="number_density", mu=0.75))
+    rho.convert_to_units('cm**-3', 'number_density', mu=0.75)
+    assert rho.units == (1/u.cm**3).units
+    assert n.units == (1/u.m**3).units
     assert_allclose_units(n, rho)
-    rho.convert_to_units('g/cm**3', 'number_density', mu=0.75)
-    assert_allclose_units(u.mp/u.cm**3, rho)
-    assert rho.units == (u.g/u.cm**3).units
+    rho.convert_to_units('kg/m**3', 'number_density', mu=0.75)
+    assert_allclose_units(u.mp/u.m**3, rho)
+    assert rho.units == (u.kg/u.m**3).units
 
     # Effective temperature
 
     T = 1e4*u.K
-    F = T.in_units("erg/s/cm**2", equivalence="effective_temperature")
-    assert_allclose_units(F, u.stefan_boltzmann_constant_cgs*T**4)
+    F = T.in_units("W/m**2", equivalence="effective_temperature")
+    assert_allclose_units(F, u.stefan_boltzmann_constant*T**4)
     assert_allclose_units(
         T, F.in_units("K", equivalence="effective_temperature"))
-    T.convert_to_units('W/m**2', 'effective_temperature')
+    T.convert_to_units('erg/s/cm**2', 'effective_temperature')
     assert_allclose_units(T, F)
-    assert T.units == (u.W/u.m**2).units
-    assert F.units == (u.erg/u.s/u.cm**2).units
+    assert T.units == unyt_quantity(1.0, "erg/cm**2/s").units
+    assert F.units == u.W/u.m**2
     T.convert_to_units('K', 'effective_temperature')
     assert_almost_equal(T.value, 1e4)
     assert T.units == u.K.units
@@ -1422,10 +1422,10 @@ def test_equivalencies():
 
     assert_allclose_units(
         F.value,
-        T.to_value("erg/s/cm**2", equivalence="effective_temperature"))
+        T.to_value("W/m**2", equivalence="effective_temperature"))
     assert_allclose_units(
         n.value,
-        rho.to_value("cm**-3", equivalence="number_density", mu=0.75))
+        rho.to_value("m**-3", equivalence="number_density", mu=0.75))
 
 
 def test_electromagnetic():
@@ -1433,29 +1433,29 @@ def test_electromagnetic():
 
     # Various tests of SI and CGS electromagnetic units
 
-    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 = 1*u.qp
+    qp_cgs = u.qp_cgs.in_units("C", "SI")
+    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")
     assert_equal(qp.units.dimensions, dimensions.charge_mks)
     assert_almost_equal(qp.v, 10*u.qp.v/u.clight.v)
 
-    qp_cgs = qp_mks.in_units("esu", "CGS")
-    assert_array_almost_equal(qp_cgs, u.qp)
-    assert_equal(qp_cgs.units.dimensions, u.qp.units.dimensions)
-    qp = qp_mks.copy()
+    qp_cgs = u.qp.in_units("esu", "CGS")
+    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')
-    assert_almost_equal(qp_cgs, qp)
+    assert_almost_equal(qp_cgs, qp_cgs)
     assert qp.units == u.esu.units
     qp.convert_to_units('C', 'MKS')
-    assert_almost_equal(qp_mks, qp)
+    assert_almost_equal(u.qp, qp)
     assert qp.units == u.C.units
 
-    qp_mks_k = u.qp.in_units("kC", "SI")
+    qp_mks_k = u.qp_cgs.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 = 1*u.qp
+        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')
     assert_almost_equal(qp, qp_mks_k)
 
@@ -1464,8 +1464,13 @@ def test_electromagnetic():
     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.convert_to_cgs(equivalence="cgs")
+    B.convert_to_cgs()
     assert_almost_equal(B, B_cgs)
+    B.convert_to_mks()
+    B_cgs2 = B.to("gauss")
+    assert_almost_equal(B_cgs, B_cgs2)
+    B_mks2 = B_cgs2.to("T")
+    assert_almost_equal(B, B_mks2)
 
     B = 1.0*u.T
     u_mks = B*B/(2*u.mu_0)
@@ -1485,6 +1490,10 @@ def test_electromagnetic():
     assert current.units == u.statA.units
     current.convert_to_units('A', equivalence='MKS')
     assert current.units == u.A.units
+    I_cgs2 = current.to('statA')
+    assert I_cgs2.units == u.statA.units
+    assert_array_almost_equal(
+        I_cgs2, unyt_quantity(0.1*speed_of_light_cm_per_s, "statA"))
 
     current = 1.0*u.A
     R = unyt_quantity(1.0, "ohm")
@@ -1544,15 +1553,13 @@ def test_modified_unit_division():
     reg1 = UnitRegistry()
     reg2 = UnitRegistry()
 
-    # this mocks comoving coordinates without going through the trouble
-    # of setting up a fake cosmological dataset
-    reg1.modify('m', 50)
+    reg1.modify('g', 50)
 
-    a = unyt_quantity(3, 'm', registry=reg1)
-    b = unyt_quantity(3, 'm', registry=reg2)
+    a = unyt_quantity(3, 'g', registry=reg1)
+    b = unyt_quantity(3, 'g', registry=reg2)
 
     ret = a/b
-    assert ret == 0.5
+    assert ret == 50000.
     assert ret.units.is_dimensionless
     assert ret.units.base_value == 1.0
 
@@ -1602,8 +1609,8 @@ def test_initialization_different_registries():
     l1 = unyt_quantity(1.0, 'code_length', registry=reg1)
     l2 = unyt_quantity(1.0, 'code_length', registry=reg2)
 
-    assert_almost_equal(float(l1.in_cgs()), 1.0)
-    assert_almost_equal(float(l2.in_cgs()), 3.0)
+    assert_almost_equal(float(l1.in_mks()), 1.0)
+    assert_almost_equal(float(l2.in_mks()), 3.0)
 
 
 def test_ones_and_zeros_like():