Added Model Entry Points (#23)

+ Implemented `Model Entry Points` to create Model objects within the
template and initialise them through entry points.
+ Added the `SPM` model which can be initialised through the model entry
points.
+ A wrapper method to load a model object called
`models("modelname/authorname")` is added.
Example - 
To load the `SPM` model, after installing the `pybamm_cookiecutter`
project, it can be accessed by calling,
`pybamm_cookiecutter.Model("SPM")`. This would return an initialised
model object of the `SPM` model.
+ Added two basic tests for model entry points.

---------

Co-authored-by: Agriya Khetarpal <[email protected]>
Co-authored-by: Ferran Brosa Planella <[email protected]>

pyproject.toml

@@ -61,9 +61,12 @@ Homepage = "https://github.com/pybamm-team/pybamm-cookiecutter"
 Discussions = "https://github.com/pybamm-team/pybamm-cookiecutter/discussions"
 Changelog = "https://github.com/pybamm-team/pybamm-cookiecutter/releases"
 
-[project.entry-points."cookie_parameter_sets"]
+[project.entry-points."parameter_sets"]
 Chen2020 = "pybamm_cookiecutter.parameters.input.Chen2020:get_parameter_values"
 
+[project.entry-points."models"]
+SPM = "pybamm_cookiecutter.models.input.SPM:SPM"
+
 [tool.hatch]
 version.source = "vcs"
 build.hooks.vcs.version-file = "src/pybamm_cookiecutter/_version.py"

src/pybamm_cookiecutter/__init__.py

@@ -8,10 +8,12 @@
 import pybamm
 
 from ._version import version as __version__
-from .parameters.parameter_sets import parameter_sets
+from .entry_point import Model, parameter_sets, models
 
 __all__ : list[str] = [
     "__version__",
     "pybamm",
     "parameter_sets",
+    "Model",
+    "models",
 ]

src/pybamm_cookiecutter/parameters/parameter_sets.py → src/pybamm_cookiecutter/entry_point.py

@@ -37,19 +37,22 @@
 from collections.abc import Mapping
 from typing import Callable
 
-class ParameterSets(Mapping):
+class EntryPoint(Mapping):
     """
-    Dict-like interface for accessing parameter sets through entry points in cookiecutter template.
-    Access via :py:data:`pybamm_cookiecutter.parameter_sets`
+    Dict-like interface for accessing parameter sets and models through entry points in cookiecutter template.
+    Access via :py:data:`pybamm_cookiecutter.parameter_sets` for parameter_sets
+    Access via :py:data:`pybamm_cookiecutter.Model` for Models
 
     Examples
     --------
     Listing available parameter sets:
         >>> import pybamm_cookiecutter
         >>> list(pybamm_cookiecutter.parameter_sets)
         ['Chen2020', ...]
+        >>> list(pybamm_cookiecutter.models)
+        ['SPM', ...]
 
-    Get the docstring for a parameter set:
+    Get the docstring for a parameter set/model:
 
 
         >>> print(pybamm_cookiecutter.parameter_sets.get_docstring("Ai2020"))
@@ -58,15 +61,23 @@ class ParameterSets(Mapping):
         :footcite:t:`rieger2016new` and references therein.
         ...
 
-    See also: :ref:`adding-parameter-sets`
-
+        >>> print(pybamm_cookiecutter.models.get_docstring("SPM"))
+        <BLANKLINE>
+        Single Particle Model (SPM) model of a lithium-ion battery, from :footcite:t:`Marquis2019`. This class differs from the :class:`pybamm.lithium_ion.SPM` model class in that it shows the whole model in a single class. This comes at the cost of flexibility in combining different physical effects, and in general the main SPM class should be used instead.
+        ...
+        See also: :ref:`adding-parameter-sets`
     """
 
-    def __init__(self):
-        """Dict of entry points for parameter sets, lazily load entry points as"""
-        self.__all_parameter_sets = dict()
-        for entry_point in self.get_entries("cookie_parameter_sets"):
-            self.__all_parameter_sets[entry_point.name] = entry_point
+    _instances = 0
+    def __init__(self, group):
+        """Dict of entry points for parameter sets or models, lazily load entry points as"""
+        if not hasattr(self, 'initialized'):    # Ensure __init__ is called once per instance
+            self.initialized = True
+            EntryPoint._instances += 1
+            self._all_entries = dict()
+            self.group = group
+            for entry_point in self.get_entries(self.group):
+                self._all_entries[entry_point.name] = entry_point
 
     @staticmethod
     def get_entries(group_name):
@@ -76,35 +87,35 @@ def get_entries(group_name):
         else:
             return importlib.metadata.entry_points(group=group_name)
 
-    def __new__(cls):
-        """Ensure only one instance of ParameterSets exists"""
-        if not hasattr(cls, "instance"):
+    def __new__(cls, group):
+        """Ensure only two instances of entry points exist, one for parameter sets and the other for models"""
+        if EntryPoint._instances < 2:
             cls.instance = super().__new__(cls)
         return cls.instance
 
     def __getitem__(self, key) -> dict:
         return self._load_entry_point(key)()
 
     def _load_entry_point(self, key) -> Callable:
-        """Check that ``key`` is a registered ``cookie_parameter_sets``,
-        and return the entry point for the parameter set, loading it needed."""
-        if key not in self.__all_parameter_sets:
-            raise KeyError(f"Unknown parameter set: {key}")
-        ps = self.__all_parameter_sets[key]
+        """Check that ``key`` is a registered ``parameter_sets`` or ``models` ,
+        and return the entry point for the parameter set/model, loading it needed."""
+        if key not in self._all_entries:
+            raise KeyError(f"Unknown parameter set or model: {key}")
+        ps = self._all_entries[key]
         try:
-            ps = self.__all_parameter_sets[key] = ps.load()
+            ps = self._all_entries[key] = ps.load()
         except AttributeError:
             pass
         return ps
 
     def __iter__(self):
-        return self.__all_parameter_sets.__iter__()
+        return self._all_entries.__iter__()
 
     def __len__(self) -> int:
-        return len(self.__all_parameter_sets)
+        return len(self._all_entries)
 
     def get_docstring(self, key):
-        """Return the docstring for the ``key`` parameter set"""
+        """Return the docstring for the ``key`` parameter set or model"""
         return textwrap.dedent(self._load_entry_point(key).__doc__)
 
     def __getattribute__(self, name):
@@ -114,4 +125,30 @@ def __getattribute__(self, name):
               raise error
 
 #: Singleton Instance of :class:ParameterSets """
-parameter_sets = ParameterSets()
+parameter_sets = EntryPoint(group="parameter_sets")
+
+#: Singleton Instance of :class:ModelEntryPoints"""
+models = EntryPoint(group="models")
+
+def Model(model:str):
+    """
+    Returns the loaded model object
+
+    Parameters
+    ----------
+    model : str
+        The model name or author name of the model mentioned at the model entry point.
+    Returns
+    -------
+    pybamm.model
+        Model object of the initialised model.
+    Examples
+    --------
+    Listing available models:
+        >>> import pybamm_cookiecutter
+        >>> list(pybamm_cookiecutter.models)
+        ['SPM', ...]
+        >>> pybamm_cookiecutter.Model('Author/Year')
+        <pybamm_cookiecutter.models.input.SPM.SPM object>
+    """
+    return models[model]

src/pybamm_cookiecutter/models/input/SPM.py

@@ -0,0 +1,214 @@
+"""
+This code is adopted from the PyBaMM project under the BSD-3-Clause
+
+Copyright (c) 2018-2024, the PyBaMM team.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are met:
+
+* Redistributions of source code must retain the above copyright notice, this
+  list of conditions and the following disclaimer.
+
+* Redistributions in binary form must reproduce the above copyright notice,
+  this list of conditions and the following disclaimer in the documentation
+  and/or other materials provided with the distribution.
+
+* Neither the name of the copyright holder nor the names of its
+  contributors may be used to endorse or promote products derived from
+  this software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
+FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+"""
+
+
+#
+# Basic Single Particle Model (SPM)
+#
+import pybamm
+
+class SPM(pybamm.lithium_ion.BaseModel):
+    """Single Particle Model (SPM) model of a lithium-ion battery, from
+    :footcite:t:`Marquis2019`.
+
+    This class differs from the :class:`pybamm.lithium_ion.SPM` model class in that it
+    shows the whole model in a single class. This comes at the cost of flexibility in
+    combining different physical effects, and in general the main SPM class should be
+    used instead.
+
+    Parameters
+    ----------
+    name : str, optional
+        The name of the model.
+    """
+
+    def __init__(self, name="Single Particle Model"):
+        super().__init__({}, name)
+        pybamm.citations.register("Marquis2019")
+        # `param` is a class containing all the relevant parameters and functions for
+        # this model. These are purely symbolic at this stage, and will be set by the
+        # `ParameterValues` class when the model is processed.
+        param = self.param
+
+        ######################
+        # Variables
+        ######################
+        # Variables that depend on time only are created without a domain
+        Q = pybamm.Variable("Discharge capacity [A.h]")
+        # Variables that vary spatially are created with a domain
+        c_s_n = pybamm.Variable(
+            "X-averaged negative particle concentration [mol.m-3]",
+            domain="negative particle",
+        )
+        c_s_p = pybamm.Variable(
+            "X-averaged positive particle concentration [mol.m-3]",
+            domain="positive particle",
+        )
+
+        # Constant temperature
+        T = param.T_init
+
+        ######################
+        # Other set-up
+        ######################
+
+        # Current density
+        i_cell = param.current_density_with_time
+        a_n = 3 * param.n.prim.epsilon_s_av / param.n.prim.R_typ
+        a_p = 3 * param.p.prim.epsilon_s_av / param.p.prim.R_typ
+        j_n = i_cell / (param.n.L * a_n)
+        j_p = -i_cell / (param.p.L * a_p)
+
+        ######################
+        # State of Charge
+        ######################
+        I = param.current_with_time
+        # The `rhs` dictionary contains differential equations, with the key being the
+        # variable in the d/dt
+        self.rhs[Q] = I / 3600
+        # Initial conditions must be provided for the ODEs
+        self.initial_conditions[Q] = pybamm.Scalar(0)
+
+        ######################
+        # Particles
+        ######################
+
+        # The div and grad operators will be converted to the appropriate matrix
+        # multiplication at the discretisation stage
+        N_s_n = -param.n.prim.D(c_s_n, T) * pybamm.grad(c_s_n)
+        N_s_p = -param.p.prim.D(c_s_p, T) * pybamm.grad(c_s_p)
+        self.rhs[c_s_n] = -pybamm.div(N_s_n)
+        self.rhs[c_s_p] = -pybamm.div(N_s_p)
+        # Surf takes the surface value of a variable, i.e. its boundary value on the
+        # right side. This is also accessible via `boundary_value(x, "right")`, with
+        # "left" providing the boundary value of the left side
+        c_s_surf_n = pybamm.surf(c_s_n)
+        c_s_surf_p = pybamm.surf(c_s_p)
+        # Boundary conditions must be provided for equations with spatial derivatives
+        self.boundary_conditions[c_s_n] = {
+            "left": (pybamm.Scalar(0), "Neumann"),
+            "right": (
+                -j_n / (param.F * pybamm.surf(param.n.prim.D(c_s_n, T))),
+                "Neumann",
+            ),
+        }
+        self.boundary_conditions[c_s_p] = {
+            "left": (pybamm.Scalar(0), "Neumann"),
+            "right": (
+                -j_p / (param.F * pybamm.surf(param.p.prim.D(c_s_p, T))),
+                "Neumann",
+            ),
+        }
+        # c_n_init and c_p_init are functions of r and x, but for the SPM we
+        # take the x-averaged value since there is no x-dependence in the particles
+        self.initial_conditions[c_s_n] = pybamm.x_average(param.n.prim.c_init)
+        self.initial_conditions[c_s_p] = pybamm.x_average(param.p.prim.c_init)
+        # Events specify points at which a solution should terminate
+        sto_surf_n = c_s_surf_n / param.n.prim.c_max
+        sto_surf_p = c_s_surf_p / param.p.prim.c_max
+        self.events += [
+            pybamm.Event(
+                "Minimum negative particle surface stoichiometry",
+                pybamm.min(sto_surf_n) - 0.01,
+            ),
+            pybamm.Event(
+                "Maximum negative particle surface stoichiometry",
+                (1 - 0.01) - pybamm.max(sto_surf_n),
+            ),
+            pybamm.Event(
+                "Minimum positive particle surface stoichiometry",
+                pybamm.min(sto_surf_p) - 0.01,
+            ),
+            pybamm.Event(
+                "Maximum positive particle surface stoichiometry",
+                (1 - 0.01) - pybamm.max(sto_surf_p),
+            ),
+        ]
+
+        # Note that the SPM does not have any algebraic equations, so the `algebraic`
+        # dictionary remains empty
+
+        ######################
+        # (Some) variables
+        ######################
+        # Interfacial reactions
+        RT_F = param.R * T / param.F
+        j0_n = param.n.prim.j0(param.c_e_init_av, c_s_surf_n, T)
+        j0_p = param.p.prim.j0(param.c_e_init_av, c_s_surf_p, T)
+        eta_n = (2 / param.n.prim.ne) * RT_F * pybamm.arcsinh(j_n / (2 * j0_n))
+        eta_p = (2 / param.p.prim.ne) * RT_F * pybamm.arcsinh(j_p / (2 * j0_p))
+        phi_s_n = 0
+        phi_e = -eta_n - param.n.prim.U(sto_surf_n, T)
+        phi_s_p = eta_p + phi_e + param.p.prim.U(sto_surf_p, T)
+        V = phi_s_p
+        num_cells = pybamm.Parameter(
+            "Number of cells connected in series to make a battery"
+        )
+
+        whole_cell = ["negative electrode", "separator", "positive electrode"]
+        # The `variables` dictionary contains all variables that might be useful for
+        # visualising the solution of the model
+        # Primary broadcasts are used to broadcast scalar quantities across a domain
+        # into a vector of the right shape, for multiplying with other vectors
+        self.variables = {
+            "Time [s]": pybamm.t,
+            "Discharge capacity [A.h]": Q,
+            "X-averaged negative particle concentration [mol.m-3]": c_s_n,
+            "Negative particle surface "
+            "concentration [mol.m-3]": pybamm.PrimaryBroadcast(
+                c_s_surf_n, "negative electrode"
+            ),
+            "Electrolyte concentration [mol.m-3]": pybamm.PrimaryBroadcast(
+                param.c_e_init_av, whole_cell
+            ),
+            "X-averaged positive particle concentration [mol.m-3]": c_s_p,
+            "Positive particle surface "
+            "concentration [mol.m-3]": pybamm.PrimaryBroadcast(
+                c_s_surf_p, "positive electrode"
+            ),
+            "Current [A]": I,
+            "Current variable [A]": I,  # for compatibility with pybamm.Experiment
+            "Negative electrode potential [V]": pybamm.PrimaryBroadcast(
+                phi_s_n, "negative electrode"
+            ),
+            "Electrolyte potential [V]": pybamm.PrimaryBroadcast(phi_e, whole_cell),
+            "Positive electrode potential [V]": pybamm.PrimaryBroadcast(
+                phi_s_p, "positive electrode"
+            ),
+            "Voltage [V]": V,
+            "Battery voltage [V]": V * num_cells,
+        }
+        # Events specify points at which a solution should terminate
+        self.events += [
+            pybamm.Event("Minimum voltage [V]", V - param.voltage_low_cut),
+            pybamm.Event("Maximum voltage [V]", param.voltage_high_cut - V),
+        ]