add ocp as a separate submodel

CHANGELOG.md

@@ -1,8 +1,8 @@
 # [Unreleased](https://github.com/pybamm-team/PyBaMM/)
 ## Features
 
+-   Added open-circuit potential as a separate submodel ([#2094](https://github.com/pybamm-team/PyBaMM/pull/2094))
 -   Added partially reversible lithium plating model and new `OKane2022` parameter set to go with it ([#2043](https://github.com/pybamm-team/PyBaMM/pull/2043))
-
 -   Added `__eq__` and `__hash__` methods for `Symbol` objects, using `.id` ([#1978](https://github.com/pybamm-team/PyBaMM/pull/1978))
 
 ## Optimizations

docs/requirements.txt

@@ -1,5 +1,5 @@
 # Requirements for readthedocs.io
-numpy >= 1.16
+numpy <= 1.22 # change back to numpy>=1.16 once scikit.odes is fixed
 scipy >= 1.3
 pandas >= 0.24
 anytree >= 2.4.3

examples/notebooks/models/using-submodels.ipynb

@@ -417,7 +417,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "In the Single Particle Model, the overpotential can be obtianed by inverting the Butler-Volmer relation, so we choose the `InverseButlerVolmer` submodel for the interface, with the \"main\" lithium-ion reaction (and default lithium ion options). Because of how the current is implemented, we also need to separately specify the `CurrentForInverseButlerVolmer` submodel"
+    "In the Single Particle Model, the overpotential can be obtained by inverting the Butler-Volmer relation, so we choose the `InverseButlerVolmer` submodel for the interface, with the \"main\" lithium-ion reaction (and default lithium ion options). Because of how the current is implemented, we also need to separately specify the `CurrentForInverseButlerVolmer` submodel. We also need to specify the submodel for open-circuit potential."
    ]
   },
   {
@@ -426,6 +426,16 @@
    "metadata": {},
    "outputs": [],
    "source": [
+    "model.submodels[\n",
+    "    \"negative open circuit potential\"\n",
+    "] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(\n",
+    "    model.param, \"Negative\", \"lithium-ion main\", options=model.options\n",
+    ")\n",
+    "model.submodels[\n",
+    "    \"positive open circuit potential\"\n",
+    "] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(\n",
+    "    model.param, \"Positive\", \"lithium-ion main\", options=model.options\n",
+    ")\n",
     "model.submodels[\n",
     "    \"negative interface\"\n",
     "] = pybamm.kinetics.InverseButlerVolmer(\n",

examples/scripts/custom_model.py

@@ -37,6 +37,17 @@
     model.param, "Positive", "uniform profile", options=model.options
 )
 model.submodels["positive particle"] = particle_p
+
+model.submodels[
+    "negative open circuit potential"
+] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+    model.param, "Negative", "lithium-ion main", options=model.options
+)
+model.submodels[
+    "positive open circuit potential"
+] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+    model.param, "Positive", "lithium-ion main", options=model.options
+)
 model.submodels["negative interface"] = pybamm.kinetics.InverseButlerVolmer(
     model.param, "Negative", "lithium-ion main", options=model.options
 )

pybamm/__init__.py

@@ -136,6 +136,7 @@
 from .models.submodels.interface import sei
 from .models.submodels.interface import lithium_plating
 from .models.submodels.interface import interface_utilisation
+from .models.submodels.interface import open_circuit_potential
 
 #
 # Geometry

pybamm/models/full_battery_models/lead_acid/base_lead_acid_model.py

@@ -81,6 +81,19 @@ def set_soc_variables(self):
             self.rhs[fci] = -self.variables["Total current density"] * 100
             self.initial_conditions[fci] = self.param.q_init * 100
 
+    def set_open_circuit_potential_submodel(self):
+        for domain in ["Negative", "Positive"]:
+            self.submodels[
+                f"{domain.lower()} open circuit potential"
+            ] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+                self.param, domain, "lead-acid main", self.options
+            )
+            self.submodels[
+                f"{domain.lower()} oxygen open circuit potential"
+            ] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+                self.param, domain, "lead-acid oxygen", self.options
+            )
+
     def set_active_material_submodel(self):
         self.submodels["negative active material"] = pybamm.active_material.Constant(
             self.param, "Negative", self.options

pybamm/models/full_battery_models/lead_acid/full.py

@@ -35,6 +35,7 @@ def __init__(self, options=None, name="Full model", build=True):
         super().__init__(options, name)
 
         self.set_external_circuit_submodel()
+        self.set_open_circuit_potential_submodel()
         self.set_intercalation_kinetics_submodel()
         self.set_interface_utilisation_submodel()
         self.set_porosity_submodel()

pybamm/models/full_battery_models/lead_acid/higher_order.py

@@ -41,6 +41,7 @@ def __init__(self, options=None, name="Composite model", build=True):
         self.set_electrolyte_diffusion_submodel()
         self.set_other_species_diffusion_submodels()
         # Average interface submodel to get average first-order potential differences
+        self.set_open_circuit_potential_submodel()
         self.set_average_interfacial_submodel()
         # Electrolyte and solid submodels to get full first-order potentials
         self.set_negative_electrode_submodel()

pybamm/models/full_battery_models/lead_acid/loqs.py

@@ -34,6 +34,7 @@ def __init__(self, options=None, name="LOQS model", build=True):
         super().__init__(options, name)
 
         self.set_external_circuit_submodel()
+        self.set_open_circuit_potential_submodel()
         self.set_intercalation_kinetics_submodel()
         self.set_interface_utilisation_submodel()
         self.set_convection_submodel()

pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py

@@ -44,6 +44,33 @@ def __init__(self, options=None, name="Unnamed lithium-ion model", build=False):
             )
         self.set_standard_output_variables()
 
+    def set_submodels(self, build):
+        self.set_external_circuit_submodel()
+        self.set_porosity_submodel()
+        self.set_interface_utilisation_submodel()
+        self.set_crack_submodel()
+        self.set_active_material_submodel()
+        self.set_transport_efficiency_submodels()
+        self.set_convection_submodel()
+        self.set_open_circuit_potential_submodel()
+        self.set_intercalation_kinetics_submodel()
+        self.set_other_reaction_submodels_to_zero()
+        self.set_particle_submodel()
+        self.set_solid_submodel()
+        self.set_electrolyte_submodel()
+        self.set_thermal_submodel()
+        self.set_current_collector_submodel()
+
+        self.set_sei_submodel()
+        self.set_lithium_plating_submodel()
+
+        if self.half_cell:
+            # This also removes "negative electrode" submodels, so should be done last
+            self.set_li_metal_counter_electrode_submodels()
+
+        if build:
+            self.build_model()
+
     @property
     def default_parameter_values(self):
         if self.half_cell:
@@ -189,6 +216,14 @@ def set_summary_variables(self):
 
         self.summary_variables = summary_variables
 
+    def set_open_circuit_potential_submodel(self):
+        for domain in ["Negative", "Positive"]:
+            self.submodels[
+                f"{domain.lower()} open circuit potential"
+            ] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+                self.param, domain, "lithium-ion main", self.options
+            )
+
     def set_sei_submodel(self):
         if self.half_cell:
             reaction_loc = "interface"
@@ -217,12 +252,15 @@ def set_lithium_plating_submodel(self):
             )
 
     def set_other_reaction_submodels_to_zero(self):
-        self.submodels["negative oxygen interface"] = pybamm.kinetics.NoReaction(
-            self.param, "Negative", "lithium-ion oxygen"
-        )
-        self.submodels["positive oxygen interface"] = pybamm.kinetics.NoReaction(
-            self.param, "Positive", "lithium-ion oxygen"
-        )
+        for domain in ["Negative", "Positive"]:
+            self.submodels[
+                f"{domain.lower()} oxygen interface"
+            ] = pybamm.kinetics.NoReaction(self.param, domain, "lithium-ion oxygen")
+            self.submodels[
+                f"{domain.lower()} oxygen open circuit potential"
+            ] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+                self.param, domain, "lithium-ion oxygen", self.options
+            )
 
     def set_crack_submodel(self):
         for domain in ["Negative", "Positive"]:
@@ -271,6 +309,12 @@ def set_porosity_submodel(self):
             )
 
     def set_li_metal_counter_electrode_submodels(self):
+        self.submodels[
+            "counter electrode open circuit potential"
+        ] = pybamm.open_circuit_potential.SingleOpenCircuitPotential(
+            self.param, "Negative", "lithium metal plating", self.options
+        )
+
         if (
             self.options["SEI"] in ["none", "constant"]
             and self.options["intercalation kinetics"] == "symmetric Butler-Volmer"

pybamm/models/full_battery_models/lithium_ion/dfn.py

@@ -37,33 +37,11 @@ class DFN(BaseModel):
     """
 
     def __init__(self, options=None, name="Doyle-Fuller-Newman model", build=True):
-        super().__init__(options, name)
         # For degradation models we use the full form since this is a full-order model
         self.x_average = False
+        super().__init__(options, name)
 
-        self.set_external_circuit_submodel()
-        self.set_porosity_submodel()
-        self.set_interface_utilisation_submodel()
-        self.set_crack_submodel()
-        self.set_active_material_submodel()
-        self.set_transport_efficiency_submodels()
-        self.set_convection_submodel()
-        self.set_intercalation_kinetics_submodel()
-        self.set_other_reaction_submodels_to_zero()
-        self.set_particle_submodel()
-        self.set_solid_submodel()
-        self.set_electrolyte_submodel()
-        self.set_thermal_submodel()
-        self.set_current_collector_submodel()
-        self.set_sei_submodel()
-        self.set_lithium_plating_submodel()
-
-        if self.half_cell:
-            # This also removes "negative electrode" submodels, so should be done last
-            self.set_li_metal_counter_electrode_submodels()
-
-        if build:
-            self.build_model()
+        self.set_submodels(build)
 
         pybamm.citations.register("Doyle1993")
 

pybamm/models/full_battery_models/lithium_ion/spm.py

@@ -42,35 +42,14 @@ def __init__(self, options=None, name="Single Particle Model", build=True):
         surface_form = options.get("surface form")
         if kinetics is not None and surface_form is None:
             options["surface form"] = "algebraic"
-        super().__init__(options, name)
+
         # For degradation models we use the "x-average" form since this is a
         # reduced-order model with uniform current density in the electrodes
         self.x_average = True
 
-        self.set_external_circuit_submodel()
-        self.set_porosity_submodel()
-        self.set_interface_utilisation_submodel()
-        self.set_crack_submodel()
-        self.set_active_material_submodel()
-        self.set_transport_efficiency_submodels()
-        self.set_convection_submodel()
-        self.set_intercalation_kinetics_submodel()
-        self.set_other_reaction_submodels_to_zero()
-        self.set_particle_submodel()
-        self.set_solid_submodel()
-        self.set_electrolyte_submodel()
-        self.set_thermal_submodel()
-        self.set_current_collector_submodel()
-
-        self.set_sei_submodel()
-        self.set_lithium_plating_submodel()
-
-        if self.half_cell:
-            # This also removes "negative electrode" submodels, so should be done last
-            self.set_li_metal_counter_electrode_submodels()
-
-        if build:
-            self.build_model()
+        super().__init__(options, name)
+
+        self.set_submodels(build)
 
         if self.__class__ != "MPM":
             pybamm.citations.register("Marquis2019")

pybamm/models/full_battery_models/lithium_ion/spme.py

@@ -39,36 +39,12 @@ class SPMe(SPM):
     def __init__(
         self, options=None, name="Single Particle Model with electrolyte", build=True
     ):
-        super().__init__(options, name, build=False)
         # For degradation models we use the "x-average" form since this is a
         # reduced-order model with uniform current density in the electrodes
         self.x_average = True
 
-        self.set_external_circuit_submodel()
-        self.set_porosity_submodel()
-        self.set_interface_utilisation_submodel()
-        self.set_crack_submodel()
-        self.set_active_material_submodel()
-        self.set_transport_efficiency_submodels()
-        self.set_convection_submodel()
-        self.set_intercalation_kinetics_submodel()
-        self.set_other_reaction_submodels_to_zero()
-        self.set_particle_submodel()
-        self.set_solid_submodel()
-        self.set_electrolyte_submodel()
-        self.set_thermal_submodel()
-        self.set_current_collector_submodel()
-        self.set_sei_submodel()
-        self.set_lithium_plating_submodel()
-
-        if self.half_cell:
-            # This also removes "negative electrode" submodels, so should be done last
-            self.set_li_metal_counter_electrode_submodels()
-
-        if build:
-            self.build_model()
-
-        pybamm.citations.register("Marquis2019")
+        # Initialize with the SPM
+        super().__init__(options, name)
 
     def set_convection_submodel(self):