Repository: pybamm-team/PyBaMM

Issue:
check that all variables in the model appear in `disc._y_slices`
Otherwise can lead to slightly obscure `KeyError` when discretising

Discussion:

Patch:
diff --git a/pybamm/models/base_models.py b/pybamm/models/base_models.py
index a5e1d9999..233799298 100644
--- a/pybamm/models/base_models.py
+++ b/pybamm/models/base_models.py
@@ -345,7 +345,11 @@ class BaseModel(object):
                     )
 
     def check_variables(self):
-        """ Chec that the right variables are provided. """
+        """
+        Check that the right variables are provided.
+        Check that all variables in the variable list appear in the rhs or algebraic
+        keys
+        """
         missing_vars = []
         for output, expression in self._variables.items():
             if expression is None:
@@ -362,6 +366,31 @@ class BaseModel(object):
             for output in missing_vars:
                 del self._variables[output]
 
+        # Create list of all Variable nodes that appear in the model's list of variables
+        all_vars = {}
+        for eqn in self.variables.values():
+            # Add all variables in the equation to the list of variables
+            all_vars.update(
+                {x.id: x for x in eqn.pre_order() if isinstance(x, pybamm.Variable)}
+            )
+        var_ids_in_keys = set()
+        for var in {**self.rhs, **self.algebraic}.keys():
+            if isinstance(var, pybamm.Variable):
+                var_ids_in_keys.add(var.id)
+            elif isinstance(var, pybamm.Concatenation):
+                var_ids_in_keys.update([x.id for x in var.children])
+        for var_id, var in all_vars.items():
+            if var_id not in var_ids_in_keys:
+                raise pybamm.ModelError(
+                    """
+                    No key set for variable '{}'. Make sure it is included in either
+                    model.rhs or model.algebraic in an unmodified form (e.g. not
+                    Broadcasted)
+                    """.format(
+                        var
+                    )
+                )
+
 
 class StandardBatteryBaseModel(BaseModel):
     """
@@ -417,7 +446,7 @@ class StandardBatteryBaseModel(BaseModel):
                 ),
                 "Positive electrode OCV entropic change": os.path.join(
                     input_path, "lico2_entropic_change_Moura.py"
-                )
+                ),
             },
         )
 
diff --git a/pybamm/solvers/dae_solver.py b/pybamm/solvers/dae_solver.py
index 7223c4a79..4f66b4dcf 100644
--- a/pybamm/solvers/dae_solver.py
+++ b/pybamm/solvers/dae_solver.py
@@ -246,7 +246,7 @@ class DaeSolver(pybamm.BaseSolver):
         # Return full set of consistent initial conditions (y0_diff unchanged)
         y0_consistent = np.concatenate([y0_diff, sol.x])
 
-        if sol.success and np.all(sol.fun < self.root_tol):
+        if sol.success and np.all(sol.fun < self.root_tol * len(sol.x)):
             pybamm.logger.info("Finish calculating consistent initial conditions")
             return y0_consistent
         elif not sol.success:
@@ -255,8 +255,12 @@ class DaeSolver(pybamm.BaseSolver):
             )
         else:
             raise pybamm.SolverError(
-                "Could not find consistent initial conditions: "
-                + "solver terminated successfully, but solution above tolerance"
+                """
+                Could not find consistent initial conditions: solver terminated
+                successfully, but maximum solution error ({}) above tolerance ({})
+                """.format(
+                    np.max(sol.fun), self.root_tol * len(sol.x)
+                )
             )
 
     def integrate(


Test patch:
diff --git a/tests/integration/test_models/standard_output_tests.py b/tests/integration/test_models/standard_output_tests.py
index 22c30be13..1ba521df0 100644
--- a/tests/integration/test_models/standard_output_tests.py
+++ b/tests/integration/test_models/standard_output_tests.py
@@ -363,7 +363,7 @@ class ElectrolyteConcentrationTests(BaseOutputTest):
         during a discharge."""
 
         # TODO: uncomment when have average concentrations
-        # # small number so that can use array less
+        # small number so that can use array less
         # epsilon = 0.001
 
         # if self.operating_condition == "discharge":
diff --git a/tests/integration/test_models/test_reaction_diffusion.py b/tests/integration/test_models/test_reaction_diffusion.py
index 0826f5203..f53e8ab73 100644
--- a/tests/integration/test_models/test_reaction_diffusion.py
+++ b/tests/integration/test_models/test_reaction_diffusion.py
@@ -32,7 +32,7 @@ class TestReactionDiffusionModel(unittest.TestCase):
     def test_convergence(self):
         # Convergence of c at x=0.5
         model = pybamm.ReactionDiffusionModel()
-        # # Make ln and lp nicer for testing
+        # Make ln and lp nicer for testing
         param = model.default_parameter_values
         # Process model and geometry
         param.process_model(model)
diff --git a/tests/unit/test_discretisations/test_discretisation.py b/tests/unit/test_discretisations/test_discretisation.py
index af915a779..964a34010 100644
--- a/tests/unit/test_discretisations/test_discretisation.py
+++ b/tests/unit/test_discretisations/test_discretisation.py
@@ -163,6 +163,7 @@ class TestDiscretise(unittest.TestCase):
         # function of one variable
         def myfun(x):
             return np.exp(x)
+
         func = pybamm.Function(myfun, var)
         func_disc = disc.process_symbol(func)
         self.assertIsInstance(func_disc, pybamm.Function)
@@ -176,6 +177,7 @@ class TestDiscretise(unittest.TestCase):
         # function of multiple variables
         def myfun(x, y):
             return np.exp(x) * y
+
         func = pybamm.Function(myfun, var, scal)
         func_disc = disc.process_symbol(func)
         self.assertIsInstance(func_disc, pybamm.Function)
@@ -752,17 +754,17 @@ class TestDiscretise(unittest.TestCase):
         # check raises error if different sized key and output var
         model.variables = {c_n.name: c_s}
         with self.assertRaisesRegex(pybamm.ModelError, "variable and its eqn"):
-            disc.process_model(model, inplace=False)
+            disc.process_model(model)
 
         # check doesn't raise if concatenation
         model.variables = {c_n.name: pybamm.Concatenation(c_n, c_s)}
-        disc.process_model(model, inplace=False)
+        disc.process_model(model)
 
         # check doesn't raise if broadcast
         model.variables = {
             c_n.name: pybamm.Broadcast(pybamm.Scalar(2), ["negative electrode"])
         }
-        disc.process_model(model, inplace=False)
+        disc.process_model(model)
 
 
 if __name__ == "__main__":
diff --git a/tests/unit/test_expression_tree/test_binary_operators.py b/tests/unit/test_expression_tree/test_binary_operators.py
index ca2390d7d..4b9659228 100644
--- a/tests/unit/test_expression_tree/test_binary_operators.py
+++ b/tests/unit/test_expression_tree/test_binary_operators.py
@@ -116,19 +116,6 @@ class TestBinaryOperators(unittest.TestCase):
             expr2.evaluate(known_evals={})[1][(a @ b).id], (a @ b).evaluate()
         )
 
-        # Expect using known evals to be faster than not
-        timer = pybamm.Timer()
-        start = timer.time()
-        for _ in range(2000):
-            expr2.evaluate()
-        end = timer.time()
-        start_known_evals = timer.time()
-        for _ in range(2000):
-            expr2.evaluate(known_evals={})
-        end_known_evals = timer.time()
-        self.assertLess(end_known_evals - start_known_evals, 1.2 * (end - start))
-        self.assertGreater(end - start, (end_known_evals - start_known_evals))
-
     def test_diff(self):
         a = pybamm.StateVector(slice(0, 1))
         b = pybamm.StateVector(slice(1, 2))
diff --git a/tests/unit/test_expression_tree/test_simplify.py b/tests/unit/test_expression_tree/test_simplify.py
index f58322162..61ebe1657 100644
--- a/tests/unit/test_expression_tree/test_simplify.py
+++ b/tests/unit/test_expression_tree/test_simplify.py
@@ -420,19 +420,6 @@ class TestSimplify(unittest.TestCase):
             expr3.evaluate(y=np.ones(300)), expr3simp.evaluate(y=np.ones(300))
         )
 
-        # we expect simplified solution to be much faster
-        timer = pybamm.Timer()
-        start = timer.time()
-        for _ in range(200):
-            expr3.evaluate(y=np.ones(300))
-        end = timer.time()
-        start_simp = timer.time()
-        for _ in range(200):
-            expr3simp.evaluate(y=np.ones(300))
-        end_simp = timer.time()
-        self.assertLess(end_simp - start_simp, 1.5 * (end - start))
-        self.assertGreater(end - start, (end_simp - start_simp))
-
         m1 = pybamm.Matrix(np.ones((300, 300)))
         m2 = pybamm.Matrix(np.ones((300, 300)))
         m3 = pybamm.Matrix(np.ones((300, 300)))
diff --git a/tests/unit/test_expression_tree/test_symbol.py b/tests/unit/test_expression_tree/test_symbol.py
index dbea2d478..028f0e631 100644
--- a/tests/unit/test_expression_tree/test_symbol.py
+++ b/tests/unit/test_expression_tree/test_symbol.py
@@ -322,8 +322,6 @@ class TestSymbol(unittest.TestCase):
         self.assertEqual(
             (var + broadcast).shape_for_testing, broadcast.shape_for_testing
         )
-        broadcast2 = pybamm.Broadcast(0, domain="separator")
-        self.assertNotEqual(var.shape_for_testing, broadcast2.shape_for_testing)
 
         sym = pybamm.Symbol("sym")
         with self.assertRaises(NotImplementedError):
diff --git a/tests/unit/test_models/test_base_model.py b/tests/unit/test_models/test_base_model.py
index 92361c5fe..98c907122 100644
--- a/tests/unit/test_models/test_base_model.py
+++ b/tests/unit/test_models/test_base_model.py
@@ -335,6 +335,26 @@ class TestBaseModel(unittest.TestCase):
         for key, item in model._variables.items():
             self.assertIsNotNone(item)
 
+        # check error raised if undefined variable in list of Variables
+        model = pybamm.BaseModel()
+        model.rhs = {c: -c}
+        model.initial_conditions = {c: 1}
+        model.variables = {"d": d}
+        with self.assertRaisesRegex(pybamm.ModelError, "No key set for variable"):
+            model.check_well_posedness()
+
+        # check error is raised even if some modified form of d is in model.rhs
+        two_d = 2 * d
+        model.rhs[two_d] = -d
+        model.initial_conditions[two_d] = 1
+        with self.assertRaisesRegex(pybamm.ModelError, "No key set for variable"):
+            model.check_well_posedness()
+
+        # add d to rhs, fine
+        model.rhs[d] = -d
+        model.initial_conditions[d] = 1
+        model.check_well_posedness()
+
 
 class TestStandardBatteryBaseModel(unittest.TestCase):
     def test_default_solver(self):
diff --git a/tests/unit/test_solvers/test_dae_solver.py b/tests/unit/test_solvers/test_dae_solver.py
index fc2e73315..b0e1e2781 100644
--- a/tests/unit/test_solvers/test_dae_solver.py
+++ b/tests/unit/test_solvers/test_dae_solver.py
@@ -43,10 +43,17 @@ class TestDaeSolver(unittest.TestCase):
 
         solver = pybamm.DaeSolver(root_method="hybr")
         y0 = np.array([2])
-        with self.assertRaises(pybamm.SolverError):
+
+        with self.assertRaisesRegex(
+            pybamm.SolverError,
+            "Could not find consistent initial conditions: The iteration is not making",
+        ):
             solver.calculate_consistent_initial_conditions(rhs, algebraic, y0)
         solver = pybamm.DaeSolver()
-        with self.assertRaises(pybamm.SolverError):
+        with self.assertRaisesRegex(
+            pybamm.SolverError,
+            "Could not find consistent initial conditions: solver terminated",
+        ):
             solver.calculate_consistent_initial_conditions(rhs, algebraic, y0)