diff --git a/interfaces/cython/cantera/base.pyx b/interfaces/cython/cantera/base.pyx
index 735dccde2..f6397b0f9 100644
--- a/interfaces/cython/cantera/base.pyx
+++ b/interfaces/cython/cantera/base.pyx
@@ -120,7 +120,7 @@ cdef class _SolutionBase:
         def __get__(self):
             return list(self._selected_species)
         def __set__(self, species):
-            if isinstance(species, (str, unicode, int)):
+            if isinstance(species, (str, int)):
                 species = (species,)
             self._selected_species.resize(len(species))
             for i,spec in enumerate(species):
diff --git a/interfaces/cython/cantera/composite.py b/interfaces/cython/cantera/composite.py
index 6d11e765e..938a3eb6c 100644
--- a/interfaces/cython/cantera/composite.py
+++ b/interfaces/cython/cantera/composite.py
@@ -93,7 +93,7 @@ class DustyGas(ThermoPhase, Kinetics, DustyGasTransport):
     __slots__ = ()
 
 
-class Quantity(object):
+class Quantity:
     """
     A class representing a specific quantity of a `Solution`. In addition to the
     properties which can be computed for class `Solution`, class `Quantity`
@@ -286,7 +286,7 @@ for _attr in dir(Solution):
         setattr(Quantity, _attr, _prop(_attr))
 
 
-class SolutionArray(object):
+class SolutionArray:
     """
     A class providing a convenient interface for representing many thermodynamic
     states using the same `Solution` object and computing properties for that
diff --git a/interfaces/cython/cantera/examples/onedim/flamespeed_sensitivity.py b/interfaces/cython/cantera/examples/onedim/flamespeed_sensitivity.py
index f5ebae78b..032b9d54c 100644
--- a/interfaces/cython/cantera/examples/onedim/flamespeed_sensitivity.py
+++ b/interfaces/cython/cantera/examples/onedim/flamespeed_sensitivity.py
@@ -4,8 +4,6 @@ flame. Computes the sensitivity of the laminar flame speed with respect
 to each reaction rate constant.
 """
 
-from __future__ import print_function
-
 import cantera as ct
 import numpy as np
 
diff --git a/interfaces/cython/cantera/examples/reactors/NonIdealShockTube.py b/interfaces/cython/cantera/examples/reactors/NonIdealShockTube.py
index 996a55332..06448b369 100644
--- a/interfaces/cython/cantera/examples/reactors/NonIdealShockTube.py
+++ b/interfaces/cython/cantera/examples/reactors/NonIdealShockTube.py
@@ -19,9 +19,6 @@ Other than the typical Cantera dependencies, plotting functions require that
 you have matplotlib (https://matplotlib.org/) installed.
 """
 
-from __future__ import division
-from __future__ import print_function
-
 # Dependencies: numpy, and matplotlib
 import numpy as np
 import matplotlib.pyplot as plt
diff --git a/interfaces/cython/cantera/examples/reactors/custom.py b/interfaces/cython/cantera/examples/reactors/custom.py
index 10e24bda5..7c982f28e 100644
--- a/interfaces/cython/cantera/examples/reactors/custom.py
+++ b/interfaces/cython/cantera/examples/reactors/custom.py
@@ -14,7 +14,7 @@ import cantera as ct
 import numpy as np
 import scipy.integrate
 
-class ReactorOde(object):
+class ReactorOde:
     def __init__(self, gas):
         # Parameters of the ODE system and auxiliary data are stored in the
         # ReactorOde object.
diff --git a/interfaces/cython/cantera/func1.pyx b/interfaces/cython/cantera/func1.pyx
index 07202130b..44dfa74ba 100644
--- a/interfaces/cython/cantera/func1.pyx
+++ b/interfaces/cython/cantera/func1.pyx
@@ -37,7 +37,7 @@ cdef class Func1:
         >>> f2(3)
         10
 
-        >>> class Multiplier(object):
+        >>> class Multiplier:
         ...     def __init__(self, factor):
         ...         self.factor = factor
         ...     def __call__(self, t):
diff --git a/interfaces/cython/cantera/kinetics.pyx b/interfaces/cython/cantera/kinetics.pyx
index 271369918..65aa7794d 100644
--- a/interfaces/cython/cantera/kinetics.pyx
+++ b/interfaces/cython/cantera/kinetics.pyx
@@ -68,7 +68,7 @@ cdef class Kinetics(_SolutionBase):
         argument is unused.
         """
         cdef int k
-        if isinstance(species, (str, unicode, bytes)):
+        if isinstance(species, (str, bytes)):
             return self.kinetics.kineticsSpeciesIndex(stringify(species))
         else:
             k = species
@@ -178,7 +178,7 @@ cdef class Kinetics(_SolutionBase):
         reaction *i_reaction*.
         """
         cdef int k
-        if isinstance(k_spec, (str, unicode, bytes)):
+        if isinstance(k_spec, (str, bytes)):
             k = self.kinetics_species_index(k_spec)
         else:
             k = k_spec
@@ -193,7 +193,7 @@ cdef class Kinetics(_SolutionBase):
         reaction *i_reaction*.
         """
         cdef int k
-        if isinstance(k_spec, (str, unicode, bytes)):
+        if isinstance(k_spec, (str, bytes)):
             k = self.kinetics_species_index(k_spec)
         else:
             k = k_spec
diff --git a/interfaces/cython/cantera/mixture.pyx b/interfaces/cython/cantera/mixture.pyx
index e4e52aa3d..bcc20c5b2 100644
--- a/interfaces/cython/cantera/mixture.pyx
+++ b/interfaces/cython/cantera/mixture.pyx
@@ -4,7 +4,7 @@
 import warnings
 
 # Need a pure-python class to store weakrefs to
-class _WeakrefProxy(object):
+class _WeakrefProxy:
     pass
 
 cdef class Mixture:
@@ -91,7 +91,7 @@ cdef class Mixture:
             >>> mix.element_index('H')
             2
         """
-        if isinstance(element, (str, unicode, bytes)):
+        if isinstance(element, (str, bytes)):
             index = self.mix.elementIndex(stringify(element))
         elif isinstance(element, (int, float)):
             index = <int>element
@@ -128,7 +128,7 @@ cdef class Mixture:
         """
         p = self.phase_index(phase)
 
-        if isinstance(species, (str, unicode, bytes)):
+        if isinstance(species, (str, bytes)):
             k = self.phase(p).species_index(species)
         elif isinstance(species, (int, float)):
             k = <int?>species
@@ -169,7 +169,7 @@ cdef class Mixture:
                 return int(p)
             else:
                 raise IndexError("Phase index '{0}' out of range.".format(p))
-        elif isinstance(p, (str, unicode, bytes)):
+        elif isinstance(p, (str, bytes)):
             for i, phase in enumerate(self._phases):
                 if phase.name == p:
                     return i
@@ -258,7 +258,7 @@ cdef class Mixture:
             return data
 
         def __set__(self, moles):
-            if isinstance(moles, (str, unicode, bytes)):
+            if isinstance(moles, (str, bytes)):
                 self.mix.setMolesByName(stringify(moles))
                 return
 
diff --git a/interfaces/cython/cantera/onedim.py b/interfaces/cython/cantera/onedim.py
index daa35ee0a..a46ba6bd0 100644
--- a/interfaces/cython/cantera/onedim.py
+++ b/interfaces/cython/cantera/onedim.py
@@ -5,12 +5,7 @@ import numpy as np
 from ._cantera import *
 from .composite import Solution
 import csv as _csv
-
-try:
-    # Python 2.7 or 3.2+
-    from math import erf
-except ImportError:
-    from scipy.special import erf
+from math import erf
 
 
 class FlameBase(Sim1D):
@@ -27,7 +22,7 @@ class FlameBase(Sim1D):
         if grid is None:
             grid = np.linspace(0.0, 0.1, 6)
         self.flame.grid = grid
-        super(FlameBase, self).__init__(domains)
+        super().__init__(domains)
         self.gas = gas
         self.flame.P = gas.P
 
@@ -54,8 +49,7 @@ class FlameBase(Sim1D):
 
         >>> f.set_refine_criteria(ratio=3.0, slope=0.1, curve=0.2, prune=0)
         """
-        super(FlameBase, self).set_refine_criteria(self.flame, ratio, slope,
-                                                   curve, prune)
+        super().set_refine_criteria(self.flame, ratio, slope, curve, prune)
 
     def get_refine_criteria(self):
         """
@@ -67,7 +61,7 @@ class FlameBase(Sim1D):
         >>> f.get_refine_criteria()
         {'ratio': 3.0, 'slope': 0.1, 'curve': 0.2, 'prune': 0.0}
         """
-        return super(FlameBase, self).get_refine_criteria(self.flame)
+        return super().get_refine_criteria(self.flame)
 
     def set_profile(self, component, locations, values):
         """
@@ -82,8 +76,7 @@ class FlameBase(Sim1D):
 
         >>> f.set_profile('T', [0.0, 0.2, 1.0], [400.0, 800.0, 1500.0])
         """
-        super(FlameBase, self).set_profile(self.flame, component, locations,
-                                           values)
+        super().set_profile(self.flame, component, locations, values)
 
     @property
     def max_grid_points(self):
@@ -91,11 +84,11 @@ class FlameBase(Sim1D):
         Get/Set the maximum number of grid points used in the solution of
         this flame.
         """
-        return super(FlameBase, self).get_max_grid_points(self.flame)
+        return super().get_max_grid_points(self.flame)
 
     @max_grid_points.setter
     def max_grid_points(self, npmax):
-        super(FlameBase, self).set_max_grid_points(self.flame, npmax)
+        super().set_max_grid_points(self.flame, npmax)
 
     @property
     def transport_model(self):
@@ -416,8 +409,7 @@ class FreeFlame(FlameBase):
         if width is not None:
             grid = np.array([0.0, 0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1.0]) * width
 
-        super(FreeFlame, self).__init__((self.inlet, self.flame, self.outlet),
-                                        gas, grid)
+        super().__init__((self.inlet, self.flame, self.outlet), gas, grid)
 
         # Setting X needs to be deferred until linked to the flow domain
         self.inlet.T = gas.T
@@ -434,7 +426,7 @@ class FreeFlame(FlameBase):
             Profiles rise linearly between the second and third location.
             Locations are given as a fraction of the entire domain
         """
-        super(FreeFlame, self).set_initial_guess()
+        super().set_initial_guess()
         self.gas.TPY = self.inlet.T, self.P, self.inlet.Y
 
         if not self.inlet.mdot:
@@ -489,7 +481,7 @@ class FreeFlame(FlameBase):
             will be calculated.
         """
         if not auto:
-            return super(FreeFlame, self).solve(loglevel, refine_grid, auto)
+            return super().solve(loglevel, refine_grid, auto)
 
         # Use a callback function to check that the domain is actually wide
         # enough to contain the flame after each steady-state solve. If the user
@@ -521,7 +513,7 @@ class FreeFlame(FlameBase):
 
         for _ in range(12):
             try:
-                return super(FreeFlame, self).solve(loglevel, refine_grid, auto)
+                return super().solve(loglevel, refine_grid, auto)
             except DomainTooNarrow:
                 self.flame.grid *= 2
                 if loglevel > 0:
@@ -611,10 +603,10 @@ class IonFlameBase(FlameBase):
     def solve(self, loglevel=1, refine_grid=True, auto=False, stage=1, enable_energy=True):
         self.flame.set_solving_stage(stage)
         if stage == 1:
-            super(IonFlameBase, self).solve(loglevel, refine_grid, auto)
+            super().solve(loglevel, refine_grid, auto)
         if stage == 2:
             self.poisson_enabled = True
-            super(IonFlameBase, self).solve(loglevel, refine_grid, auto)
+            super().solve(loglevel, refine_grid, auto)
 
 
 class IonFreeFlame(IonFlameBase, FreeFlame):
@@ -627,7 +619,7 @@ class IonFreeFlame(IonFlameBase, FreeFlame):
             self.flame = IonFlow(gas, name='flame')
             self.flame.set_free_flow()
 
-        super(IonFreeFlame, self).__init__(gas, grid, width)
+        super().__init__(gas, grid, width)
 
 
 class BurnerFlame(FlameBase):
@@ -662,8 +654,7 @@ class BurnerFlame(FlameBase):
         if width is not None:
             grid = np.array([0.0, 0.1, 0.2, 0.3, 0.5, 0.7, 1.0]) * width
 
-        super(BurnerFlame, self).__init__((self.burner, self.flame, self.outlet),
-                                          gas, grid)
+        super().__init__((self.burner, self.flame, self.outlet), gas, grid)
 
         # Setting X needs to be deferred until linked to the flow domain
         self.burner.T = gas.T
@@ -677,7 +668,7 @@ class BurnerFlame(FlameBase):
         20% of the flame to Tad, then is flat. The mass fraction profiles are
         set similarly.
         """
-        super(BurnerFlame, self).set_initial_guess()
+        super().set_initial_guess()
 
         self.gas.TPY = self.burner.T, self.P, self.burner.Y
         Y0 = self.burner.Y
@@ -740,7 +731,7 @@ class BurnerFlame(FlameBase):
             self.set_steady_callback(check_blowoff)
 
         try:
-            return super(BurnerFlame, self).solve(loglevel, refine_grid, auto)
+            return super().solve(loglevel, refine_grid, auto)
         except FlameBlowoff:
             # The eventual solution for a blown off flame is the non-reacting
             # solution, so just set the state to this now
@@ -749,7 +740,7 @@ class BurnerFlame(FlameBase):
                 self.set_flat_profile(self.flame, spec, self.burner.Y[k])
 
             self.set_steady_callback(original_callback)
-            super(BurnerFlame, self).solve(loglevel, False, False)
+            super().solve(loglevel, False, False)
             if loglevel > 0:
                 print('Flame has blown off of burner (non-reacting solution)')
 
@@ -766,7 +757,7 @@ class IonBurnerFlame(IonFlameBase, BurnerFlame):
             self.flame = IonFlow(gas, name='flame')
             self.flame.set_axisymmetric_flow()
 
-        super(IonBurnerFlame, self).__init__(gas, grid, width)
+        super().__init__(gas, grid, width)
 
 
 class CounterflowDiffusionFlame(FlameBase):
@@ -803,8 +794,7 @@ class CounterflowDiffusionFlame(FlameBase):
         if width is not None:
             grid = np.array([0.0, 0.2, 0.4, 0.6, 0.8, 1.0]) * width
 
-        super(CounterflowDiffusionFlame, self).__init__(
-                (self.fuel_inlet, self.flame, self.oxidizer_inlet), gas, grid)
+        super().__init__((self.fuel_inlet, self.flame, self.oxidizer_inlet), gas, grid)
 
     def set_initial_guess(self):
         """
@@ -812,7 +802,7 @@ class CounterflowDiffusionFlame(FlameBase):
         by assuming infinitely-fast chemistry.
         """
 
-        super(CounterflowDiffusionFlame, self).set_initial_guess()
+        super().set_initial_guess()
 
         moles = lambda el: (self.gas.elemental_mass_fraction(el) /
                             self.gas.atomic_weight(el))
@@ -907,7 +897,7 @@ class CounterflowDiffusionFlame(FlameBase):
             transport is enabled, an additional solution using these options
             will be calculated.
         """
-        super(CounterflowDiffusionFlame, self).solve(loglevel, refine_grid, auto)
+        super().solve(loglevel, refine_grid, auto)
         # Do some checks if loglevel is set
         if loglevel > 0:
             if self.extinct():
@@ -1087,8 +1077,7 @@ class ImpingingJet(FlameBase):
             self.surface.set_kinetics(surface)
             self.surface.T = surface.T
 
-        super(ImpingingJet, self).__init__(
-                (self.inlet, self.flame, self.surface), gas, grid)
+        super().__init__((self.inlet, self.flame, self.surface), gas, grid)
 
         # Setting X needs to be deferred until linked to the flow domain
         self.inlet.T = gas.T
@@ -1101,7 +1090,7 @@ class ImpingingJet(FlameBase):
         used to form the initial guess. Otherwise the inlet composition will
         be used.
         """
-        super(ImpingingJet, self).set_initial_guess(products=products)
+        super().set_initial_guess(products=products)
 
         Y0 = self.inlet.Y
         T0 = self.inlet.T
@@ -1163,8 +1152,7 @@ class CounterflowPremixedFlame(FlameBase):
             # Create grid points aligned with initial guess profile
             grid = np.array([0.0, 0.3, 0.5, 0.7, 1.0]) * width
 
-        super(CounterflowPremixedFlame, self).__init__(
-                (self.reactants, self.flame, self.products), gas, grid)
+        super().__init__((self.reactants, self.flame, self.products), gas, grid)
 
         # Setting X needs to be deferred until linked to the flow domain
         self.reactants.X = gas.X
@@ -1177,7 +1165,7 @@ class CounterflowPremixedFlame(FlameBase):
         will be set to the equilibrium state of the reactants mixture.
         """
 
-        super(CounterflowPremixedFlame, self).set_initial_guess()
+        super().set_initial_guess()
 
         Yu = self.reactants.Y
         Tu = self.reactants.T
@@ -1257,8 +1245,7 @@ class CounterflowTwinPremixedFlame(FlameBase):
             # Create grid points aligned with initial guess profile
             grid = np.array([0.0, 0.2, 0.4, 0.5, 0.6, 0.8, 1.0]) * width
 
-        super(CounterflowTwinPremixedFlame, self).__init__(
-                (self.reactants, self.flame, self.products), gas, grid)
+        super().__init__((self.reactants, self.flame, self.products), gas, grid)
 
         # Setting X needs to be deferred until linked to the flow domain
         self.reactants.X = gas.X
@@ -1267,7 +1254,7 @@ class CounterflowTwinPremixedFlame(FlameBase):
         """
         Set the initial guess for the solution.
         """
-        super(CounterflowTwinPremixedFlame, self).set_initial_guess()
+        super().set_initial_guess()
 
         Yu = self.reactants.Y
         Tu = self.reactants.T
diff --git a/interfaces/cython/cantera/onedim.pyx b/interfaces/cython/cantera/onedim.pyx
index 485d59f82..f9832e4a9 100644
--- a/interfaces/cython/cantera/onedim.pyx
+++ b/interfaces/cython/cantera/onedim.pyx
@@ -5,7 +5,7 @@ import interrupts
 import warnings
 
 # Need a pure-python class to store weakrefs to
-class _WeakrefProxy(object):
+class _WeakrefProxy:
     pass
 
 cdef class Domain1D:
@@ -695,7 +695,7 @@ cdef class Sim1D:
     def _get_indices(self, dom, comp):
         idom = self.domain_index(dom)
         dom = self.domains[idom]
-        if isinstance(comp, (str, unicode, bytes)):
+        if isinstance(comp, (str, bytes)):
             kcomp = dom.component_index(comp)
         else:
             kcomp = comp
diff --git a/interfaces/cython/cantera/reactor.pyx b/interfaces/cython/cantera/reactor.pyx
index eaf842e04..283983de5 100644
--- a/interfaces/cython/cantera/reactor.pyx
+++ b/interfaces/cython/cantera/reactor.pyx
@@ -7,7 +7,7 @@ import numbers as _numbers
 _reactor_counts = _defaultdict(int)
 
 # Need a pure-python class to store weakrefs to
-class _WeakrefProxy(object):
+class _WeakrefProxy:
     pass
 
 cdef class ReactorBase:
@@ -431,7 +431,7 @@ cdef class ReactorSurface:
             if self._kinetics is None:
                 raise CanteraError("Can't set coverages before assigning kinetics manager.")
 
-            if isinstance(coverages, (dict, str, unicode, bytes)):
+            if isinstance(coverages, (dict, str, bytes)):
                 self.surface.setCoverages(comp_map(coverages))
                 return
 
@@ -964,7 +964,7 @@ cdef class ReactorNet:
         """
         if isinstance(component, int):
             return self.net.sensitivity(component, p)
-        elif isinstance(component, (str, unicode, bytes)):
+        elif isinstance(component, (str, bytes)):
             return self.net.sensitivity(stringify(component), p, r)
 
     def sensitivities(self):
diff --git a/interfaces/cython/cantera/test/test_equilibrium.py b/interfaces/cython/cantera/test/test_equilibrium.py
index fa5555bbb..0ff4a0dc3 100644
--- a/interfaces/cython/cantera/test/test_equilibrium.py
+++ b/interfaces/cython/cantera/test/test_equilibrium.py
@@ -1,5 +1,3 @@
-from __future__ import division
-
 import unittest
 from os.path import join as pjoin
 
@@ -9,7 +7,7 @@ import cantera as ct
 from . import utilities
 
 
-class EquilTestCases(object):
+class EquilTestCases:
     def __init__(self, solver):
         self.solver = solver
 
diff --git a/interfaces/cython/cantera/test/test_func1.py b/interfaces/cython/cantera/test/test_func1.py
index d1d53d762..b8672d45d 100644
--- a/interfaces/cython/cantera/test/test_func1.py
+++ b/interfaces/cython/cantera/test/test_func1.py
@@ -17,7 +17,7 @@ class TestFunc1(utilities.CanteraTest):
             self.assertNear(f(t), np.sin(t)*np.sqrt(t))
 
     def test_callable(self):
-        class Multiplier(object):
+        class Multiplier:
             def __init__(self, factor):
                 self.factor = factor
             def __call__(self, t):
diff --git a/interfaces/cython/cantera/test/test_purefluid.py b/interfaces/cython/cantera/test/test_purefluid.py
index 919fe7276..a5c8502b6 100644
--- a/interfaces/cython/cantera/test/test_purefluid.py
+++ b/interfaces/cython/cantera/test/test_purefluid.py
@@ -1,5 +1,3 @@
-from __future__ import division
-
 import itertools
 import numpy as np
 
@@ -161,7 +159,7 @@ class TestPureFluid(utilities.CanteraTest):
 # To minimize errors when transcribing tabulated data, the input units here are:
 # T: K, P: MPa, rho: kg/m3, v: m3/kg, (u,h): kJ/kg, s: kJ/kg-K
 # Which are then converted to SI
-class StateData(object):
+class StateData:
     def __init__(self, phase, T, p, rho=None, v=None, u=None, h=None, s=None, relax=False):
         self.phase = phase
         self.T = T
@@ -172,7 +170,7 @@ class StateData(object):
         self.tolMod = 10.0 if relax else 1.0
 
 
-class Tolerances(object):
+class Tolerances:
     def __init__(self, p=None, u=None, s=None,
                  dUdS=None, dAdV=None, dPdT=None, hTs=None):
         self.p = p or 2e-5
@@ -184,7 +182,7 @@ class Tolerances(object):
         self.hTs = hTs or 2e-4
 
 
-class PureFluidTestCases(object):
+class PureFluidTestCases:
     """
     Test the results of pure fluid phase calculations against tabulated
     references and for consistency with basic thermodynamic relations.
diff --git a/interfaces/cython/cantera/test/test_reactor.py b/interfaces/cython/cantera/test/test_reactor.py
index e3681275a..7a8a5068a 100644
--- a/interfaces/cython/cantera/test/test_reactor.py
+++ b/interfaces/cython/cantera/test/test_reactor.py
@@ -370,8 +370,7 @@ class TestReactor(utilities.CanteraTest):
         mfc = ct.MassFlowController(self.r1, self.r2)
         mfc.set_mass_flow_rate(lambda t: eggs)
 
-        # TODO: replace with 'assertRasesRegex' after dropping Python 2.7 support
-        with self.assertRaisesRegexp(Exception, 'eggs'):
+        with self.assertRaisesRegex(Exception, 'eggs'):
             self.net.step()
 
     def test_valve1(self):
@@ -1283,7 +1282,7 @@ class TestReactorSensitivities(utilities.CanteraTest):
             self.assertNear(dtigdh_cvodes[i], dtigdh, atol=1e-14, rtol=5e-2)
 
 
-class CombustorTestImplementation(object):
+class CombustorTestImplementation:
     """
     These tests are based on the sample:
 
@@ -1371,7 +1370,7 @@ class CombustorTestImplementation(object):
             self.assertFalse(bad, bad)
 
 
-class WallTestImplementation(object):
+class WallTestImplementation:
     """
     These tests are based on the sample:
 
diff --git a/interfaces/cython/cantera/test/utilities.py b/interfaces/cython/cantera/test/utilities.py
index 5596ff2aa..a1afdd5ac 100644
--- a/interfaces/cython/cantera/test/utilities.py
+++ b/interfaces/cython/cantera/test/utilities.py
@@ -4,21 +4,10 @@ import os
 import warnings
 import shutil
 import tempfile
+import unittest
 import errno
 import cantera
 
-_ver = sys.version_info[:2]
-python_version = str(_ver[0])
-
-if _ver < (2,7) or (3,0) <= _ver < (3,2):
-    # unittest2 is a backport of the new features added to the unittest
-    # testing framework in Python 2.7 and Python 3.2. See
-    # https://pypi.python.org/pypi/unittest2 (for Python 2.x)
-    # https://pypi.python.org/pypi/unittest2py3k (for Python 3.x)
-    import unittest2 as unittest
-else:
-    import unittest
-
 
 class CanteraTest(unittest.TestCase):
     @classmethod
diff --git a/interfaces/cython/cantera/thermo.pyx b/interfaces/cython/cantera/thermo.pyx
index 57ec5e811..33c2c2db9 100644
--- a/interfaces/cython/cantera/thermo.pyx
+++ b/interfaces/cython/cantera/thermo.pyx
@@ -356,7 +356,7 @@ cdef class ThermoPhase(_SolutionBase):
         an integer. In the latter case, the index is checked for validity and
         returned. If no such element is present, an exception is thrown.
         """
-        if isinstance(element, (str, unicode, bytes)):
+        if isinstance(element, (str, bytes)):
             index = self.thermo.elementIndex(stringify(element))
         elif isinstance(element, (int, float)):
             index = <int>element
@@ -418,7 +418,7 @@ cdef class ThermoPhase(_SolutionBase):
         an integer. In the latter case, the index is checked for validity and
         returned. If no such species is present, an exception is thrown.
         """
-        if isinstance(species, (str, unicode, bytes)):
+        if isinstance(species, (str, bytes)):
             index = self.thermo.speciesIndex(stringify(species))
         elif isinstance(species, (int, float)):
             index = <int>species
@@ -443,7 +443,7 @@ cdef class ThermoPhase(_SolutionBase):
             return [self.species(i) for i in range(self.n_species)]
 
         s = Species(init=False)
-        if isinstance(k, (str, unicode, bytes)):
+        if isinstance(k, (str, bytes)):
             s._assign(self.thermo.species(stringify(k)))
         elif isinstance(k, (int, float)):
             s._assign(self.thermo.species(<int>k))
@@ -534,7 +534,7 @@ cdef class ThermoPhase(_SolutionBase):
         def __get__(self):
             return self._getArray1(thermo_getMassFractions)
         def __set__(self, Y):
-            if isinstance(Y, (str, unicode, bytes)):
+            if isinstance(Y, (str, bytes)):
                 self.thermo.setMassFractionsByName(stringify(Y))
             elif isinstance(Y, dict):
                 self.thermo.setMassFractionsByName(comp_map(Y))
@@ -555,7 +555,7 @@ cdef class ThermoPhase(_SolutionBase):
         def __get__(self):
             return self._getArray1(thermo_getMoleFractions)
         def __set__(self, X):
-            if isinstance(X, (str, unicode, bytes)):
+            if isinstance(X, (str, bytes)):
                 self.thermo.setMoleFractionsByName(stringify(X))
             elif isinstance(X, dict):
                 self.thermo.setMoleFractionsByName(comp_map(X))
@@ -1371,7 +1371,7 @@ cdef class InterfacePhase(ThermoPhase):
                 return data
 
         def __set__(self, theta):
-            if isinstance(theta, (dict, str, unicode, bytes)):
+            if isinstance(theta, (dict, str, bytes)):
                 self.surf.setCoveragesByName(comp_map(theta))
                 return
 
@@ -1569,7 +1569,7 @@ cdef class PureFluid(ThermoPhase):
             return self.s, self.v, self.X
 
 
-class Element(object):
+class Element:
     """
     An element or a named isotope defined in Cantera.
 
@@ -1624,7 +1624,7 @@ class Element(object):
                      for m in range(num_elements_defined)]
 
     def __init__(self, arg):
-        if isinstance(arg, (str, unicode, bytes)):
+        if isinstance(arg, (str, bytes)):
             try:
                 # Assume the argument is the element symbol and try to get the name
                 self._name = pystr(getElementName(stringify(arg)))
diff --git a/interfaces/cython/cantera/utils.pyx b/interfaces/cython/cantera/utils.pyx
index cb3a6adfb..4d927615d 100644
--- a/interfaces/cython/cantera/utils.pyx
+++ b/interfaces/cython/cantera/utils.pyx
@@ -5,14 +5,11 @@ import sys
 import os
 from cpython.ref cimport PyObject
 
-cdef int _pythonMajorVersion = sys.version_info[0]
-
 cdef CxxPythonLogger* _logger = new CxxPythonLogger()
 CxxSetLogger(_logger)
 
 cdef string stringify(x) except *:
     """ Converts Python strings to std::string. """
-    # This method works with both Python 2.x and 3.x.
     if isinstance(x, bytes):
         return string(<bytes>x)
     else:
@@ -20,13 +17,7 @@ cdef string stringify(x) except *:
         return string(tmp)
 
 cdef pystr(string x):
-    cdef bytes s = x.c_str()
-    if _pythonMajorVersion == 2:
-        # Python 2.x
-        return s
-    else:
-        # Python 3.x
-        return s.decode()
+    return x.decode()
 
 def add_directory(directory):
     """ Add a directory to search for Cantera data files. """
@@ -53,7 +44,7 @@ def suppress_thermo_warnings(pybool suppress=True):
     Cxx_suppress_thermo_warnings(suppress)
 
 cdef Composition comp_map(X) except *:
-    if isinstance(X, (str, unicode, bytes)):
+    if isinstance(X, (str, bytes)):
         return parseCompString(stringify(X))
 
     # assume X is dict-like