#1011 first attempt. quite hacky

pybamm/models/base_model.py

@@ -326,11 +326,9 @@ def _find_input_parameters(self):
     def __getitem__(self, key):
         return self.rhs[key]
 
-    def new_copy(self, options=None):
+    def new_copy(self):
         "Create an empty copy with identical options, or new options if specified"
-        options = options or self.options
-        new_model = self.__class__(options)
-        new_model.name = self.name
+        new_model = self.__class__(name=self.name)
         new_model.use_jacobian = self.use_jacobian
         new_model.use_simplify = self.use_simplify
         new_model.convert_to_format = self.convert_to_format

pybamm/models/full_battery_models/base_battery_model.py

@@ -467,6 +467,8 @@ def build_model_equations(self):
                 pybamm.logger.debug(
                     "Setting rhs for {} submodel ({})".format(submodel_name, self.name)
                 )
+                if submodel_name == "external circuit":
+                    n = 1
 
                 submodel.set_rhs(self.variables)
                 pybamm.logger.debug(
@@ -493,7 +495,7 @@ def build_model_equations(self):
                 )
                 self.update(submodel)
 
-    def build_model(self):
+    def build_model(self, build_equations=True):
 
         # Check if already built
         if self._built:
@@ -509,7 +511,10 @@ def build_model(self):
 
         self.build_coupled_variables()
 
-        self.build_model_equations()
+        if build_equations:
+            self.build_model_equations()
+        else:
+            self.update(*self.submodels.values())
 
         pybamm.logger.debug("Setting voltage variables ({})".format(self.name))
         self.set_voltage_variables()
@@ -530,6 +535,30 @@ def build_model(self):
 
         self._built = True
 
+    def new_copy(self, options=None, build=True):
+        """
+        Create a copy of the model. Overwrites the functionality of
+        :class:`pybamm.BaseModel` to make sure that the submodels are updated correctly
+        """
+        options = options or self.options
+        # create without building
+        # 'build' is not a keyword argument for the BaseBatteryModel class, but it
+        # should be for all of the subclasses
+        new_model = self.__class__(name=self.name, options=options, build=False)
+        # update submodels
+        new_model.submodels = self.submodels
+        # now build
+        if build:
+            if self._built is True:
+                new_model.build_model(build_equations=False)
+            else:
+                new_model.build_model(build_equations=True)
+        new_model.use_jacobian = self.use_jacobian
+        new_model.use_simplify = self.use_simplify
+        new_model.convert_to_format = self.convert_to_format
+        new_model.timescale = self.timescale
+        return new_model
+
     def set_external_circuit_submodel(self):
         """
         Define how the external circuit defines the boundary conditions for the model,

pybamm/simulation.py

@@ -124,12 +124,14 @@ def set_up_experiment(self, model, experiment):
         time.
         """
         self.operating_mode = "with experiment"
-        self.model = model.new_copy(
-            options={
-                **model.options,
-                "operating mode": constant_current_constant_voltage_constant_power,
-            }
+        new_model = model.new_copy(build=False)
+        new_model.submodels[
+            "external circuit"
+        ] = pybamm.external_circuit.FunctionControl(
+            new_model.param, constant_current_constant_voltage_constant_power
         )
+        new_model.build_model()
+        self.model = new_model
         if not isinstance(experiment, pybamm.Experiment):
             raise TypeError("experiment must be a pybamm `Experiment` instance")
         # Save the experiment

pybamm/solvers/base_solver.py

@@ -444,8 +444,9 @@ def calculate_consistent_state(self, model, time=0, inputs=None):
             root_sol = self.root_method._integrate(model, [time], inputs)
         except pybamm.SolverError as e:
             raise pybamm.SolverError(
-                "Could not find consistent initial conditions: {}".format(e.args[0])
+                "Could not find consistent states: {}".format(e.args[0])
             )
+        pybamm.logger.info("Found consistent states")
         return root_sol.y.flatten()
 
     def solve(self, model, t_eval=None, external_variables=None, inputs=None):

tests/unit/test_models/test_base_model.py

@@ -229,6 +229,28 @@ def test_update(self):
         self.assertEqual(model.initial_conditions[e], submodel2.initial_conditions[e])
         self.assertEqual(model.boundary_conditions[e], submodel2.boundary_conditions[e])
 
+    def test_new_copy(self):
+        model = pybamm.BaseModel(name="a model")
+        whole_cell = ["negative electrode", "separator", "positive electrode"]
+        c = pybamm.Variable("c", domain=whole_cell)
+        d = pybamm.Variable("d", domain=whole_cell)
+        model.rhs = {c: 5 * pybamm.div(pybamm.grad(d)) - 1, d: -c}
+        model.initial_conditions = {c: 1, d: 2}
+        model.boundary_conditions = {
+            c: {"left": (0, "Dirichlet"), "right": (0, "Dirichlet")},
+            d: {"left": (0, "Dirichlet"), "right": (0, "Dirichlet")},
+        }
+        model.use_jacobian = False
+        model.use_simplify = False
+        model.convert_to_format = "python"
+
+        new_model = model.new_copy()
+        self.assertEqual(new_model.name, model.name)
+        self.assertEqual(new_model.use_jacobian, model.use_jacobian)
+        self.assertEqual(new_model.use_simplify, model.use_simplify)
+        self.assertEqual(new_model.convert_to_format, model.convert_to_format)
+        self.assertEqual(new_model.timescale, model.timescale)
+
     def test_check_well_posedness_variables(self):
         # Well-posed ODE model
         model = pybamm.BaseModel()

tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_spm.py

@@ -104,6 +104,25 @@ def test_surface_form_algebraic(self):
         model = pybamm.lithium_ion.SPM(options)
         model.check_well_posedness()
 
+    def test_new_model(self):
+        model = pybamm.lithium_ion.SPM({"thermal": "x-full"})
+        new_model = model.new_copy()
+        self.assertEqual(new_model.submodels, model.submodels)
+        self.assertEqual(new_model.name, model.name)
+        self.assertEqual(new_model.use_jacobian, model.use_jacobian)
+        self.assertEqual(new_model.use_simplify, model.use_simplify)
+        self.assertEqual(new_model.convert_to_format, model.convert_to_format)
+        self.assertEqual(new_model.timescale, model.timescale)
+
+        # with custom submodels
+        model = pybamm.lithium_ion.SPM({"thermal": "x-full"}, build=False)
+        model.submodels["negative particle"] = pybamm.particle.FastSingleParticle(
+            model.param, "Negative"
+        )
+        model.build_model()
+        new_model = model.new_copy()
+        self.assertEqual(new_model.submodels, model.submodels)
+
 
 class TestSPMExternalCircuits(unittest.TestCase):
     def test_well_posed_voltage(self):