diff --git a/moldesign/_tests/conftest.py b/moldesign/_tests/conftest.py
new file mode 100644
index 0000000..108f711
--- /dev/null
+++ b/moldesign/_tests/conftest.py
@@ -0,0 +1,14 @@
+
+
+def pytest_itemcollected(item):
+    if hasattr(item.module, '__PYTEST_MARK__'):
+        item.add_marker(item.module.__PYTEST_MARK__)
+
+
+# TODO: nicer output strings for git commit status
+# see https://docs.pytest.org/en/latest/example/simple.html#post-process-test-reports-failures
+#@pytest.hookimpl(tryfirst=True, hookwrapper=True)
+#def pytest_runtest_makereport(item, call):i
+#    pass
+# Also possibly useful:  item.add_report_section
+
diff --git a/moldesign/_tests/helpers.py b/moldesign/_tests/helpers.py
index d634479..e8bb3ea 100644
--- a/moldesign/_tests/helpers.py
+++ b/moldesign/_tests/helpers.py
@@ -53,24 +53,6 @@ def num_grad(mol, fn, step=DEFSTEP, atoms=None, fnargs=None, fnkwargs=None):
 
 ENERGY_TOLERANCE = 1e-8 * u.hartree
 
-def minimization_tester(mol):
-    assert 'potential_energy' not in mol.properties
-
-    e1 = mol.calculate_potential_energy()
-    p1 = mol.positions.copy()
-
-    traj = mol.minimize()
-
-    # check trajectory has correct initial and final positions, energies
-    assert (p1 == traj.frames[0].positions).all()
-    assert (mol.positions == traj.frames[-1].positions).all()
-    assert abs(e1 - traj.frames[0].potential_energy) < ENERGY_TOLERANCE
-    assert abs(mol.potential_energy - traj.frames[-1].potential_energy) < ENERGY_TOLERANCE
-
-    # Force recalculation of energy to check that it's correct
-    mol.calculate(use_cache=False)
-    assert abs(mol.potential_energy - traj.frames[-1].potential_energy) < ENERGY_TOLERANCE
-
 
 def _make_mol_with_n_hydrogens(n):
     return mdt.Molecule([mdt.Atom('H') for i in range(n)])
@@ -94,7 +76,6 @@ def calculate(self):
                     forces=np.zeros(self.mol.positions.shape)*u.default.force)
 
 
-
 def _internet(host="8.8.8.8", port=53, timeout=3):
     """
     Host: 8.8.8.8 (google-public-dns-a.google.com)
@@ -129,21 +110,42 @@ def assert_something_resembling_minimization_happened(p0, e0, traj, mol):
     Returns:
 
     """
+    # Use assert_almost_equal to account for small numerical differences in
+    # convergence noise, serialization methods, etc.
+
     import scipy.optimize.optimize
 
     assert traj.num_frames > 1
 
-    assert traj.potential_energy[0] == e0
+    assert_almost_equal(traj.potential_energy[0],
+                        e0, decimal=9)
+    assert_almost_equal(traj.potential_energy[-1],
+                        mol.potential_energy, decimal=9)
     assert traj.potential_energy[-1] < e0
-    assert traj.potential_energy[-1] == mol.potential_energy
 
     assert (traj.positions[0] == p0).all()
     assert (traj.positions[-1] != p0).any()
     assert (traj.positions[-1] == mol.positions).all()
 
+    recalc = mol.calculate_potential_energy(usecache=False)
+    assert (recalc - traj.potential_energy[-1]) < 1.0e-8 * u.hartree
+
     scipyresult = getattr(traj, 'info', None)
     if isinstance(scipyresult, scipy.optimize.optimize.OptimizeResult):
         np.testing.assert_allclose(scipyresult.x,
                                    mol.positions.defunits_value().flat)
-        np.testing.assert_almost_equal(scipyresult.fun,
-                                       mol.potential_energy.defunits_value())
+        assert_almost_equal(scipyresult.fun,
+                            mol.potential_energy.defunits_value())
+
+
+
+def assert_almost_equal(actual, desired, **kwargs):
+    units = mdt.units.get_units(actual)
+    desired = units.value_of(desired)
+
+    np.testing.assert_almost_equal(actual.value_in(units),
+                                   desired,
+                                   **kwargs)
+
+
+
diff --git a/moldesign/_tests/test_atoms_and_bonds.py b/moldesign/_tests/test_atom_bond_computed_properties.py
similarity index 100%
rename from moldesign/_tests/test_atoms_and_bonds.py
rename to moldesign/_tests/test_atom_bond_computed_properties.py
diff --git a/moldesign/_tests/test_atom_containers.py b/moldesign/_tests/test_atom_containers.py
index c1131da..c1084fb 100644
--- a/moldesign/_tests/test_atom_containers.py
+++ b/moldesign/_tests/test_atom_containers.py
@@ -11,6 +11,7 @@
 from . import helpers
 from .molecule_fixtures import pdb3aid, ethylene_waterbox_2na_2cl
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
 
 registered_types = {}
 
diff --git a/moldesign/_tests/test_cli.py b/moldesign/_tests/test_cli.py
index 6c38dcb..2d85539 100644
--- a/moldesign/_tests/test_cli.py
+++ b/moldesign/_tests/test_cli.py
@@ -5,6 +5,9 @@
 from moldesign.external import pathlib
 
 
+__PYTEST_MARK__ = 'io'
+
+
 @pytest.fixture
 def example_path(tmpdir):
     path = pathlib.Path(str(tmpdir))
diff --git a/moldesign/_tests/test_copies.py b/moldesign/_tests/test_copies.py
index 8b4b085..4a68822 100644
--- a/moldesign/_tests/test_copies.py
+++ b/moldesign/_tests/test_copies.py
@@ -9,6 +9,9 @@
 from .object_fixtures import *
 
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 def test_carbon_copy(carbon_copy, carbon_atom):
     atom = carbon_copy
     assert atom.symbol == carbon_atom.symbol
diff --git a/moldesign/_tests/test_dna_primary_structure.py b/moldesign/_tests/test_dna_primary_structure.py
index beaabab..51e8c82 100644
--- a/moldesign/_tests/test_dna_primary_structure.py
+++ b/moldesign/_tests/test_dna_primary_structure.py
@@ -7,6 +7,9 @@
 from moldesign import units as u
 
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 fixture_types = {}
 
 
diff --git a/moldesign/_tests/test_geometry.py b/moldesign/_tests/test_geometry.py
index 714cfca..d19b7e9 100644
--- a/moldesign/_tests/test_geometry.py
+++ b/moldesign/_tests/test_geometry.py
@@ -14,6 +14,10 @@
 
 registered_types = {}
 
+
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 # TODO: automated method testing based on its metadata - i.e. test to make sure parameters are
 #       honored, test that it calcultes what it says it does, test that properties have the right
 #       units and array shapes, etc.
diff --git a/moldesign/_tests/test_imports.py b/moldesign/_tests/test_imports.py
index 77082d8..636c58b 100644
--- a/moldesign/_tests/test_imports.py
+++ b/moldesign/_tests/test_imports.py
@@ -2,6 +2,9 @@
 import sys
 
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 @pytest.mark.tryfirst
 def test_lazy_imports():
     import moldesign
diff --git a/moldesign/_tests/test_io.py b/moldesign/_tests/test_io.py
index 50974d0..c09782a 100644
--- a/moldesign/_tests/test_io.py
+++ b/moldesign/_tests/test_io.py
@@ -23,6 +23,9 @@
 from .object_fixtures import h2_trajectory, h2_harmonic, h2
 
 
+__PYTEST_MARK__ = 'io'
+
+
 @pytest.fixture
 def bipyridine_sdf():
     return mdt.read(get_data_path('bipyridine.sdf'))
diff --git a/moldesign/_tests/test_molecules.py b/moldesign/_tests/test_molecules.py
index b2aa866..af90080 100644
--- a/moldesign/_tests/test_molecules.py
+++ b/moldesign/_tests/test_molecules.py
@@ -12,6 +12,9 @@
 from .test_qm_xfaces import h2_rhfwfn
 
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 def test_h2_protected_atom_arrays(h2):
     atom1, atom2 = h2.atoms
     with pytest.raises(TypeError):
diff --git a/moldesign/_tests/test_objects.py b/moldesign/_tests/test_objects.py
index b48bfa5..e242aee 100644
--- a/moldesign/_tests/test_objects.py
+++ b/moldesign/_tests/test_objects.py
@@ -4,6 +4,9 @@
 from .object_fixtures import *
 
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 class ComposedClass(object):
     delegated = Alias('s.lower')
 
diff --git a/moldesign/_tests/test_qm_xfaces.py b/moldesign/_tests/test_qm_xfaces.py
index 28f7f66..de20d75 100644
--- a/moldesign/_tests/test_qm_xfaces.py
+++ b/moldesign/_tests/test_qm_xfaces.py
@@ -15,6 +15,7 @@
 #       honored, test that it calcultes what it says it does, test that properties have the right
 #       units and array shapes, etc.
 
+
 def typedfixture(*types, **kwargs):
     """This is a decorator that lets us associate fixtures with one or more arbitrary types.
     We'll later use this type to determine what tests to run on the result"""
@@ -54,12 +55,26 @@ def heh_plus():
 @pytest.fixture(params=models_to_test, ids=model_ids, scope='function')
 def h2_with_model(request, h2):
     model, basis, theory = request.param
+
+    if model is mdt.models.NWChemQM and theory == 'mp2':
+        pytest.xfail('Not implemented')
+
     h2.set_energy_model(model, basis=basis, theory=theory)
     return h2
 
 
 def test_minimization_trajectory(h2_with_model):
-    helpers.minimization_tester(h2_with_model)
+    mol = h2_with_model
+    if mol.energy_model.params.theory == 'mp2':
+        pytest.skip('Not testing mp2 minimizations at this time')
+
+    assert 'potential_energy' not in mol.properties
+
+    e1 = mol.calculate_potential_energy()
+    p1 = mol.positions.copy()
+
+    traj = mol.minimize()
+    helpers.assert_something_resembling_minimization_happened(p1, e1, traj, mol)
 
 
 @typedfixture('model')
diff --git a/moldesign/_tests/test_trajectory.py b/moldesign/_tests/test_trajectory.py
index e0a616a..fb99065 100644
--- a/moldesign/_tests/test_trajectory.py
+++ b/moldesign/_tests/test_trajectory.py
@@ -6,6 +6,8 @@
 
 from .object_fixtures import h2, h2_harmonic, h2_trajectory
 
+
+@pytest.mark.internal
 def test_frames_synched_with_trajectory(h2_trajectory):
     traj = h2_trajectory
 
@@ -44,6 +46,7 @@ def precanned_trajectory():
     return traj
 
 
+@pytest.mark.internal
 def test_geometry_analysis_precanned(precanned_trajectory):
     traj = precanned_trajectory
     a1, a2, a3 = traj.mol.atoms
@@ -72,6 +75,7 @@ def test_geometry_analysis_precanned(precanned_trajectory):
     assert traj.someletter == list('abc')
 
 
+@pytest.mark.internal
 def test_frame_to_molecule_conversion(precanned_trajectory):
     traj = precanned_trajectory
 
@@ -101,6 +105,7 @@ def test_frame_to_molecule_conversion(precanned_trajectory):
     assert m2.time == 2.0 * u.fs
 
 
+@pytest.mark.internal
 def test_property_backfill(precanned_trajectory):
     traj = precanned_trajectory
     oldnumframes = len(traj)
@@ -111,6 +116,7 @@ def test_property_backfill(precanned_trajectory):
     assert traj.somenewthing == [None] * oldnumframes + [5]
 
 
+@pytest.mark.internal
 def test_add_traj(precanned_trajectory):
     newtraj = precanned_trajectory + precanned_trajectory
 
diff --git a/moldesign/_tests/test_units.py b/moldesign/_tests/test_units.py
index 99ad49a..465d0b6 100644
--- a/moldesign/_tests/test_units.py
+++ b/moldesign/_tests/test_units.py
@@ -6,6 +6,9 @@
 from moldesign import units
 
 
+__PYTEST_MARK__ = 'internal'  # mark all tests in this module with this label (see ./conftest.py)
+
+
 def test_scalar_comparison_dimensionality_errors():
     with pytest.raises(units.DimensionalityError):
         x = 1.0 * units.angstrom == 1.0*units.ureg.kilograms
@@ -123,9 +126,17 @@ def test_array_unit_checks():
     np.testing.assert_allclose(arr.magnitude,
                                np.arange(10, 15))
 
+
 def test_default_unit_conversions():
     assert abs(10.0 - (1.0*units.nm).defunits_value()) < 1e-10
     assert abs(1000.0 - (1.0*units.ps).defunits_value()) < 1e-10
     assert abs(1.0 - 6.022140857e23/((1.0*units.ureg.grams).defunits_value())) < 1e-6
     assert abs(103.642685491 - (1.0*units.angstrom**2*units.dalton/units.fs**2).defunits_value()
                ) < 1e-7
+
+
+def test_getunits_doctests():
+    assert units.get_units(1.0*units.angstrom) == units.MdtUnit('ang')
+    assert units.get_units(np.array([1.0, 2, 3.0])) == units.MdtUnit('dimensionless')
+    assert units.get_units([[1.0*units.dalton, 3.0*units.eV],
+                            ['a'], 'gorilla']) == units.MdtUnit('amu')
diff --git a/pytest.ini b/pytest.ini
index c6846c3..37e738b 100644
--- a/pytest.ini
+++ b/pytest.ini
@@ -2,5 +2,7 @@
 testpaths = moldesign/_tests
 
 markers =
-   base: tests data structures, molecule construction
-   slow: heavy duty simulation and numerical tests
+   internal: tests MDT data structures, consistency of molecular data
+   io: just molecular reading and writing from various formats
+   slow: tests simulations, numerical computatoins, and external interfaces
+