diff --git a/tests/__init__.py b/tests/__init__.py
index f23a008ce0..9ca19981bb 100644
--- a/tests/__init__.py
+++ b/tests/__init__.py
@@ -44,5 +44,5 @@
     get_optional_distribution_deps,
     get_present_optional_import_deps,
     no_internet_connection,
+    assert_domain_equal,
 )
-from .testcase import TestCase
diff --git a/tests/integration/test_experiments.py b/tests/integration/test_experiments.py
index 30db3e446c..f43c7293d2 100644
--- a/tests/integration/test_experiments.py
+++ b/tests/integration/test_experiments.py
@@ -1,13 +1,13 @@
 #
 # Test some experiments
 #
-from tests import TestCase
+
 import pybamm
 import numpy as np
 import unittest
 
 
-class TestExperiments(TestCase):
+class TestExperiments(unittest.TestCase):
     def test_discharge_rest_charge(self):
         experiment = pybamm.Experiment(
             [
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py
index 136350fe65..57067d6e3b 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py
@@ -1,14 +1,14 @@
 #
 # Tests for the surface formulation
 #
-from tests import TestCase
+
 import pybamm
 import numpy as np
 import unittest
 from tests import StandardOutputComparison
 
 
-class TestCompareOutputs(TestCase):
+class TestCompareOutputs(unittest.TestCase):
     def test_compare_outputs_surface_form(self):
         # load models
         options = [
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs_two_phase.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs_two_phase.py
index b0c0fe5898..ed6d707d77 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs_two_phase.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs_two_phase.py
@@ -4,10 +4,9 @@
 import pybamm
 import numpy as np
 import unittest
-from tests import TestCase
 
 
-class TestCompareOutputsTwoPhase(TestCase):
+class TestCompareOutputsTwoPhase(unittest.TestCase):
     def compare_outputs_two_phase_graphite_graphite(self, model_class):
         """
         Check that a two-phase graphite-graphite model gives the same results as a
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py
index 083e3b648b..1b6b29ea10 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py
@@ -1,7 +1,7 @@
 #
 # Tests for the lithium-ion DFN model
 #
-from tests import TestCase
+
 import pybamm
 import tests
 import numpy as np
@@ -9,7 +9,7 @@
 from tests import BaseIntegrationTestLithiumIon
 
 
-class TestDFN(BaseIntegrationTestLithiumIon, TestCase):
+class TestDFN(BaseIntegrationTestLithiumIon, unittest.TestCase):
     def setUp(self):
         self.model = pybamm.lithium_ion.DFN
 
@@ -35,7 +35,7 @@ def positive_radius(x):
         self.run_basic_processing_test({}, parameter_values=param)
 
 
-class TestDFNWithSizeDistribution(TestCase):
+class TestDFNWithSizeDistribution(unittest.TestCase):
     def setUp(self):
         params = pybamm.ParameterValues("Marquis2019")
         self.params = pybamm.get_size_distribution_parameters(params)
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py
index 82d228badb..7f5a5941c8 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py
@@ -1,14 +1,14 @@
 #
 # Tests for the lithium-ion MPM model
 #
-from tests import TestCase
+
 import pybamm
 import tests
 import numpy as np
 import unittest
 
 
-class TestMPM(TestCase):
+class TestMPM(unittest.TestCase):
     def test_basic_processing(self):
         options = {"thermal": "isothermal"}
         model = pybamm.lithium_ion.MPM(options)
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_thermal_models.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_thermal_models.py
index 7de36da042..2b74366652 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_thermal_models.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_thermal_models.py
@@ -5,10 +5,9 @@
 import pybamm
 import numpy as np
 import unittest
-from tests import TestCase
 
 
-class TestThermal(TestCase):
+class TestThermal(unittest.TestCase):
     def test_consistent_cooling(self):
         "Test the cooling is consistent between the 1D, 1+1D and 2+1D SPMe models"
 
diff --git a/tests/integration/test_models/test_submodels/test_external_circuit/test_function_control.py b/tests/integration/test_models/test_submodels/test_external_circuit/test_function_control.py
index 8909367892..c67604412c 100644
--- a/tests/integration/test_models/test_submodels/test_external_circuit/test_function_control.py
+++ b/tests/integration/test_models/test_submodels/test_external_circuit/test_function_control.py
@@ -1,13 +1,13 @@
 #
 # Test function control submodel
 #
-from tests import TestCase
+
 import numpy as np
 import pybamm
 import unittest
 
 
-class TestFunctionControl(TestCase):
+class TestFunctionControl(unittest.TestCase):
     def test_constant_current(self):
         def constant_current(variables):
             I = variables["Current [A]"]
diff --git a/tests/integration/test_models/test_submodels/test_interface/test_butler_volmer.py b/tests/integration/test_models/test_submodels/test_interface/test_butler_volmer.py
index 536573bad1..f1f02350cd 100644
--- a/tests/integration/test_models/test_submodels/test_interface/test_butler_volmer.py
+++ b/tests/integration/test_models/test_submodels/test_interface/test_butler_volmer.py
@@ -1,7 +1,7 @@
 #
 # Tests for the electrode-electrolyte interface equations
 #
-from tests import TestCase
+
 import pybamm
 from tests import get_discretisation_for_testing
 
@@ -9,7 +9,7 @@
 import numpy as np
 
 
-class TestButlerVolmer(TestCase):
+class TestButlerVolmer(unittest.TestCase):
     def setUp(self):
         self.delta_phi_s_n = pybamm.Variable(
             "surface potential difference [V]",
diff --git a/tests/integration/test_models/test_submodels/test_interface/test_lead_acid.py b/tests/integration/test_models/test_submodels/test_interface/test_lead_acid.py
index daf05ff6ca..ca41414c70 100644
--- a/tests/integration/test_models/test_submodels/test_interface/test_lead_acid.py
+++ b/tests/integration/test_models/test_submodels/test_interface/test_lead_acid.py
@@ -1,13 +1,13 @@
 #
 # Tests for the electrode-electrolyte interface equations for lead-acid models
 #
-from tests import TestCase
+
 import pybamm
 from tests import get_discretisation_for_testing
 import unittest
 
 
-class TestMainReaction(TestCase):
+class TestMainReaction(unittest.TestCase):
     def setUp(self):
         c_e_n = pybamm.Variable(
             "Negative electrolyte concentration [mol.m-3]",
diff --git a/tests/integration/test_models/test_submodels/test_interface/test_lithium_ion.py b/tests/integration/test_models/test_submodels/test_interface/test_lithium_ion.py
index 90d27e2de7..135e29ad8b 100644
--- a/tests/integration/test_models/test_submodels/test_interface/test_lithium_ion.py
+++ b/tests/integration/test_models/test_submodels/test_interface/test_lithium_ion.py
@@ -1,7 +1,7 @@
 #
 # Tests for the electrode-electrolyte interface equations for lithium-ion models
 #
-from tests import TestCase
+
 import pybamm
 from tests import get_discretisation_for_testing
 
@@ -9,7 +9,7 @@
 import numpy as np
 
 
-class TestExchangeCurrentDensity(TestCase):
+class TestExchangeCurrentDensity(unittest.TestCase):
     def setUp(self):
         c_e_n = pybamm.Variable("concentration", domain=["negative electrode"])
         c_e_s = pybamm.Variable("concentration", domain=["separator"])
diff --git a/tests/integration/test_solvers/test_idaklu.py b/tests/integration/test_solvers/test_idaklu.py
index eecc5dfe2b..693ec849a9 100644
--- a/tests/integration/test_solvers/test_idaklu.py
+++ b/tests/integration/test_solvers/test_idaklu.py
@@ -1,12 +1,12 @@
 import pybamm
 import numpy as np
 import sys
-from tests import TestCase
+
 import unittest
 
 
 @unittest.skipIf(not pybamm.have_idaklu(), "idaklu solver is not installed")
-class TestIDAKLUSolver(TestCase):
+class TestIDAKLUSolver(unittest.TestCase):
     def test_on_spme(self):
         model = pybamm.lithium_ion.SPMe()
         geometry = model.default_geometry
diff --git a/tests/integration/test_solvers/test_solution.py b/tests/integration/test_solvers/test_solution.py
index d45df9fd4c..bb276a62a2 100644
--- a/tests/integration/test_solvers/test_solution.py
+++ b/tests/integration/test_solvers/test_solution.py
@@ -1,13 +1,13 @@
 #
 # Tests for the Solution class
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 import numpy as np
 
 
-class TestSolution(TestCase):
+class TestSolution(unittest.TestCase):
     def test_append(self):
         model = pybamm.lithium_ion.SPMe()
         # create geometry
diff --git a/tests/integration/test_spatial_methods/test_finite_volume.py b/tests/integration/test_spatial_methods/test_finite_volume.py
index ce2d07c2de..d331554853 100644
--- a/tests/integration/test_spatial_methods/test_finite_volume.py
+++ b/tests/integration/test_spatial_methods/test_finite_volume.py
@@ -1,7 +1,7 @@
 #
 # Test for the operator class
 #
-from tests import TestCase
+
 import pybamm
 from tests import (
     get_mesh_for_testing,
@@ -13,7 +13,7 @@
 import unittest
 
 
-class TestFiniteVolumeConvergence(TestCase):
+class TestFiniteVolumeConvergence(unittest.TestCase):
     def test_grad_div_broadcast(self):
         # create mesh and discretisation
         spatial_methods = {"macroscale": pybamm.FiniteVolume()}
@@ -319,7 +319,7 @@ def solve_laplace_equation(coord_sys="cartesian"):
     return solver.solve(model)
 
 
-class TestFiniteVolumeLaplacian(TestCase):
+class TestFiniteVolumeLaplacian(unittest.TestCase):
     def test_laplacian_cartesian(self):
         solution = solve_laplace_equation(coord_sys="cartesian")
         np.testing.assert_array_almost_equal(
@@ -374,7 +374,7 @@ def solve_advection_equation(direction="upwind", source=1, bc=0):
     return solver.solve(model, [0, 1])
 
 
-class TestUpwindDownwind(TestCase):
+class TestUpwindDownwind(unittest.TestCase):
     def test_upwind(self):
         solution = solve_advection_equation("upwind")
         np.testing.assert_array_almost_equal(
diff --git a/tests/integration/test_spatial_methods/test_spectral_volume.py b/tests/integration/test_spatial_methods/test_spectral_volume.py
index 4b998cb0ed..deba95ebac 100644
--- a/tests/integration/test_spatial_methods/test_spectral_volume.py
+++ b/tests/integration/test_spatial_methods/test_spectral_volume.py
@@ -1,7 +1,7 @@
 #
 # Test for the operator class
 #
-from tests import TestCase
+
 import pybamm
 
 import numpy as np
@@ -76,7 +76,7 @@ def get_p2d_mesh_for_testing(xpts=None, rpts=10):
     return get_mesh_for_testing(xpts=xpts, rpts=rpts, geometry=geometry)
 
 
-class TestSpectralVolumeConvergence(TestCase):
+class TestSpectralVolumeConvergence(unittest.TestCase):
     def test_grad_div_broadcast(self):
         # create mesh and discretisation
         spatial_methods = {"macroscale": pybamm.SpectralVolume()}
diff --git a/tests/shared.py b/tests/shared.py
index 2fa9c24960..48e54e19d8 100644
--- a/tests/shared.py
+++ b/tests/shared.py
@@ -338,3 +338,10 @@ def no_internet_connection():
         return False
     except socket.gaierror:
         return True
+
+
+def assert_domain_equal(a, b):
+    "Check that two domains are equal, ignoring empty domains"
+    a_dict = {k: v for k, v in a.items() if v != []}
+    b_dict = {k: v for k, v in b.items() if v != []}
+    assert a_dict == b_dict
diff --git a/tests/unit/test_batch_study.py b/tests/unit/test_batch_study.py
index 2713a3f35b..9d13d71d9d 100644
--- a/tests/unit/test_batch_study.py
+++ b/tests/unit/test_batch_study.py
@@ -2,14 +2,13 @@
 Tests for the batch_study.py
 """
 
-from tests import TestCase
 import os
 import pybamm
 import unittest
 from tempfile import TemporaryDirectory
 
 
-class TestBatchStudy(TestCase):
+class TestBatchStudy(unittest.TestCase):
     def test_solve(self):
         spm = pybamm.lithium_ion.SPM()
         spm_uniform = pybamm.lithium_ion.SPM({"particle": "uniform profile"})
diff --git a/tests/unit/test_callbacks.py b/tests/unit/test_callbacks.py
index 649c7d9ec8..bc523d404f 100644
--- a/tests/unit/test_callbacks.py
+++ b/tests/unit/test_callbacks.py
@@ -1,7 +1,7 @@
 #
 # Tests the citations class.
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 import os
@@ -18,7 +18,7 @@ def on_experiment_end(self, logs):
             print(self.name, file=f)
 
 
-class TestCallbacks(TestCase):
+class TestCallbacks(unittest.TestCase):
     def tearDown(self):
         # Remove any test log files that were created, even if the test fails
         for logfile in ["test_callback.log", "test_callback_2.log"]:
diff --git a/tests/unit/test_discretisations/test_discretisation.py b/tests/unit/test_discretisations/test_discretisation.py
index 61e3c2bc6e..bf698dd39c 100644
--- a/tests/unit/test_discretisations/test_discretisation.py
+++ b/tests/unit/test_discretisations/test_discretisation.py
@@ -1,7 +1,7 @@
 #
 # Tests for the base model class
 #
-from tests import TestCase
+
 import pybamm
 
 import numpy as np
@@ -18,7 +18,7 @@
 from scipy.sparse.linalg import inv
 
 
-class TestDiscretise(TestCase):
+class TestDiscretise(unittest.TestCase):
     def test_concatenate_in_order(self):
         a = pybamm.Variable("a")
         b = pybamm.Variable("b")
diff --git a/tests/unit/test_experiments/test_experiment.py b/tests/unit/test_experiments/test_experiment.py
index 6c342bd269..87672ee0aa 100644
--- a/tests/unit/test_experiments/test_experiment.py
+++ b/tests/unit/test_experiments/test_experiment.py
@@ -1,13 +1,13 @@
 #
 # Test the base experiment class
 #
-from tests import TestCase
+
 from datetime import datetime
 import pybamm
 import unittest
 
 
-class TestExperiment(TestCase):
+class TestExperiment(unittest.TestCase):
     def test_cycle_unpacking(self):
         experiment = pybamm.Experiment(
             [
diff --git a/tests/unit/test_experiments/test_simulation_with_experiment.py b/tests/unit/test_experiments/test_simulation_with_experiment.py
index 7784726fc2..5dfa6b9aaf 100644
--- a/tests/unit/test_experiments/test_simulation_with_experiment.py
+++ b/tests/unit/test_experiments/test_simulation_with_experiment.py
@@ -1,7 +1,6 @@
 #
 # Test setting up a simulation with an experiment
 #
-from tests import TestCase
 import casadi
 import pybamm
 import numpy as np
@@ -16,7 +15,7 @@ def default_duration(self, value):
         return 1
 
 
-class TestSimulationExperiment(TestCase):
+class TestSimulationExperiment(unittest.TestCase):
     def test_set_up(self):
         experiment = pybamm.Experiment(
             [
diff --git a/tests/unit/test_expression_tree/test_averages.py b/tests/unit/test_expression_tree/test_averages.py
index 1351ee3c05..1f1385db65 100644
--- a/tests/unit/test_expression_tree/test_averages.py
+++ b/tests/unit/test_expression_tree/test_averages.py
@@ -2,12 +2,12 @@
 # Tests for the Unary Operator classes
 #
 import unittest
-from tests import TestCase
 import numpy as np
 import pybamm
+from tests import assert_domain_equal
 
 
-class TestUnaryOperators(TestCase):
+class TestUnaryOperators(unittest.TestCase):
     def test_x_average(self):
         a = pybamm.Scalar(4)
         average_a = pybamm.x_average(a)
@@ -112,7 +112,7 @@ def test_x_average(self):
         )
         average_conc_broad = pybamm.x_average(conc_broad)
         self.assertIsInstance(average_conc_broad, pybamm.FullBroadcast)
-        self.assertDomainEqual(
+        assert_domain_equal(
             average_conc_broad.domains,
             {"primary": ["current collector"], "secondary": ["test"]},
         )
diff --git a/tests/unit/test_expression_tree/test_binary_operators.py b/tests/unit/test_expression_tree/test_binary_operators.py
index 0d843c6db2..de58cdab28 100644
--- a/tests/unit/test_expression_tree/test_binary_operators.py
+++ b/tests/unit/test_expression_tree/test_binary_operators.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Binary Operator classes
 #
-from tests import TestCase
+
 import unittest
 import unittest.mock as mock
 
@@ -19,7 +19,7 @@
 }
 
 
-class TestBinaryOperators(TestCase):
+class TestBinaryOperators(unittest.TestCase):
     def test_binary_operator(self):
         a = pybamm.Symbol("a")
         b = pybamm.Symbol("b")
diff --git a/tests/unit/test_expression_tree/test_broadcasts.py b/tests/unit/test_expression_tree/test_broadcasts.py
index be8fe1a677..e5516395e9 100644
--- a/tests/unit/test_expression_tree/test_broadcasts.py
+++ b/tests/unit/test_expression_tree/test_broadcasts.py
@@ -2,13 +2,13 @@
 # Tests for the Broadcast class
 #
 import unittest
-from tests import TestCase
+from tests import assert_domain_equal
 import numpy as np
 
 import pybamm
 
 
-class TestBroadcasts(TestCase):
+class TestBroadcasts(unittest.TestCase):
     def test_primary_broadcast(self):
         a = pybamm.Symbol("a")
         broad_a = pybamm.PrimaryBroadcast(a, ["negative electrode"])
@@ -24,7 +24,7 @@ def test_primary_broadcast(self):
             auxiliary_domains={"secondary": "current collector"},
         )
         broad_a = pybamm.PrimaryBroadcast(a, ["negative particle"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             broad_a.domains,
             {
                 "primary": ["negative particle"],
@@ -41,7 +41,7 @@ def test_primary_broadcast(self):
             },
         )
         broad_a = pybamm.PrimaryBroadcast(a, ["negative particle"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             broad_a.domains,
             {
                 "primary": ["negative particle"],
@@ -83,7 +83,7 @@ def test_secondary_broadcast(self):
             auxiliary_domains={"secondary": "current collector"},
         )
         broad_a = pybamm.SecondaryBroadcast(a, ["negative electrode"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             broad_a.domains,
             {
                 "primary": ["negative particle"],
@@ -93,7 +93,7 @@ def test_secondary_broadcast(self):
         )
         self.assertTrue(broad_a.broadcasts_to_nodes)
         broadbroad_a = pybamm.SecondaryBroadcast(broad_a, ["negative particle size"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             broadbroad_a.domains,
             {
                 "primary": ["negative particle"],
@@ -140,7 +140,7 @@ def test_tertiary_broadcast(self):
             },
         )
         broad_a = pybamm.TertiaryBroadcast(a, "negative electrode")
-        self.assertDomainEqual(
+        assert_domain_equal(
             broad_a.domains,
             {
                 "primary": ["negative particle"],
@@ -269,7 +269,7 @@ def test_broadcast_to_edges(self):
             auxiliary_domains={"secondary": "current collector"},
         )
         broad_a = pybamm.SecondaryBroadcastToEdges(a, ["negative electrode"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             broad_a.domains,
             {
                 "primary": ["negative particle"],
@@ -290,7 +290,7 @@ def test_broadcast_to_edges(self):
             },
         )
         broad_a = pybamm.TertiaryBroadcastToEdges(a, ["negative electrode"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             broad_a.domains,
             {
                 "primary": ["negative particle"],
diff --git a/tests/unit/test_expression_tree/test_concatenations.py b/tests/unit/test_expression_tree/test_concatenations.py
index 9c0ee88964..b8536017d9 100644
--- a/tests/unit/test_expression_tree/test_concatenations.py
+++ b/tests/unit/test_expression_tree/test_concatenations.py
@@ -3,7 +3,8 @@
 #
 import unittest
 import unittest.mock as mock
-from tests import TestCase
+from tests import assert_domain_equal
+
 
 import numpy as np
 
@@ -12,7 +13,7 @@
 from tests import get_discretisation_for_testing, get_mesh_for_testing
 
 
-class TestConcatenations(TestCase):
+class TestConcatenations(unittest.TestCase):
     def test_base_concatenation(self):
         a = pybamm.Symbol("a", domain="test a")
         b = pybamm.Symbol("b", domain="test b")
@@ -81,7 +82,7 @@ def test_concatenation_auxiliary_domains(self):
             auxiliary_domains={"secondary": "current collector"},
         )
         conc = pybamm.concatenation(a, b)
-        self.assertDomainEqual(
+        assert_domain_equal(
             conc.domains,
             {
                 "primary": ["negative electrode", "separator", "positive electrode"],
@@ -158,7 +159,7 @@ def test_concatenation_simplify(self):
         concat = pybamm.concatenation(a, b, c)
         self.assertIsInstance(concat, pybamm.FullBroadcast)
         self.assertEqual(concat.orphans[0], pybamm.Scalar(0))
-        self.assertDomainEqual(
+        assert_domain_equal(
             concat.domains,
             {
                 "primary": ["negative electrode", "separator", "positive electrode"],
diff --git a/tests/unit/test_expression_tree/test_functions.py b/tests/unit/test_expression_tree/test_functions.py
index b619331577..003a3eae9a 100644
--- a/tests/unit/test_expression_tree/test_functions.py
+++ b/tests/unit/test_expression_tree/test_functions.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Function classes
 #
-from tests import TestCase
+
 import unittest
 import unittest.mock as mock
 
@@ -16,7 +16,7 @@
 )
 
 
-class TestFunction(TestCase):
+class TestFunction(unittest.TestCase):
     def test_number_input(self):
         # with numbers
         log = pybamm.Function(np.log, 10)
@@ -107,7 +107,7 @@ def test_to_from_json_error(self):
             pybamm.Function._from_json({})
 
 
-class TestSpecificFunctions(TestCase):
+class TestSpecificFunctions(unittest.TestCase):
     def test_to_json(self):
         a = pybamm.InputParameter("a")
         fun = pybamm.cos(a)
@@ -453,7 +453,7 @@ def test_erfc(self):
         )
 
 
-class TestNonObjectFunctions(TestCase):
+class TestNonObjectFunctions(unittest.TestCase):
     def test_normal_pdf(self):
         x = pybamm.InputParameter("x")
         mu = pybamm.InputParameter("mu")
diff --git a/tests/unit/test_expression_tree/test_interpolant.py b/tests/unit/test_expression_tree/test_interpolant.py
index f5ded9cf8e..e40b9daab0 100644
--- a/tests/unit/test_expression_tree/test_interpolant.py
+++ b/tests/unit/test_expression_tree/test_interpolant.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Function classes
 #
-from tests import TestCase
+
 import pybamm
 
 import unittest
@@ -9,7 +9,7 @@
 import numpy as np
 
 
-class TestInterpolant(TestCase):
+class TestInterpolant(unittest.TestCase):
     def test_errors(self):
         with self.assertRaisesRegex(ValueError, "x1"):
             pybamm.Interpolant(np.ones(10), np.ones(11), pybamm.Symbol("a"))
diff --git a/tests/unit/test_expression_tree/test_operations/test_convert_to_casadi.py b/tests/unit/test_expression_tree/test_operations/test_convert_to_casadi.py
index 20c2e40db0..e3301fdcc3 100644
--- a/tests/unit/test_expression_tree/test_operations/test_convert_to_casadi.py
+++ b/tests/unit/test_expression_tree/test_operations/test_convert_to_casadi.py
@@ -1,7 +1,7 @@
 #
 # Test for the Simplify class
 #
-from tests import TestCase
+
 import casadi
 import numpy as np
 import pybamm
@@ -10,7 +10,7 @@
 from scipy import special
 
 
-class TestCasadiConverter(TestCase):
+class TestCasadiConverter(unittest.TestCase):
     def assert_casadi_equal(self, a, b, evalf=False):
         if evalf is True:
             self.assertTrue((casadi.evalf(a) - casadi.evalf(b)).is_zero())
diff --git a/tests/unit/test_expression_tree/test_operations/test_copy.py b/tests/unit/test_expression_tree/test_operations/test_copy.py
index ca6ac6c448..5c1e3b29ca 100644
--- a/tests/unit/test_expression_tree/test_operations/test_copy.py
+++ b/tests/unit/test_expression_tree/test_operations/test_copy.py
@@ -1,14 +1,14 @@
 #
 # Test for making copies
 #
-from tests import TestCase
+
 import numpy as np
 import pybamm
 import unittest
 from tests import get_mesh_for_testing
 
 
-class TestCopy(TestCase):
+class TestCopy(unittest.TestCase):
     def test_symbol_create_copy(self):
         a = pybamm.Parameter("a")
         b = pybamm.Parameter("b")
diff --git a/tests/unit/test_expression_tree/test_operations/test_evaluate_python.py b/tests/unit/test_expression_tree/test_operations/test_evaluate_python.py
index 6e1b155eca..667522c286 100644
--- a/tests/unit/test_expression_tree/test_operations/test_evaluate_python.py
+++ b/tests/unit/test_expression_tree/test_operations/test_evaluate_python.py
@@ -1,7 +1,7 @@
 #
 # Test for the evaluate-to-python functions
 #
-from tests import TestCase
+
 import pybamm
 
 from tests import get_discretisation_for_testing, get_1p1d_discretisation_for_testing
@@ -18,7 +18,7 @@
 )
 
 
-class TestEvaluate(TestCase):
+class TestEvaluate(unittest.TestCase):
     def test_find_symbols(self):
         a = pybamm.StateVector(slice(0, 1))
         b = pybamm.StateVector(slice(1, 2))
diff --git a/tests/unit/test_expression_tree/test_operations/test_jac.py b/tests/unit/test_expression_tree/test_operations/test_jac.py
index bec38d7243..34a65f007b 100644
--- a/tests/unit/test_expression_tree/test_operations/test_jac.py
+++ b/tests/unit/test_expression_tree/test_operations/test_jac.py
@@ -1,7 +1,7 @@
 #
 # Tests for the jacobian methods
 #
-from tests import TestCase
+
 import pybamm
 
 import numpy as np
@@ -10,7 +10,7 @@
 from tests import get_mesh_for_testing
 
 
-class TestJacobian(TestCase):
+class TestJacobian(unittest.TestCase):
     def test_variable_is_statevector(self):
         a = pybamm.Symbol("a")
         with self.assertRaisesRegex(
diff --git a/tests/unit/test_expression_tree/test_operations/test_jac_2D.py b/tests/unit/test_expression_tree/test_operations/test_jac_2D.py
index 2be71d85e4..e7001b184b 100644
--- a/tests/unit/test_expression_tree/test_operations/test_jac_2D.py
+++ b/tests/unit/test_expression_tree/test_operations/test_jac_2D.py
@@ -1,7 +1,7 @@
 #
 # Tests for the jacobian methods for two-dimensional objects
 #
-from tests import TestCase
+
 import pybamm
 
 import numpy as np
@@ -12,7 +12,7 @@
 )
 
 
-class TestJacobian(TestCase):
+class TestJacobian(unittest.TestCase):
     def test_linear(self):
         y = pybamm.StateVector(slice(0, 8))
         u = pybamm.StateVector(slice(0, 2), slice(4, 6))
diff --git a/tests/unit/test_expression_tree/test_operations/test_latexify.py b/tests/unit/test_expression_tree/test_operations/test_latexify.py
index 0340a1d53f..746efb4269 100644
--- a/tests/unit/test_expression_tree/test_operations/test_latexify.py
+++ b/tests/unit/test_expression_tree/test_operations/test_latexify.py
@@ -2,7 +2,6 @@
 Tests for the latexify.py
 """
 
-from tests import TestCase
 import os
 import platform
 import unittest
@@ -11,7 +10,7 @@
 import pybamm
 
 
-class TestLatexify(TestCase):
+class TestLatexify(unittest.TestCase):
     def test_latexify(self):
         model_dfn = pybamm.lithium_ion.DFN()
         func_dfn = str(model_dfn.latexify())
diff --git a/tests/unit/test_expression_tree/test_parameter.py b/tests/unit/test_expression_tree/test_parameter.py
index ba85203bbd..15ced916c3 100644
--- a/tests/unit/test_expression_tree/test_parameter.py
+++ b/tests/unit/test_expression_tree/test_parameter.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Parameter class
 #
-from tests import TestCase
+
 import numbers
 import unittest
 
@@ -9,7 +9,7 @@
 import sympy
 
 
-class TestParameter(TestCase):
+class TestParameter(unittest.TestCase):
     def test_parameter_init(self):
         a = pybamm.Parameter("a")
         self.assertEqual(a.name, "a")
@@ -41,7 +41,7 @@ def test_to_json_error(self):
             pybamm.Parameter._from_json({})
 
 
-class TestFunctionParameter(TestCase):
+class TestFunctionParameter(unittest.TestCase):
     def test_function_parameter_init(self):
         var = pybamm.Variable("var")
         func = pybamm.FunctionParameter("func", {"var": var})
diff --git a/tests/unit/test_expression_tree/test_state_vector.py b/tests/unit/test_expression_tree/test_state_vector.py
index 18025c0aa3..0cbc6faa4d 100644
--- a/tests/unit/test_expression_tree/test_state_vector.py
+++ b/tests/unit/test_expression_tree/test_state_vector.py
@@ -1,7 +1,7 @@
 #
 # Tests for the State Vector class
 #
-from tests import TestCase
+
 import pybamm
 import numpy as np
 
@@ -9,7 +9,7 @@
 import unittest.mock as mock
 
 
-class TestStateVector(TestCase):
+class TestStateVector(unittest.TestCase):
     def test_evaluate(self):
         sv = pybamm.StateVector(slice(0, 10))
         y = np.linspace(0, 2, 19)
@@ -97,7 +97,7 @@ def test_to_from_json(self):
         pybamm.settings.debug_mode = original_debug_mode
 
 
-class TestStateVectorDot(TestCase):
+class TestStateVectorDot(unittest.TestCase):
     def test_evaluate(self):
         sv = pybamm.StateVectorDot(slice(0, 10))
         y_dot = np.linspace(0, 2, 19)
diff --git a/tests/unit/test_expression_tree/test_symbol.py b/tests/unit/test_expression_tree/test_symbol.py
index 6a3a0c16f8..f829b9122e 100644
--- a/tests/unit/test_expression_tree/test_symbol.py
+++ b/tests/unit/test_expression_tree/test_symbol.py
@@ -1,7 +1,7 @@
 #
 # Test for the Symbol class
 #
-from tests import TestCase
+
 import os
 import unittest
 from tempfile import TemporaryDirectory
@@ -14,7 +14,7 @@
 import sympy
 
 
-class TestSymbol(TestCase):
+class TestSymbol(unittest.TestCase):
     def test_symbol_init(self):
         sym = pybamm.Symbol("a symbol")
         with self.assertRaises(TypeError):
@@ -520,7 +520,7 @@ def test_to_from_json(self):
         self.assertEqual(pybamm.Symbol._from_json(json_dict), symp)
 
 
-class TestIsZero(TestCase):
+class TestIsZero(unittest.TestCase):
     def test_is_scalar_zero(self):
         a = pybamm.Scalar(0)
         b = pybamm.Scalar(2)
diff --git a/tests/unit/test_expression_tree/test_symbolic_diff.py b/tests/unit/test_expression_tree/test_symbolic_diff.py
index a6493b66ce..fb08740305 100644
--- a/tests/unit/test_expression_tree/test_symbolic_diff.py
+++ b/tests/unit/test_expression_tree/test_symbolic_diff.py
@@ -1,14 +1,14 @@
 #
 # Tests for the symbolic differentiation methods
 #
-from tests import TestCase
+
 import numpy as np
 import pybamm
 import unittest
 from numpy import testing
 
 
-class TestSymbolicDifferentiation(TestCase):
+class TestSymbolicDifferentiation(unittest.TestCase):
     def test_advanced(self):
         a = pybamm.StateVector(slice(0, 1))
         b = pybamm.StateVector(slice(1, 2))
diff --git a/tests/unit/test_expression_tree/test_unary_operators.py b/tests/unit/test_expression_tree/test_unary_operators.py
index 094c6c5661..216d7e736d 100644
--- a/tests/unit/test_expression_tree/test_unary_operators.py
+++ b/tests/unit/test_expression_tree/test_unary_operators.py
@@ -2,7 +2,7 @@
 # Tests for the Unary Operator classes
 #
 import unittest
-from tests import TestCase
+
 import unittest.mock as mock
 
 import numpy as np
@@ -10,12 +10,13 @@
 import sympy
 from sympy.vector.operators import Divergence as sympy_Divergence
 from sympy.vector.operators import Gradient as sympy_Gradient
+from tests import assert_domain_equal
 
 import pybamm
 from pybamm.util import have_optional_dependency
 
 
-class TestUnaryOperators(TestCase):
+class TestUnaryOperators(unittest.TestCase):
     def test_unary_operator(self):
         a = pybamm.Symbol("a", domain=["test"])
         un = pybamm.UnaryOperator("unary test", a)
@@ -291,7 +292,7 @@ def test_integral(self):
         self.assertEqual(inta.name, "integral dx ['negative electrode']")
         self.assertEqual(inta.children[0].name, a.name)
         self.assertEqual(inta.integration_variable[0], x)
-        self.assertDomainEqual(inta.domains, {})
+        assert_domain_equal(inta.domains, {})
         # space integral with secondary domain
         a_sec = pybamm.Symbol(
             "a",
@@ -300,7 +301,7 @@ def test_integral(self):
         )
         x = pybamm.SpatialVariable("x", ["negative electrode"])
         inta_sec = pybamm.Integral(a_sec, x)
-        self.assertDomainEqual(inta_sec.domains, {"primary": ["current collector"]})
+        assert_domain_equal(inta_sec.domains, {"primary": ["current collector"]})
         # space integral with tertiary domain
         a_tert = pybamm.Symbol(
             "a",
@@ -312,7 +313,7 @@ def test_integral(self):
         )
         x = pybamm.SpatialVariable("x", ["negative electrode"])
         inta_tert = pybamm.Integral(a_tert, x)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_tert.domains,
             {"primary": ["current collector"], "secondary": ["some extra domain"]},
         )
@@ -327,7 +328,7 @@ def test_integral(self):
             },
         )
         inta_quat = pybamm.Integral(a_quat, x)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_quat.domains,
             {
                 "primary": ["current collector"],
@@ -340,16 +341,16 @@ def test_integral(self):
         y = pybamm.SpatialVariable("y", ["current collector"])
         # without a tertiary domain
         inta_sec_y = pybamm.Integral(a_sec, y)
-        self.assertDomainEqual(inta_sec_y.domains, {"primary": ["negative electrode"]})
+        assert_domain_equal(inta_sec_y.domains, {"primary": ["negative electrode"]})
         # with a tertiary domain
         inta_tert_y = pybamm.Integral(a_tert, y)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_tert_y.domains,
             {"primary": ["negative electrode"], "secondary": ["some extra domain"]},
         )
         # with a quaternary domain
         inta_quat_y = pybamm.Integral(a_quat, y)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_quat_y.domains,
             {
                 "primary": ["negative electrode"],
@@ -361,13 +362,13 @@ def test_integral(self):
         # space integral *in* tertiary domain
         z = pybamm.SpatialVariable("z", ["some extra domain"])
         inta_tert_z = pybamm.Integral(a_tert, z)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_tert_z.domains,
             {"primary": ["negative electrode"], "secondary": ["current collector"]},
         )
         # with a quaternary domain
         inta_quat_z = pybamm.Integral(a_quat, z)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_quat_z.domains,
             {
                 "primary": ["negative electrode"],
@@ -379,7 +380,7 @@ def test_integral(self):
         # space integral *in* quaternary domain
         Z = pybamm.SpatialVariable("Z", ["another extra domain"])
         inta_quat_Z = pybamm.Integral(a_quat, Z)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_quat_Z.domains,
             {
                 "primary": ["negative electrode"],
@@ -406,7 +407,7 @@ def test_integral(self):
         self.assertEqual(inta.integration_variable[0], x)
         self.assertEqual(inta.domain, ["negative electrode"])
         inta_sec = pybamm.IndefiniteIntegral(a_sec, x)
-        self.assertDomainEqual(
+        assert_domain_equal(
             inta_sec.domains,
             {"primary": ["negative electrode"], "secondary": ["current collector"]},
         )
@@ -557,7 +558,7 @@ def test_delta_function(self):
         a = pybamm.Symbol("a", domain="some domain")
         delta_a = pybamm.DeltaFunction(a, "left", "another domain")
         self.assertEqual(delta_a.side, "left")
-        self.assertDomainEqual(
+        assert_domain_equal(
             delta_a.domains,
             {"primary": ["another domain"], "secondary": ["some domain"]},
         )
@@ -596,7 +597,7 @@ def test_boundary_value(self):
         boundary_a = pybamm.boundary_value(a, "right")
         self.assertIsInstance(boundary_a, pybamm.BoundaryValue)
         self.assertEqual(boundary_a.side, "right")
-        self.assertDomainEqual(boundary_a.domains, {})
+        assert_domain_equal(boundary_a.domains, {})
         # test with secondary domain
         a_sec = pybamm.Symbol(
             "a",
@@ -604,9 +605,7 @@ def test_boundary_value(self):
             auxiliary_domains={"secondary": "current collector"},
         )
         boundary_a_sec = pybamm.boundary_value(a_sec, "right")
-        self.assertDomainEqual(
-            boundary_a_sec.domains, {"primary": ["current collector"]}
-        )
+        assert_domain_equal(boundary_a_sec.domains, {"primary": ["current collector"]})
         # test with secondary domain and tertiary domain
         a_tert = pybamm.Symbol(
             "a",
@@ -614,7 +613,7 @@ def test_boundary_value(self):
             auxiliary_domains={"secondary": "current collector", "tertiary": "bla"},
         )
         boundary_a_tert = pybamm.boundary_value(a_tert, "right")
-        self.assertDomainEqual(
+        assert_domain_equal(
             boundary_a_tert.domains,
             {"primary": ["current collector"], "secondary": ["bla"]},
         )
@@ -630,7 +629,7 @@ def test_boundary_value(self):
         )
         boundary_a_quat = pybamm.boundary_value(a_quat, "right")
         self.assertEqual(boundary_a_quat.domain, ["current collector"])
-        self.assertDomainEqual(
+        assert_domain_equal(
             boundary_a_quat.domains,
             {
                 "primary": ["current collector"],
diff --git a/tests/unit/test_expression_tree/test_variable.py b/tests/unit/test_expression_tree/test_variable.py
index 900cccb78f..62c5a9cf62 100644
--- a/tests/unit/test_expression_tree/test_variable.py
+++ b/tests/unit/test_expression_tree/test_variable.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Variable class
 #
-from tests import TestCase
+
 import unittest
 
 import numpy as np
@@ -10,7 +10,7 @@
 import sympy
 
 
-class TestVariable(TestCase):
+class TestVariable(unittest.TestCase):
     def test_variable_init(self):
         a = pybamm.Variable("a")
         self.assertEqual(a.name, "a")
@@ -70,7 +70,7 @@ def test_to_json_error(self):
             func.to_json()
 
 
-class TestVariableDot(TestCase):
+class TestVariableDot(unittest.TestCase):
     def test_variable_init(self):
         a = pybamm.VariableDot("a'")
         self.assertEqual(a.name, "a'")
diff --git a/tests/unit/test_geometry/test_battery_geometry.py b/tests/unit/test_geometry/test_battery_geometry.py
index 57d638e5ef..38e1ce1908 100644
--- a/tests/unit/test_geometry/test_battery_geometry.py
+++ b/tests/unit/test_geometry/test_battery_geometry.py
@@ -1,12 +1,12 @@
 #
 # Tests for the base model class
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestBatteryGeometry(TestCase):
+class TestBatteryGeometry(unittest.TestCase):
     def test_geometry_keys(self):
         for cc_dimension in [0, 1, 2]:
             geometry = pybamm.battery_geometry(
@@ -84,7 +84,7 @@ def test_geometry_error(self):
             pybamm.battery_geometry(form_factor="triangle")
 
 
-class TestReadParameters(TestCase):
+class TestReadParameters(unittest.TestCase):
     # This is the most complicated geometry and should test the parameters are
     # all returned for the deepest dict
     def test_read_parameters(self):
diff --git a/tests/unit/test_meshes/test_meshes.py b/tests/unit/test_meshes/test_meshes.py
index 3066d14534..2f3bffddfb 100644
--- a/tests/unit/test_meshes/test_meshes.py
+++ b/tests/unit/test_meshes/test_meshes.py
@@ -1,7 +1,7 @@
 #
 # Test for the Finite Volume Mesh class
 #
-from tests import TestCase
+
 import pybamm
 import numpy as np
 import unittest
@@ -19,7 +19,7 @@ def get_param():
     )
 
 
-class TestMesh(TestCase):
+class TestMesh(unittest.TestCase):
     def test_mesh_creation_no_parameters(self):
         r = pybamm.SpatialVariable(
             "r", domain=["negative particle"], coord_sys="spherical polar"
@@ -415,7 +415,7 @@ def test_to_json(self):
         self.assertEqual(mesh_json, expected_json)
 
 
-class TestMeshGenerator(TestCase):
+class TestMeshGenerator(unittest.TestCase):
     def test_init_name(self):
         mesh_generator = pybamm.MeshGenerator(pybamm.SubMesh0D)
         self.assertEqual(mesh_generator.__repr__(), "Generator for SubMesh0D")
diff --git a/tests/unit/test_meshes/test_one_dimensional_submesh.py b/tests/unit/test_meshes/test_one_dimensional_submesh.py
index 514de4248b..82429e475c 100644
--- a/tests/unit/test_meshes/test_one_dimensional_submesh.py
+++ b/tests/unit/test_meshes/test_one_dimensional_submesh.py
@@ -1,10 +1,9 @@
 import pybamm
 import unittest
 import numpy as np
-from tests import TestCase
 
 
-class TestSubMesh1D(TestCase):
+class TestSubMesh1D(unittest.TestCase):
     def test_tabs(self):
         edges = np.linspace(0, 1, 10)
         tabs = {"negative": {"z_centre": 0}, "positive": {"z_centre": 1}}
@@ -49,7 +48,7 @@ def test_to_json(self):
         self.assertEqual(mesh.tabs, new_mesh.tabs)
 
 
-class TestUniform1DSubMesh(TestCase):
+class TestUniform1DSubMesh(unittest.TestCase):
     def test_exceptions(self):
         lims = {"a": 1, "b": 2}
         with self.assertRaises(pybamm.GeometryError):
@@ -82,7 +81,7 @@ def test_symmetric_mesh_creation_no_parameters(self):
         )
 
 
-class TestExponential1DSubMesh(TestCase):
+class TestExponential1DSubMesh(unittest.TestCase):
     def test_symmetric_mesh_creation_no_parameters_even(self):
         r = pybamm.SpatialVariable(
             "r", domain=["negative particle"], coord_sys="spherical polar"
@@ -208,7 +207,7 @@ def test_right_mesh_creation_no_parameters(self):
         )
 
 
-class TestChebyshev1DSubMesh(TestCase):
+class TestChebyshev1DSubMesh(unittest.TestCase):
     def test_mesh_creation_no_parameters(self):
         r = pybamm.SpatialVariable(
             "r", domain=["negative particle"], coord_sys="spherical polar"
@@ -236,7 +235,7 @@ def test_mesh_creation_no_parameters(self):
         )
 
 
-class TestUser1DSubMesh(TestCase):
+class TestUser1DSubMesh(unittest.TestCase):
     def test_exceptions(self):
         edges = np.array([0, 0.3, 1])
         submesh_params = {"edges": edges}
@@ -298,7 +297,7 @@ def test_mesh_creation_no_parameters(self):
         )
 
 
-class TestSpectralVolume1DSubMesh(TestCase):
+class TestSpectralVolume1DSubMesh(unittest.TestCase):
     def test_exceptions(self):
         edges = np.array([0, 0.3, 1])
         submesh_params = {"edges": edges}
diff --git a/tests/unit/test_meshes/test_scikit_fem_submesh.py b/tests/unit/test_meshes/test_scikit_fem_submesh.py
index 83c0192d30..07e7dd016a 100644
--- a/tests/unit/test_meshes/test_scikit_fem_submesh.py
+++ b/tests/unit/test_meshes/test_scikit_fem_submesh.py
@@ -1,7 +1,7 @@
 #
 # Test for the scikit-fem Finite Element Mesh class
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 import numpy as np
@@ -25,7 +25,7 @@ def get_param():
     )
 
 
-class TestScikitFiniteElement2DSubMesh(TestCase):
+class TestScikitFiniteElement2DSubMesh(unittest.TestCase):
     def test_mesh_creation(self):
         param = get_param()
         geometry = pybamm.battery_geometry(
@@ -284,7 +284,7 @@ def test_to_json(self):
             np.testing.assert_array_equal(x, y)
 
 
-class TestScikitFiniteElementChebyshev2DSubMesh(TestCase):
+class TestScikitFiniteElementChebyshev2DSubMesh(unittest.TestCase):
     def test_mesh_creation(self):
         param = get_param()
 
@@ -347,7 +347,7 @@ def test_init_failure(self):
             pybamm.ScikitChebyshev2DSubMesh(lims, None)
 
 
-class TestScikitExponential2DSubMesh(TestCase):
+class TestScikitExponential2DSubMesh(unittest.TestCase):
     def test_mesh_creation(self):
         param = get_param()
 
@@ -416,7 +416,7 @@ def test_init_failure(self):
             pybamm.ScikitExponential2DSubMesh(None, None, side="bottom")
 
 
-class TestScikitUser2DSubMesh(TestCase):
+class TestScikitUser2DSubMesh(unittest.TestCase):
     def test_mesh_creation(self):
         param = get_param()
 
diff --git a/tests/unit/test_models/test_base_model.py b/tests/unit/test_models/test_base_model.py
index 7caa4c94b8..4d5f71201a 100644
--- a/tests/unit/test_models/test_base_model.py
+++ b/tests/unit/test_models/test_base_model.py
@@ -1,7 +1,7 @@
 #
 # Tests for the base model class
 #
-from tests import TestCase
+
 import os
 import platform
 import subprocess  # nosec
@@ -16,7 +16,7 @@
 import pybamm
 
 
-class TestBaseModel(TestCase):
+class TestBaseModel(unittest.TestCase):
     def test_rhs_set_get(self):
         model = pybamm.BaseModel()
         rhs = {
diff --git a/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py b/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py
index 25d142d60e..af95a4da35 100644
--- a/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py
+++ b/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py
@@ -1,7 +1,7 @@
 #
 # Tests for the base battery model class
 #
-from tests import TestCase
+
 from pybamm.models.full_battery_models.base_battery_model import BatteryModelOptions
 import pybamm
 import unittest
@@ -55,7 +55,7 @@
 """
 
 
-class TestBaseBatteryModel(TestCase):
+class TestBaseBatteryModel(unittest.TestCase):
     def test_process_parameters_and_discretise(self):
         model = pybamm.lithium_ion.SPM()
         # Set up geometry and parameters
@@ -489,7 +489,7 @@ def test_save_load_model(self):
         os.remove("test_base_battery_model.json")
 
 
-class TestOptions(TestCase):
+class TestOptions(unittest.TestCase):
     def test_print_options(self):
         with io.StringIO() as buffer, redirect_stdout(buffer):
             BatteryModelOptions(OPTIONS_DICT).print_options()
diff --git a/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py b/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py
index 95365fb42b..f6abf592f8 100644
--- a/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py
+++ b/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py
@@ -1,12 +1,12 @@
 #
 # Tests for the Thevenin equivalant circuit model
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestThevenin(TestCase):
+class TestThevenin(unittest.TestCase):
     def test_standard_model(self):
         model = pybamm.equivalent_circuit.Thevenin()
         model.check_well_posedness()
diff --git a/tests/unit/test_models/test_full_battery_models/test_lead_acid/test_loqs.py b/tests/unit/test_models/test_full_battery_models/test_lead_acid/test_loqs.py
index c20626856b..d68686936c 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lead_acid/test_loqs.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lead_acid/test_loqs.py
@@ -1,12 +1,12 @@
 #
 # Tests for the lead-acid LOQS model
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestLeadAcidLOQS(TestCase):
+class TestLeadAcidLOQS(unittest.TestCase):
     def test_well_posed(self):
         options = {"thermal": "isothermal"}
         model = pybamm.lead_acid.LOQS(options)
@@ -80,7 +80,7 @@ def test_well_posed_2plus1D(self):
         )
 
 
-class TestLeadAcidLOQSWithSideReactions(TestCase):
+class TestLeadAcidLOQSWithSideReactions(unittest.TestCase):
     def test_well_posed_differential(self):
         options = {"surface form": "differential", "hydrolysis": "true"}
         model = pybamm.lead_acid.LOQS(options)
@@ -92,7 +92,7 @@ def test_well_posed_algebraic(self):
         model.check_well_posedness()
 
 
-class TestLeadAcidLOQSSurfaceForm(TestCase):
+class TestLeadAcidLOQSSurfaceForm(unittest.TestCase):
     def test_well_posed_differential(self):
         options = {"surface form": "differential"}
         model = pybamm.lead_acid.LOQS(options)
@@ -121,7 +121,7 @@ def test_default_geometry(self):
         self.assertIn("current collector", model.default_geometry)
 
 
-class TestLeadAcidLOQSExternalCircuits(TestCase):
+class TestLeadAcidLOQSExternalCircuits(unittest.TestCase):
     def test_well_posed_voltage(self):
         options = {"operating mode": "voltage"}
         model = pybamm.lead_acid.LOQS(options)
diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py
index dd7d35b683..9f044b0566 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_electrode_soh.py
@@ -1,12 +1,12 @@
 #
 # Tests for the lithium-ion electrode-specific SOH model
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestElectrodeSOH(TestCase):
+class TestElectrodeSOH(unittest.TestCase):
     def test_known_solution(self):
         param = pybamm.LithiumIonParameters()
         parameter_values = pybamm.ParameterValues("Mohtat2020")
@@ -146,7 +146,7 @@ def test_error(self):
             esoh_solver.solve(inputs)
 
 
-class TestElectrodeSOHMSMR(TestCase):
+class TestElectrodeSOHMSMR(unittest.TestCase):
     def test_known_solution(self):
         options = {
             "open-circuit potential": "MSMR",
@@ -232,7 +232,7 @@ def test_error(self):
             esoh_solver._get_electrode_soh_sims_split()
 
 
-class TestElectrodeSOHHalfCell(TestCase):
+class TestElectrodeSOHHalfCell(unittest.TestCase):
     def test_known_solution(self):
         model = pybamm.lithium_ion.ElectrodeSOHHalfCell()
         param = pybamm.LithiumIonParameters({"working electrode": "positive"})
@@ -247,7 +247,7 @@ def test_known_solution(self):
         self.assertAlmostEqual(sol["Uw(x_0)"].data[0], V_min, places=5)
 
 
-class TestCalculateTheoreticalEnergy(TestCase):
+class TestCalculateTheoreticalEnergy(unittest.TestCase):
     def test_efficiency(self):
         model = pybamm.lithium_ion.DFN(options={"calculate discharge energy": "true"})
         parameter_values = pybamm.ParameterValues("Chen2020")
@@ -266,7 +266,7 @@ def test_efficiency(self):
         self.assertLess(0, theoretical_energy)
 
 
-class TestGetInitialSOC(TestCase):
+class TestGetInitialSOC(unittest.TestCase):
     def test_initial_soc(self):
         param = pybamm.LithiumIonParameters()
         parameter_values = pybamm.ParameterValues("Mohtat2020")
@@ -390,7 +390,7 @@ def test_error(self):
             )
 
 
-class TestGetInitialOCP(TestCase):
+class TestGetInitialOCP(unittest.TestCase):
     def test_get_initial_ocp(self):
         param = pybamm.LithiumIonParameters()
         parameter_values = pybamm.ParameterValues("Mohtat2020")
@@ -412,7 +412,7 @@ def test_min_max_ocp(self):
         self.assertAlmostEqual(Up_0 - Un_0, 2.8)
 
 
-class TestGetInitialOCPMSMR(TestCase):
+class TestGetInitialOCPMSMR(unittest.TestCase):
     def test_get_initial_ocp(self):
         options = {
             "open-circuit potential": "MSMR",
diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py
index 389aa55849..e5147f01e2 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_mpm.py
@@ -1,12 +1,12 @@
 #
 # Tests for the lithium-ion MPM model
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestMPM(TestCase):
+class TestMPM(unittest.TestCase):
     def test_well_posed(self):
         options = {"thermal": "isothermal"}
         model = pybamm.lithium_ion.MPM(options)
@@ -124,7 +124,7 @@ def test_wycisk_ocp(self):
         model.check_well_posedness()
 
 
-class TestMPMExternalCircuits(TestCase):
+class TestMPMExternalCircuits(unittest.TestCase):
     def test_well_posed_voltage(self):
         options = {"operating mode": "voltage"}
         model = pybamm.lithium_ion.MPM(options)
@@ -146,7 +146,7 @@ def external_circuit_function(variables):
         model.check_well_posedness()
 
 
-class TestMPMWithSEI(TestCase):
+class TestMPMWithSEI(unittest.TestCase):
     def test_reaction_limited_not_implemented(self):
         options = {"SEI": "reaction limited"}
         with self.assertRaises(NotImplementedError):
@@ -176,7 +176,7 @@ def test_ec_reaction_limited_not_implemented(self):
             pybamm.lithium_ion.MPM(options)
 
 
-class TestMPMWithMechanics(TestCase):
+class TestMPMWithMechanics(unittest.TestCase):
     def test_well_posed_negative_cracking_not_implemented(self):
         options = {"particle mechanics": ("swelling and cracking", "none")}
         with self.assertRaises(NotImplementedError):
@@ -198,7 +198,7 @@ def test_well_posed_both_swelling_only_not_implemented(self):
             pybamm.lithium_ion.MPM(options)
 
 
-class TestMPMWithPlating(TestCase):
+class TestMPMWithPlating(unittest.TestCase):
     def test_well_posed_reversible_plating_not_implemented(self):
         options = {"lithium plating": "reversible"}
         with self.assertRaises(NotImplementedError):
diff --git a/tests/unit/test_parameters/test_bpx.py b/tests/unit/test_parameters/test_bpx.py
index 916eb8d161..f57fe8f7fa 100644
--- a/tests/unit/test_parameters/test_bpx.py
+++ b/tests/unit/test_parameters/test_bpx.py
@@ -1,7 +1,7 @@
 #
 # Tests for the create_from_bpx function
 #
-from tests import TestCase
+
 
 import tempfile
 import unittest
@@ -12,7 +12,7 @@
 import pytest
 
 
-class TestBPX(TestCase):
+class TestBPX(unittest.TestCase):
     def setUp(self):
         self.base = {
             "Header": {
diff --git a/tests/unit/test_parameters/test_current_functions.py b/tests/unit/test_parameters/test_current_functions.py
index 8a8cc266ce..b00cba0b89 100644
--- a/tests/unit/test_parameters/test_current_functions.py
+++ b/tests/unit/test_parameters/test_current_functions.py
@@ -1,7 +1,7 @@
 #
 # Tests for current input functions
 #
-from tests import TestCase
+
 import pybamm
 import numbers
 import unittest
@@ -11,7 +11,7 @@
 from tests import no_internet_connection
 
 
-class TestCurrentFunctions(TestCase):
+class TestCurrentFunctions(unittest.TestCase):
     def test_constant_current(self):
         # test simplify
         param = pybamm.electrical_parameters
diff --git a/tests/unit/test_parameters/test_ecm_parameters.py b/tests/unit/test_parameters/test_ecm_parameters.py
index 4764409e4e..543b4f4e5b 100644
--- a/tests/unit/test_parameters/test_ecm_parameters.py
+++ b/tests/unit/test_parameters/test_ecm_parameters.py
@@ -1,7 +1,7 @@
 #
 # Tests for the equivalent circuit parameters
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
@@ -33,7 +33,7 @@
 parameter_values = pybamm.ParameterValues(values)
 
 
-class TestEcmParameters(TestCase):
+class TestEcmParameters(unittest.TestCase):
     def test_init_parameters(self):
         param = pybamm.EcmParameters()
 
diff --git a/tests/unit/test_parameters/test_lead_acid_parameters.py b/tests/unit/test_parameters/test_lead_acid_parameters.py
index ddc73f61ee..3fc62fde93 100644
--- a/tests/unit/test_parameters/test_lead_acid_parameters.py
+++ b/tests/unit/test_parameters/test_lead_acid_parameters.py
@@ -2,14 +2,14 @@
 # Test for the standard lead acid parameters
 #
 import os
-from tests import TestCase
+
 import pybamm
 from tests import get_discretisation_for_testing
 from tempfile import TemporaryDirectory
 import unittest
 
 
-class TestStandardParametersLeadAcid(TestCase):
+class TestStandardParametersLeadAcid(unittest.TestCase):
     def test_scipy_constants(self):
         constants = pybamm.constants
         self.assertAlmostEqual(constants.R.evaluate(), 8.314, places=3)
diff --git a/tests/unit/test_parameters/test_lithium_ion_parameters.py b/tests/unit/test_parameters/test_lithium_ion_parameters.py
index 0c46eec16e..66c4ea398e 100644
--- a/tests/unit/test_parameters/test_lithium_ion_parameters.py
+++ b/tests/unit/test_parameters/test_lithium_ion_parameters.py
@@ -2,14 +2,14 @@
 # Tests lithium-ion parameters load and give expected values
 #
 import os
-from tests import TestCase
+
 import pybamm
 from tempfile import TemporaryDirectory
 import unittest
 import numpy as np
 
 
-class TestLithiumIonParameterValues(TestCase):
+class TestLithiumIonParameterValues(unittest.TestCase):
     def test_print_parameters(self):
         with TemporaryDirectory() as dir_name:
             parameters = pybamm.LithiumIonParameters()
diff --git a/tests/unit/test_parameters/test_parameter_sets/test_Ecker2015.py b/tests/unit/test_parameters/test_parameter_sets/test_Ecker2015.py
index 2dc73d4484..4be67175d7 100644
--- a/tests/unit/test_parameters/test_parameter_sets/test_Ecker2015.py
+++ b/tests/unit/test_parameters/test_parameter_sets/test_Ecker2015.py
@@ -1,12 +1,12 @@
 #
 # Tests for O'Kane (2022) parameter set
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestEcker2015(TestCase):
+class TestEcker2015(unittest.TestCase):
     def test_functions(self):
         param = pybamm.ParameterValues("Ecker2015")
         sto = pybamm.Scalar(0.5)
diff --git a/tests/unit/test_parameters/test_parameter_sets/test_OKane2022.py b/tests/unit/test_parameters/test_parameter_sets/test_OKane2022.py
index 1ab7e7930e..e34f837b38 100644
--- a/tests/unit/test_parameters/test_parameter_sets/test_OKane2022.py
+++ b/tests/unit/test_parameters/test_parameter_sets/test_OKane2022.py
@@ -1,12 +1,12 @@
 #
 # Tests for O'Kane (2022) parameter set
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestOKane2022(TestCase):
+class TestOKane2022(unittest.TestCase):
     def test_functions(self):
         param = pybamm.ParameterValues("OKane2022")
         sto = pybamm.Scalar(0.9)
diff --git a/tests/unit/test_parameters/test_parameter_sets/test_parameters_with_default_models.py b/tests/unit/test_parameters/test_parameter_sets/test_parameters_with_default_models.py
index 4e18f1ef50..d7133a73e0 100644
--- a/tests/unit/test_parameters/test_parameter_sets/test_parameters_with_default_models.py
+++ b/tests/unit/test_parameters/test_parameter_sets/test_parameters_with_default_models.py
@@ -1,12 +1,12 @@
 #
 # Tests each parameter set with the standard model associated with that parameter set
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 
 
-class TestParameterValuesWithModel(TestCase):
+class TestParameterValuesWithModel(unittest.TestCase):
     def test_parameter_values_with_model(self):
         param_to_model = {
             "Ai2020": pybamm.lithium_ion.DFN(
diff --git a/tests/unit/test_parameters/test_parameter_values.py b/tests/unit/test_parameters/test_parameter_values.py
index 28b8aa2ef9..4826885f1c 100644
--- a/tests/unit/test_parameters/test_parameter_values.py
+++ b/tests/unit/test_parameters/test_parameter_values.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Base Parameter Values class
 #
-from tests import TestCase
+
 
 import os
 import unittest
@@ -19,7 +19,7 @@
 import casadi
 
 
-class TestParameterValues(TestCase):
+class TestParameterValues(unittest.TestCase):
     def test_init(self):
         # from dict
         param = pybamm.ParameterValues({"a": 1})
diff --git a/tests/unit/test_parameters/test_process_parameter_data.py b/tests/unit/test_parameters/test_process_parameter_data.py
index 22f83a6bdb..bf27f71e4f 100644
--- a/tests/unit/test_parameters/test_process_parameter_data.py
+++ b/tests/unit/test_parameters/test_process_parameter_data.py
@@ -1,7 +1,7 @@
 #
 # Tests for the parameter processing functions
 #
-from tests import TestCase
+
 
 import os
 import numpy as np
@@ -10,7 +10,7 @@
 import unittest
 
 
-class TestProcessParameterData(TestCase):
+class TestProcessParameterData(unittest.TestCase):
     def test_process_1D_data(self):
         name = "lico2_ocv_example"
         path = os.path.join(pybamm.root_dir(), "tests", "unit", "test_parameters")
diff --git a/tests/unit/test_plotting/test_quick_plot.py b/tests/unit/test_plotting/test_quick_plot.py
index e9e5dd810b..eb6c0607e3 100644
--- a/tests/unit/test_plotting/test_quick_plot.py
+++ b/tests/unit/test_plotting/test_quick_plot.py
@@ -1,12 +1,12 @@
 import os
 import pybamm
 import unittest
-from tests import TestCase
+
 import numpy as np
 from tempfile import TemporaryDirectory
 
 
-class TestQuickPlot(TestCase):
+class TestQuickPlot(unittest.TestCase):
     def test_simple_ode_model(self):
         model = pybamm.lithium_ion.BaseModel(name="Simple ODE Model")
 
diff --git a/tests/unit/test_serialisation/test_serialisation.py b/tests/unit/test_serialisation/test_serialisation.py
index e7dcba6702..a1286cad26 100644
--- a/tests/unit/test_serialisation/test_serialisation.py
+++ b/tests/unit/test_serialisation/test_serialisation.py
@@ -1,7 +1,7 @@
 #
 # Tests for the serialisation class
 #
-from tests import TestCase
+
 import json
 import os
 import unittest
@@ -86,7 +86,7 @@ def mesh_var_dict():
     return mesh, mesh_json
 
 
-class TestSerialiseModels(TestCase):
+class TestSerialiseModels(unittest.TestCase):
     def test_user_defined_model_recreaction(self):
         # Start with a base model
         model = pybamm.BaseModel()
@@ -146,7 +146,7 @@ def test_user_defined_model_recreaction(self):
         os.remove("heat_equation.json")
 
 
-class TestSerialise(TestCase):
+class TestSerialise(unittest.TestCase):
     # test the symbol encoder
 
     def test_symbol_encoder_symbol(self):
diff --git a/tests/unit/test_simulation.py b/tests/unit/test_simulation.py
index 368c84607c..744ea2457c 100644
--- a/tests/unit/test_simulation.py
+++ b/tests/unit/test_simulation.py
@@ -1,7 +1,7 @@
 import pybamm
 import numpy as np
 import pandas as pd
-from tests import TestCase
+
 import os
 import sys
 import unittest
@@ -12,7 +12,7 @@
 from tests import no_internet_connection
 
 
-class TestSimulation(TestCase):
+class TestSimulation(unittest.TestCase):
     def test_simple_model(self):
         model = pybamm.BaseModel()
         v = pybamm.Variable("v")
diff --git a/tests/unit/test_solvers/test_algebraic_solver.py b/tests/unit/test_solvers/test_algebraic_solver.py
index 22017f8cf8..6e8b3a3d80 100644
--- a/tests/unit/test_solvers/test_algebraic_solver.py
+++ b/tests/unit/test_solvers/test_algebraic_solver.py
@@ -1,14 +1,14 @@
 #
 # Tests for the Algebraic Solver class
 #
-from tests import TestCase
+
 import pybamm
 import unittest
 import numpy as np
 from tests import get_discretisation_for_testing
 
 
-class TestAlgebraicSolver(TestCase):
+class TestAlgebraicSolver(unittest.TestCase):
     def test_algebraic_solver_init(self):
         solver = pybamm.AlgebraicSolver(
             method="hybr", tol=1e-4, extra_options={"maxfev": 100}
diff --git a/tests/unit/test_solvers/test_base_solver.py b/tests/unit/test_solvers/test_base_solver.py
index c444722929..884c85f87f 100644
--- a/tests/unit/test_solvers/test_base_solver.py
+++ b/tests/unit/test_solvers/test_base_solver.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Base Solver class
 #
-from tests import TestCase
+
 import casadi
 import pybamm
 import numpy as np
@@ -10,7 +10,7 @@
 import unittest
 
 
-class TestBaseSolver(TestCase):
+class TestBaseSolver(unittest.TestCase):
     def test_base_solver_init(self):
         solver = pybamm.BaseSolver(rtol=1e-2, atol=1e-4)
         self.assertEqual(solver.rtol, 1e-2)
diff --git a/tests/unit/test_solvers/test_casadi_algebraic_solver.py b/tests/unit/test_solvers/test_casadi_algebraic_solver.py
index 5001b37d82..b85f4292b9 100644
--- a/tests/unit/test_solvers/test_casadi_algebraic_solver.py
+++ b/tests/unit/test_solvers/test_casadi_algebraic_solver.py
@@ -1,4 +1,3 @@
-from tests import TestCase
 import casadi
 import pybamm
 import unittest
@@ -7,7 +6,7 @@
 import tests
 
 
-class TestCasadiAlgebraicSolver(TestCase):
+class TestCasadiAlgebraicSolver(unittest.TestCase):
     def test_algebraic_solver_init(self):
         solver = pybamm.CasadiAlgebraicSolver(tol=1e-4)
         self.assertEqual(solver.tol, 1e-4)
@@ -171,7 +170,7 @@ def test_solve_with_input(self):
         np.testing.assert_array_equal(solution.y, -7)
 
 
-class TestCasadiAlgebraicSolverSensitivity(TestCase):
+class TestCasadiAlgebraicSolverSensitivity(unittest.TestCase):
     def test_solve_with_symbolic_input(self):
         # Simple system: a single algebraic equation
         var = pybamm.Variable("var")
diff --git a/tests/unit/test_solvers/test_casadi_solver.py b/tests/unit/test_solvers/test_casadi_solver.py
index c798024579..5ce29a365d 100644
--- a/tests/unit/test_solvers/test_casadi_solver.py
+++ b/tests/unit/test_solvers/test_casadi_solver.py
@@ -1,4 +1,3 @@
-from tests import TestCase
 import pybamm
 import unittest
 import numpy as np
@@ -6,7 +5,7 @@
 from scipy.sparse import eye
 
 
-class TestCasadiSolver(TestCase):
+class TestCasadiSolver(unittest.TestCase):
     def test_bad_mode(self):
         with self.assertRaisesRegex(ValueError, "invalid mode"):
             pybamm.CasadiSolver(mode="bad mode")
@@ -582,7 +581,7 @@ def test_modulo_non_smooth_events(self):
         )
 
 
-class TestCasadiSolverODEsWithForwardSensitivityEquations(TestCase):
+class TestCasadiSolverODEsWithForwardSensitivityEquations(unittest.TestCase):
     def test_solve_sensitivity_scalar_var_scalar_input(self):
         # Create model
         model = pybamm.BaseModel()
@@ -963,7 +962,7 @@ def test_solve_sensitivity_subset(self):
         )
 
 
-class TestCasadiSolverDAEsWithForwardSensitivityEquations(TestCase):
+class TestCasadiSolverDAEsWithForwardSensitivityEquations(unittest.TestCase):
     def test_solve_sensitivity_scalar_var_scalar_input(self):
         # Create model
         model = pybamm.BaseModel()
diff --git a/tests/unit/test_solvers/test_idaklu_jax.py b/tests/unit/test_solvers/test_idaklu_jax.py
index 6d891010d6..7bae5d74e9 100644
--- a/tests/unit/test_solvers/test_idaklu_jax.py
+++ b/tests/unit/test_solvers/test_idaklu_jax.py
@@ -1,7 +1,7 @@
 #
 # Tests for the KLU-Jax interface class
 #
-from tests import TestCase
+
 from parameterized import parameterized
 
 import pybamm
@@ -87,7 +87,7 @@ def no_jit(f):
     pybamm.have_idaklu() and pybamm.have_jax(),
     "Both IDAKLU and JAX are available",
 )
-class TestIDAKLUJax_NoJax(TestCase):
+class TestIDAKLUJax_NoJax(unittest.TestCase):
     def test_instantiate_fails(self):
         with self.assertRaises(ModuleNotFoundError):
             pybamm.IDAKLUJax([], [], [])
@@ -97,7 +97,7 @@ def test_instantiate_fails(self):
     not pybamm.have_idaklu() or not pybamm.have_jax(),
     "IDAKLU Solver and/or JAX are not available",
 )
-class TestIDAKLUJax(TestCase):
+class TestIDAKLUJax(unittest.TestCase):
     # Initialisation tests
 
     def test_initialise_twice(self):
diff --git a/tests/unit/test_solvers/test_idaklu_solver.py b/tests/unit/test_solvers/test_idaklu_solver.py
index 02303a364e..c14c29fbb6 100644
--- a/tests/unit/test_solvers/test_idaklu_solver.py
+++ b/tests/unit/test_solvers/test_idaklu_solver.py
@@ -1,7 +1,7 @@
 #
 # Tests for the KLU Solver class
 #
-from tests import TestCase
+
 from contextlib import redirect_stdout
 import io
 import unittest
@@ -13,7 +13,7 @@
 
 
 @unittest.skipIf(not pybamm.have_idaklu(), "idaklu solver is not installed")
-class TestIDAKLUSolver(TestCase):
+class TestIDAKLUSolver(unittest.TestCase):
     def test_ida_roberts_klu(self):
         # this test implements a python version of the ida Roberts
         # example provided in sundials
diff --git a/tests/unit/test_solvers/test_jax_bdf_solver.py b/tests/unit/test_solvers/test_jax_bdf_solver.py
index 854a618fba..e02bdb2510 100644
--- a/tests/unit/test_solvers/test_jax_bdf_solver.py
+++ b/tests/unit/test_solvers/test_jax_bdf_solver.py
@@ -1,7 +1,7 @@
 import pybamm
 import unittest
 from tests import get_mesh_for_testing
-from tests import TestCase
+
 import sys
 import numpy as np
 
@@ -10,7 +10,7 @@
 
 
 @unittest.skipIf(not pybamm.have_jax(), "jax or jaxlib is not installed")
-class TestJaxBDFSolver(TestCase):
+class TestJaxBDFSolver(unittest.TestCase):
     def test_solver_(self):  # Trailing _ manipulates the random seed
         # Create model
         model = pybamm.BaseModel()
diff --git a/tests/unit/test_solvers/test_jax_solver.py b/tests/unit/test_solvers/test_jax_solver.py
index 9df28e8ac2..4f34497626 100644
--- a/tests/unit/test_solvers/test_jax_solver.py
+++ b/tests/unit/test_solvers/test_jax_solver.py
@@ -1,7 +1,7 @@
 import pybamm
 import unittest
 from tests import get_mesh_for_testing
-from tests import TestCase
+
 import sys
 import numpy as np
 
@@ -10,7 +10,7 @@
 
 
 @unittest.skipIf(not pybamm.have_jax(), "jax or jaxlib is not installed")
-class TestJaxSolver(TestCase):
+class TestJaxSolver(unittest.TestCase):
     def test_model_solver(self):
         # Create model
         model = pybamm.BaseModel()
diff --git a/tests/unit/test_solvers/test_processed_variable.py b/tests/unit/test_solvers/test_processed_variable.py
index 4cf3f9392e..b6ae669878 100644
--- a/tests/unit/test_solvers/test_processed_variable.py
+++ b/tests/unit/test_solvers/test_processed_variable.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Processed Variable class
 #
-from tests import TestCase
+
 import casadi
 import pybamm
 import tests
@@ -61,7 +61,7 @@ def process_and_check_2D_variable(
     return y_sol, first_sol, second_sol, t_sol
 
 
-class TestProcessedVariable(TestCase):
+class TestProcessedVariable(unittest.TestCase):
     def test_processed_variable_0D(self):
         # without space
         t = pybamm.t
diff --git a/tests/unit/test_solvers/test_processed_variable_computed.py b/tests/unit/test_solvers/test_processed_variable_computed.py
index 7e0616c81b..407d422e4c 100644
--- a/tests/unit/test_solvers/test_processed_variable_computed.py
+++ b/tests/unit/test_solvers/test_processed_variable_computed.py
@@ -5,7 +5,7 @@
 #  by the idaklu solver, and does not possesses any capability to calculate
 #  values itself since it does not have access to the full state vector
 #
-from tests import TestCase
+
 import casadi
 import pybamm
 import tests
@@ -68,7 +68,7 @@ def process_and_check_2D_variable(
     return y_sol, first_sol, second_sol, t_sol
 
 
-class TestProcessedVariableComputed(TestCase):
+class TestProcessedVariableComputed(unittest.TestCase):
     def test_processed_variable_0D(self):
         # without space
         y = pybamm.StateVector(slice(0, 1))
diff --git a/tests/unit/test_solvers/test_scipy_solver.py b/tests/unit/test_solvers/test_scipy_solver.py
index fad6651d55..c6afd16704 100644
--- a/tests/unit/test_solvers/test_scipy_solver.py
+++ b/tests/unit/test_solvers/test_scipy_solver.py
@@ -1,7 +1,7 @@
 # Tests for the Scipy Solver class
 #
 import pybamm
-from tests import TestCase
+
 import unittest
 import numpy as np
 from tests import get_mesh_for_testing, get_discretisation_for_testing
@@ -9,7 +9,7 @@
 import sys
 
 
-class TestScipySolver(TestCase):
+class TestScipySolver(unittest.TestCase):
     def test_model_solver_python_and_jax(self):
         if pybamm.have_jax():
             formats = ["python", "jax"]
@@ -492,7 +492,7 @@ def test_scale_and_reference(self):
         )
 
 
-class TestScipySolverWithSensitivity(TestCase):
+class TestScipySolverWithSensitivity(unittest.TestCase):
     def test_solve_sensitivity_scalar_var_scalar_input(self):
         # Create model
         model = pybamm.BaseModel()
diff --git a/tests/unit/test_solvers/test_solution.py b/tests/unit/test_solvers/test_solution.py
index ecc8d1bd8e..995898e8dd 100644
--- a/tests/unit/test_solvers/test_solution.py
+++ b/tests/unit/test_solvers/test_solution.py
@@ -2,7 +2,7 @@
 # Tests for the Solution class
 #
 import os
-from tests import TestCase
+
 import json
 import pybamm
 import unittest
@@ -13,7 +13,7 @@
 from tempfile import TemporaryDirectory
 
 
-class TestSolution(TestCase):
+class TestSolution(unittest.TestCase):
     def test_init(self):
         t = np.linspace(0, 1)
         y = np.tile(t, (20, 1))
diff --git a/tests/unit/test_spatial_methods/test_base_spatial_method.py b/tests/unit/test_spatial_methods/test_base_spatial_method.py
index d48ea69a7b..647616c924 100644
--- a/tests/unit/test_spatial_methods/test_base_spatial_method.py
+++ b/tests/unit/test_spatial_methods/test_base_spatial_method.py
@@ -1,7 +1,7 @@
 #
 # Test for the base Spatial Method class
 #
-from tests import TestCase
+
 import numpy as np
 import pybamm
 import unittest
@@ -12,7 +12,7 @@
 )
 
 
-class TestSpatialMethod(TestCase):
+class TestSpatialMethod(unittest.TestCase):
     def test_basics(self):
         mesh = get_mesh_for_testing()
         spatial_method = pybamm.SpatialMethod()
diff --git a/tests/unit/test_spatial_methods/test_finite_volume/test_extrapolation.py b/tests/unit/test_spatial_methods/test_finite_volume/test_extrapolation.py
index 86597e6a1c..c51e2d9a13 100644
--- a/tests/unit/test_spatial_methods/test_finite_volume/test_extrapolation.py
+++ b/tests/unit/test_spatial_methods/test_finite_volume/test_extrapolation.py
@@ -1,7 +1,7 @@
 #
 # Test for the extrapolations in the finite volume class
 #
-from tests import TestCase
+
 import pybamm
 from tests import (
     get_mesh_for_testing,
@@ -57,7 +57,7 @@ def get_errors(function, method_options, pts, bcs=None):
     return l_errors, r_errors
 
 
-class TestExtrapolation(TestCase):
+class TestExtrapolation(unittest.TestCase):
     def test_convergence_without_bcs(self):
         # all tests are performed on x in [0, 1]
         linear = {"extrapolation": {"order": "linear"}}
diff --git a/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py b/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py
index 1d5844d7b0..de31b770ff 100644
--- a/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py
+++ b/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py
@@ -1,7 +1,7 @@
 #
 # Tests for the Finite Volume Method
 #
-from tests import TestCase
+
 import pybamm
 from tests import (
     get_mesh_for_testing,
@@ -13,7 +13,7 @@
 import unittest
 
 
-class TestFiniteVolume(TestCase):
+class TestFiniteVolume(unittest.TestCase):
     def test_node_to_edge_to_node(self):
         # Create discretisation
         mesh = get_mesh_for_testing()
diff --git a/tests/unit/test_spatial_methods/test_finite_volume/test_ghost_nodes_and_neumann.py b/tests/unit/test_spatial_methods/test_finite_volume/test_ghost_nodes_and_neumann.py
index 47e1a2fda1..ba82f2fb09 100644
--- a/tests/unit/test_spatial_methods/test_finite_volume/test_ghost_nodes_and_neumann.py
+++ b/tests/unit/test_spatial_methods/test_finite_volume/test_ghost_nodes_and_neumann.py
@@ -1,14 +1,14 @@
 #
 # Test for adding ghost nodes in finite volumes class
 #
-from tests import TestCase
+
 import pybamm
 from tests import get_mesh_for_testing, get_p2d_mesh_for_testing
 import numpy as np
 import unittest
 
 
-class TestGhostNodes(TestCase):
+class TestGhostNodes(unittest.TestCase):
     def test_add_ghost_nodes(self):
         # Set up
 
diff --git a/tests/unit/test_spatial_methods/test_finite_volume/test_grad_div_shapes.py b/tests/unit/test_spatial_methods/test_finite_volume/test_grad_div_shapes.py
index 24132f5cc1..a1dd402f56 100644
--- a/tests/unit/test_spatial_methods/test_finite_volume/test_grad_div_shapes.py
+++ b/tests/unit/test_spatial_methods/test_finite_volume/test_grad_div_shapes.py
@@ -1,7 +1,7 @@
 #
 # Test for the gradient and divergence in Finite Volumes
 #
-from tests import TestCase
+
 import pybamm
 from tests import (
     get_mesh_for_testing,
@@ -13,7 +13,7 @@
 import unittest
 
 
-class TestFiniteVolumeGradDiv(TestCase):
+class TestFiniteVolumeGradDiv(unittest.TestCase):
     def test_grad_div_shapes_Dirichlet_bcs(self):
         """
         Test grad and div with Dirichlet boundary conditions in Cartesian coordinates
diff --git a/tests/unit/test_spatial_methods/test_finite_volume/test_integration.py b/tests/unit/test_spatial_methods/test_finite_volume/test_integration.py
index 57113259c1..e9730a8eb7 100644
--- a/tests/unit/test_spatial_methods/test_finite_volume/test_integration.py
+++ b/tests/unit/test_spatial_methods/test_finite_volume/test_integration.py
@@ -1,7 +1,7 @@
 #
 # Tests for integration using Finite Volume method
 #
-from tests import TestCase
+
 import pybamm
 from tests import (
     get_mesh_for_testing,
@@ -12,7 +12,7 @@
 import unittest
 
 
-class TestFiniteVolumeIntegration(TestCase):
+class TestFiniteVolumeIntegration(unittest.TestCase):
     def test_definite_integral(self):
         # create discretisation
         mesh = get_mesh_for_testing(xpts=200, rpts=200)
diff --git a/tests/unit/test_spatial_methods/test_scikit_finite_element.py b/tests/unit/test_spatial_methods/test_scikit_finite_element.py
index 05b424e053..18c941517b 100644
--- a/tests/unit/test_spatial_methods/test_scikit_finite_element.py
+++ b/tests/unit/test_spatial_methods/test_scikit_finite_element.py
@@ -1,14 +1,14 @@
 #
 # Test for the operator class
 #
-from tests import TestCase
+
 import pybamm
 from tests import get_2p1d_mesh_for_testing, get_unit_2p1D_mesh_for_testing
 import numpy as np
 import unittest
 
 
-class TestScikitFiniteElement(TestCase):
+class TestScikitFiniteElement(unittest.TestCase):
     def test_not_implemented(self):
         mesh = get_2p1d_mesh_for_testing(include_particles=False)
         spatial_method = pybamm.ScikitFiniteElement()
diff --git a/tests/unit/test_spatial_methods/test_spectral_volume.py b/tests/unit/test_spatial_methods/test_spectral_volume.py
index 3988bdd266..f6a631e84c 100644
--- a/tests/unit/test_spatial_methods/test_spectral_volume.py
+++ b/tests/unit/test_spatial_methods/test_spectral_volume.py
@@ -1,7 +1,7 @@
 #
 # Test for the operator class
 #
-from tests import TestCase
+
 import pybamm
 import numpy as np
 import unittest
@@ -87,7 +87,7 @@ def get_1p1d_mesh_for_testing(
     )
 
 
-class TestSpectralVolume(TestCase):
+class TestSpectralVolume(unittest.TestCase):
     def test_exceptions(self):
         sp_meth = pybamm.SpectralVolume()
         with self.assertRaises(ValueError):