Add MultiFittingProblem class and example (#364)

* Add MultiFittingProblem, example and test * Remove unused n_problems property * Update multi_fitting.py * Update CHANGELOG.md * Remove unused weights * Update problem_list to problem args * Concatenate the whole list * Apply suggestions from code review * Update description * Update default init_soc * Update CHANGELOG.md * Update check_params * Add pybamm_model as default attribute * Ensure predict uses unprocessed_model * Move rebuild check to model.simulate * Align simulate output with predict * Replace init_soc with init_ocv for FittingProblem * Update notebooks * Update test_observers.py * Update descriptions and simplify * Add test_set_initial_state * Copy each model into MultiFittingProblem * Update test_problem.py * Update ecm.py * style: pre-commit fixes * Break connection between parameter_sets * Allow predict to update initial state * Fix typo * Add nbstripout pre-commit hook * Add -q and re-run all notebooks * Copy parameter sets and remove model.initial_state * Reset spm_NelderMead.py * Update CHANGELOG.md * Update CHANGELOG.md * Allow parameter_set is None * Re-run notebooks * Update bounds * Update notebooks * Update notebooks * Set numpy random seed in notebooks * Re-run with fixed seed * Update bounds * Update notebooks to initial_state * Add set_initial_state for ECMs * Add init_ocv setter * Add init_ocv values * Re-run notebooks * Add tests for ECM get_initial_state * Add ECM initial state error tests * Remove unused store_optimised_parameters * Update parameters.initial_value * Use any Initial SoC from parameter_set * Update bounds again * Update init_soc in notebooks * Move dataset check within unscented_kalman * Remove unnecessary lines from spm_UKF * Update all parameters for rebuild * Update init_ocv to _init_ocv * Ensure value updates alongside initial_value * Update multi_model_identification * Update spm_electrode_design.ipynb * Update spm_electrode_design.ipynb * Fix identation * Fix test_plots design problem * Move Changelog entry to breaking changes * Move Changelog entry * style: pre-commit fixes * Fix merge mistake * style: pre-commit fixes * Allow kwargs in MultiFitting evaluate * Add tests * Update integration tests * Update spm_weighted_cost.py * Fix tests * style: pre-commit fixes * Fix model type check * Update _parameter_set to parameter_set * style: pre-commit fixes * Update tests with parameter set * Add model build description * Revert to _parameter_set * Fix predict without pybamm test * Apply suggestions from code review Co-authored-by: Brady Planden <[email protected]> * Apply suggestions from code review Co-authored-by: Brady Planden <[email protected]> * Fix syntax * Fix variable name * Update model type check * Update parameter_set setter * style: pre-commit fixes * Add parameters.reset_initial_value * Add n_outputs property * style: pre-commit fixes * Remove public parameter_set setter * Correct integer to float * Convert initial_state to dict * Add guidance * Remove empty dictionary defaults * style: pre-commit fixes * Add warning stacklevels * Catch simulation errors in problem evaluation * Add pybamm version comment * Add set initial ocv check * Add model.clear and remove setters * Update unscented_kalman.py * Update unscented_kalman.py * Update test_models.py * Update test_set_initial_state * Use clear in model.new_copy * Reference public attributes * Move MultiFittingProblem into separate file * Update description * Add dataset property * Fix changes due to linting * Add test_multi_fitting_problem * Add problem.set_initial_state * Merge rebuild into build * Update CHANGELOG.md * Update base_model.py * Fix notebooks * Update multi_fitting with different initial SoC * Update copying * Add check for identical models * refactor: model.new_copy() args as dictionary and single construction --------- Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com> Co-authored-by: Brady Planden <[email protected]> Co-authored-by: Brady Planden <[email protected]>
pybop-team · Aug 12, 2024 · 778a72f · 778a72f
1 parent df429f7
commit 778a72f
Show file tree

Hide file tree

Showing 14 changed files with 437 additions and 41 deletions.
diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -2,6 +2,7 @@
 
 ## Features
 
+- [#364](https://github.com/pybop-team/PyBOP/pull/364) - Adds the MultiFittingProblem class and the multi_fitting example script.
 - [#444](https://github.com/pybop-team/PyBOP/issues/444) - Merge `BaseModel` `build()` and `rebuild()` functionality.
 - [#435](https://github.com/pybop-team/PyBOP/pull/435) - Adds SLF001 linting for private members.
 - [#418](https://github.com/pybop-team/PyBOP/issues/418) - Wraps the `get_parameter_info` method from PyBaMM to get a dictionary of parameter names and types.
@@ -28,6 +29,7 @@
 
 ## Bug Fixes
 
+
 ## Breaking Changes
 
 # [v24.6](https://github.com/pybop-team/PyBOP/tree/v24.6) - 2024-07-08

diff --git a/examples/scripts/multi_fitting.py b/examples/scripts/multi_fitting.py
@@ -0,0 +1,79 @@
+import numpy as np
+
+import pybop
+
+# Parameter set and model definition
+parameter_set = pybop.ParameterSet.pybamm("Chen2020")
+model = pybop.lithium_ion.SPM(parameter_set=parameter_set)
+
+# Fitting parameters
+parameters = pybop.Parameters(
+    pybop.Parameter(
+        "Negative electrode active material volume fraction",
+        prior=pybop.Gaussian(0.68, 0.05),
+        true_value=parameter_set["Negative electrode active material volume fraction"],
+    ),
+    pybop.Parameter(
+        "Positive electrode active material volume fraction",
+        prior=pybop.Gaussian(0.58, 0.05),
+        true_value=parameter_set["Positive electrode active material volume fraction"],
+    ),
+)
+
+# Generate a dataset and a fitting problem
+sigma = 0.001
+experiment = pybop.Experiment([("Discharge at 0.5C for 2 minutes (4 second period)")])
+values = model.predict(initial_state={"Initial SoC": 0.8}, experiment=experiment)
+dataset_1 = pybop.Dataset(
+    {
+        "Time [s]": values["Time [s]"].data,
+        "Current function [A]": values["Current [A]"].data,
+        "Voltage [V]": values["Voltage [V]"].data
+        + np.random.normal(0, sigma, len(values["Voltage [V]"].data)),
+    }
+)
+problem_1 = pybop.FittingProblem(model, parameters, dataset_1)
+
+# Generate a second dataset and problem
+model = model.new_copy()
+experiment = pybop.Experiment([("Discharge at 1C for 1 minutes (4 second period)")])
+values = model.predict(initial_state={"Initial SoC": 0.8}, experiment=experiment)
+dataset_2 = pybop.Dataset(
+    {
+        "Time [s]": values["Time [s]"].data,
+        "Current function [A]": values["Current [A]"].data,
+        "Voltage [V]": values["Voltage [V]"].data
+        + np.random.normal(0, sigma, len(values["Voltage [V]"].data)),
+    }
+)
+problem_2 = pybop.FittingProblem(model, parameters, dataset_2)
+
+# Combine the problems into one
+problem = pybop.MultiFittingProblem(problem_1, problem_2)
+
+# Generate the cost function and optimisation class
+cost = pybop.SumSquaredError(problem)
+optim = pybop.IRPropMin(
+    cost,
+    verbose=True,
+    max_iterations=100,
+)
+
+# Run optimisation
+x, final_cost = optim.run()
+print("True parameters:", parameters.true_value())
+print("Estimated parameters:", x)
+
+# Plot the timeseries output
+pybop.quick_plot(problem_1, problem_inputs=x, title="Optimised Comparison")
+pybop.quick_plot(problem_2, problem_inputs=x, title="Optimised Comparison")
+
+# Plot convergence
+pybop.plot_convergence(optim)
+
+# Plot the parameter traces
+pybop.plot_parameters(optim)
+
+# Plot the cost landscape with optimisation path
+bounds = np.array([[0.5, 0.8], [0.4, 0.7]])
+pybop.plot2d(optim, bounds=bounds, steps=15)
diff --git a/pybop/__init__.py b/pybop/__init__.py
@@ -87,6 +87,7 @@
 #
 from .problems.base_problem import BaseProblem
 from .problems.fitting_problem import FittingProblem
+from .problems.multi_fitting_problem import MultiFittingProblem
 from .problems.design_problem import DesignProblem
 
 #

diff --git a/pybop/costs/_weighted_cost.py b/pybop/costs/_weighted_cost.py
@@ -16,14 +16,14 @@ class WeightedCost(BaseCost):
 
     Attributes
     ---------------------
-    costs : list[pybop.BaseCost]
-        A list of PyBOP cost objects.
+    costs : pybop.BaseCost
+        The individual PyBOP cost objects.
     weights : list[float]
         A list of values with which to weight the cost values.
-    _has_identical_problems : bool
+    has_identical_problems : bool
         If True, the shared problem will be evaluated once and saved before the
         self.compute() method of each cost is called (default: False).
-    _has_separable_problem: bool
+    has_separable_problem: bool
         This attribute must be set to False for WeightedCost objects. If the
         corresponding attribute of an individual cost is True, the problem is
         separable from the cost function and will be evaluated before the

diff --git a/pybop/models/base_model.py b/pybop/models/base_model.py
@@ -302,7 +302,8 @@ def classify_parameters(
         Parameters
         ----------
         parameters : Parameters, optional
-            The optimisation parameters. Defaults to None, resulting in the internal `pybop.Parameters` object to be used.
+            The optimisation parameters. Defaults to None, resulting in the internal
+            `pybop.Parameters` object to be used.
         inputs : Inputs, optional
             The input parameters for the simulation (default: None).
         """
@@ -337,6 +338,7 @@ def classify_parameters(
             if requires_rebuild:
                 self.clear()
                 self._geometry = self.pybamm_model.default_geometry
+                # Update both the active and unprocessed parameter sets for consistency
                 self._parameter_set.update(rebuild_parameters)
                 self._unprocessed_parameter_set.update(rebuild_parameters)
 
@@ -530,7 +532,7 @@ def predict(
             if PyBaMM models are not supported by the current simulation method.
 
         """
-        if self._unprocessed_model is None:
+        if self.pybamm_model is None:
             raise ValueError(
                 "The predict method currently only supports PyBaMM models."
             )
@@ -644,6 +646,32 @@ def copy(self):
         """
         return copy.copy(self)
 
+    def new_copy(self):
+        """
+        Return a new copy of the model, explicitly copying all the mutable attributes
+        to avoid issues with shared objects.
+
+        Returns
+        -------
+        BaseModel
+            A new copy of the model.
+        """
+        model_class = type(self)
+        if self.pybamm_model is None:
+            model_args = {"parameter_set": self.parameter_set.copy()}
+        else:
+            model_args = {
+                "options": self._unprocessed_model.options,
+                "parameter_set": self._unprocessed_parameter_set.copy(),
+                "geometry": self.pybamm_model.default_geometry.copy(),
+                "submesh_types": self.pybamm_model.default_submesh_types.copy(),
+                "var_pts": self.pybamm_model.default_var_pts.copy(),
+                "spatial_methods": self.pybamm_model.default_spatial_methods.copy(),
+                "solver": self.pybamm_model.default_solver.copy(),
+            }
+
+        return model_class(**model_args)
+
     def get_parameter_info(self, print_info: bool = False):
         """
         Extracts the parameter names and types and returns them as a dictionary.

diff --git a/pybop/optimisers/base_optimiser.py b/pybop/optimisers/base_optimiser.py
@@ -244,7 +244,11 @@ def set_allow_infeasible_solutions(self, allow=True):
         self.physical_viability = allow
         self.allow_infeasible_solutions = allow
 
-        if hasattr(self.cost, "problem") and hasattr(self.cost.problem, "_model"):
+        if (
+            hasattr(self.cost, "problem")
+            and hasattr(self.cost.problem, "model")
+            and self.cost.problem.model is not None
+        ):
             self.cost.problem.model.allow_infeasible_solutions = (
                 self.allow_infeasible_solutions
             )

diff --git a/pybop/problems/base_problem.py b/pybop/problems/base_problem.py
@@ -63,12 +63,23 @@ def __init__(
         self.check_model = check_model
         self.signal = signal or ["Voltage [V]"]
         self.additional_variables = additional_variables or []
-        self.initial_state = initial_state
+        self.set_initial_state(initial_state)
         self._dataset = None
         self._time_data = None
         self._target = None
         self.verbose = False
 
+    def set_initial_state(self, initial_state: Optional[dict] = None):
+        """
+        Set the initial state to be applied to evaluations of the problem.
+
+        Parameters
+        ----------
+        initial_state : dict, optional
+            A valid initial state (default: None).
+        """
+        self.initial_state = initial_state
+
     @property
     def n_parameters(self):
         return len(self.parameters)
@@ -156,3 +167,7 @@ def time_data(self):
     @time_data.setter
     def time_data(self, time_data):
         self._time_data = time_data
+
+    @property
+    def dataset(self):
+        return self._dataset
diff --git a/pybop/problems/design_problem.py b/pybop/problems/design_problem.py
@@ -47,22 +47,6 @@ def __init__(
         additional_variables.extend(["Time [s]", "Current [A]"])
         additional_variables = list(set(additional_variables))
 
-        if initial_state is None:
-            if isinstance(model, ECircuitModel):
-                initial_state = {"Initial SoC": model.parameter_set["Initial SoC"]}
-            else:
-                initial_state = {"Initial SoC": 1.0}  # default value
-        elif "Initial open-circuit voltage [V]" in initial_state.keys():
-            warnings.warn(
-                "It is usually better to define an initial state of charge as the "
-                "initial_state for a DesignProblem because this state will scale with "
-                "design properties such as the capacity of the battery, as opposed to the "
-                "initial open-circuit voltage which may correspond to a different state "
-                "of charge for each design.",
-                UserWarning,
-                stacklevel=1,
-            )
-
         super().__init__(
             parameters, model, check_model, signal, additional_variables, initial_state
         )
@@ -78,6 +62,33 @@ def __init__(
             "Non-physical point encountered",
         ]
 
+    def set_initial_state(self, initial_state):
+        """
+        Set the initial state to be applied to evaluations of the problem.
+
+        Parameters
+        ----------
+        initial_state : dict, optional
+            A valid initial state (default: None).
+        """
+        if initial_state is None:
+            if isinstance(self.model, ECircuitModel):
+                initial_state = {"Initial SoC": self.model.parameter_set["Initial SoC"]}
+            else:
+                initial_state = {"Initial SoC": 1.0}  # default value
+        elif "Initial open-circuit voltage [V]" in initial_state.keys():
+            warnings.warn(
+                "It is usually better to define an initial state of charge as the "
+                "initial_state for a DesignProblem because this state will scale with "
+                "design properties such as the capacity of the battery, as opposed to the "
+                "initial open-circuit voltage which may correspond to a different state "
+                "of charge for each design.",
+                UserWarning,
+                stacklevel=1,
+            )
+
+        self.initial_state = initial_state
+
     def evaluate(self, inputs: Inputs, update_capacity=False):
         """
         Evaluate the model with the given parameters and return the signal.

diff --git a/pybop/problems/fitting_problem.py b/pybop/problems/fitting_problem.py
@@ -21,12 +21,25 @@ class FittingProblem(BaseProblem):
         An object or list of the parameters for the problem.
     dataset : Dataset
         Dataset object containing the data to fit the model to.
+    check_model : bool, optional
+        Flag to indicate if the model should be checked (default: True).
     signal : str, optional
         The variable used for fitting (default: "Voltage [V]").
     additional_variables : list[str], optional
         Additional variables to observe and store in the solution (default additions are: ["Time [s]"]).
     initial_state : dict, optional
         A valid initial state, e.g. the initial open-circuit voltage (default: None).
+
+    Additional Attributes
+    ---------------------
+    dataset : dictionary
+        The dictionary from a Dataset object containing the signal keys and values to fit the model to.
+    time_data : np.ndarray
+        The time points in the dataset.
+    n_time_data : int
+        The number of time points.
+    target : np.ndarray
+        The target values of the signals.
     """
 
     def __init__(
@@ -44,19 +57,6 @@ def __init__(
         additional_variables.extend(["Time [s]"])
         additional_variables = list(set(additional_variables))
 
-        if initial_state is not None and "Initial SoC" in initial_state.keys():
-            warnings.warn(
-                "It is usually better to define an initial open-circuit voltage as the "
-                "initial_state for a FittingProblem because this value can typically be "
-                "obtained from the data, unlike the intrinsic initial state of charge. "
-                "In the case where the fitting parameters do not change the OCV-SOC "
-                "relationship, the initial state of charge may be passed to the model "
-                'using, for example, `model.set_initial_state({"Initial SoC": 1.0})` '
-                "before constructing the FittingProblem.",
-                UserWarning,
-                stacklevel=1,
-            )
-
         super().__init__(
             parameters, model, check_model, signal, additional_variables, initial_state
         )
@@ -82,6 +82,30 @@ def __init__(
                 initial_state=self.initial_state,
             )
 
+    def set_initial_state(self, initial_state: Optional[dict] = None):
+        """
+        Set the initial state to be applied to evaluations of the problem.
+
+        Parameters
+        ----------
+        initial_state : dict, optional
+            A valid initial state (default: None).
+        """
+        if initial_state is not None and "Initial SoC" in initial_state.keys():
+            warnings.warn(
+                "It is usually better to define an initial open-circuit voltage as the "
+                "initial_state for a FittingProblem because this value can typically be "
+                "obtained from the data, unlike the intrinsic initial state of charge. "
+                "In the case where the fitting parameters do not change the OCV-SOC "
+                "relationship, the initial state of charge may be passed to the model "
+                'using, for example, `model.set_initial_state({"Initial SoC": 1.0})` '
+                "before constructing the FittingProblem.",
+                UserWarning,
+                stacklevel=1,
+            )
+
+        self.initial_state = initial_state
+
     def evaluate(
         self, inputs: Inputs, update_capacity=False
     ) -> dict[str, np.ndarray[np.float64]]:
@@ -130,6 +154,7 @@ def evaluateS1(self, inputs: Inputs):
             dy/dx(t) evaluated with given inputs.
         """
         inputs = self.parameters.verify(inputs)
+        self.parameters.update(values=list(inputs.values()))
 
         try:
             sol = self._model.simulateS1(