Tip

An interactive online version of this notebook is available, which can be accessed via

Alternatively, you may download this notebook and run it offline.

Attention

You are viewing this notebook on the latest version of the documentation, where these notebooks may not be compatible with the stable release of PyBaMM since they can contain features that are not yet released. We recommend viewing these notebooks from the stable version of the documentation. To install the latest version of PyBaMM that is compatible with the latest notebooks, build PyBaMM from source.

Using submodels in PyBaMM#

In this notebook we show how to modify existing models by swapping out submodels, and how to build your own model from scratch using existing submodels. To see all of the models and submodels available in PyBaMM, please take a look at the documentation here.

Changing a submodel in an existing battery model#

PyBaMM is designed to be a flexible modelling package that allows users to easily compare different models and numerical techniques within a common framework. Battery models within PyBaMM are built up using a number of submodels that describe different physics included within the model, such as mass conservation in the electrolyte or charge conservation in the solid. For ease of use, a number of popular battery models are pre-configured in PyBaMM. As an example, we look at the Single Particle Model (for more information see here).

First we import pybamm

[1]:

%pip install "pybamm[plot,cite]" -q    # install PyBaMM if it is not installed
import pybamm

Note: you may need to restart the kernel to use updated packages.

Then we load the SPM

[2]:

model = pybamm.lithium_ion.SPM()

We can look at the submodels that make up the SPM by accessing model.submodels, which is a dictionary of the submodel’s name (i.e. the physics it represents) and the submodel that is selected

[3]:

for name, submodel in model.submodels.items():
    print(name, submodel)

external circuit <pybamm.models.submodels.external_circuit.explicit_control_external_circuit.ExplicitCurrentControl object at 0x7f9cd479ea00>
porosity <pybamm.models.submodels.porosity.constant_porosity.Constant object at 0x7f9cd479ecd0>
Negative interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x7f9cd479ee80>
Positive interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x7f9cd479eee0>
negative particle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x7f9cd479ef40>
positive particle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x7f9cd479efd0>
negative primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x7f9cd47a60a0>
positive primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x7f9cd47a6160>
electrolyte transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman_transport_efficiency.Bruggeman object at 0x7f9cd479ed30>
electrode transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman_transport_efficiency.Bruggeman object at 0x7f9cd47a61f0>
transverse convection <pybamm.models.submodels.convection.transverse.no_convection.NoConvection object at 0x7f9cd47a6220>
through-cell convection <pybamm.models.submodels.convection.through_cell.no_convection.NoConvection object at 0x7f9cd47a62b0>
negative primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x7f9cd479ec40>
positive primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x7f9cd479ed60>
negative interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x7f9cd479ea90>
negative interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.CurrentForInverseButlerVolmer object at 0x7f9cd479e880>
positive interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x7f9cd479e610>
positive interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.CurrentForInverseButlerVolmer object at 0x7f9cd479e6d0>
negative primary particle <pybamm.models.submodels.particle.fickian_diffusion.FickianDiffusion object at 0x7f9cd479e3a0>
negative primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x7f9cd479e490>
positive primary particle <pybamm.models.submodels.particle.fickian_diffusion.FickianDiffusion object at 0x7f9cd479e460>
positive primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x7f9cd479e2e0>
negative electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x7f9cd479e310>
positive electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x7f9cd49119a0>
electrolyte diffusion <pybamm.models.submodels.electrolyte_diffusion.constant_concentration.ConstantConcentration object at 0x7f9cd4816d30>
leading-order electrolyte conductivity <pybamm.models.submodels.electrolyte_conductivity.leading_order_conductivity.LeadingOrder object at 0x7f9cd4911fd0>
negative surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x7f9cd4911e50>
positive surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x7f9cd4911cd0>
thermal <pybamm.models.submodels.thermal.isothermal.Isothermal object at 0x7f9cd4911bb0>
current collector <pybamm.models.submodels.current_collector.homogeneous_current_collector.Uniform object at 0x7f9cd4911ac0>
negative primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd4907220>
positive primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd49078b0>
negative primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd4819880>
positive primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd4819610>
negative lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x7f9cd48e41c0>
positive lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x7f9cd48f66a0>
total interface <pybamm.models.submodels.interface.total_interfacial_current.TotalInterfacialCurrent object at 0x7f9cd4911700>

When you load a model in PyBaMM it builds by default. Building the model sets all of the model variables and sets up any variables which are coupled between different submodels: this is the process which couples the submodels together and allows one submodel to access variables from another. If you would like to swap out a submodel in an existing battery model you need to load it without building it by passing the keyword build=False

[4]:

model = pybamm.lithium_ion.SPM(build=False)

This collects all of the submodels which make up the SPM, but doesn’t build the model. Now you are free to swap out one submodel for another. For instance, you may want to assume that diffusion within the negative particles is infinitely fast, so that the PDE describing diffusion is replaced with an ODE for the uniform particle concentration. To change a submodel you simply update the dictionary entry, in this case to the XAveragedPolynomialProfile submodel

[5]:

model.submodels["negative primary particle"] = (
    pybamm.particle.XAveragedPolynomialProfile(
        model.param,
        "negative",
        options={**model.options, "particle": "uniform profile"},
    )
)

where we pass in the model parameters, the electrode (negative or positive) the submodel corresponds to, and the name of the polynomial we want to use. In the example we assume uniform concentration within the particle, corresponding to a zero-order polynomial.

Now if we look at the submodels again we see that the model for the negative particle has been changed

[6]:

for name, submodel in model.submodels.items():
    print(name, submodel)

external circuit <pybamm.models.submodels.external_circuit.explicit_control_external_circuit.ExplicitCurrentControl object at 0x7f9cd48d94c0>
porosity <pybamm.models.submodels.porosity.constant_porosity.Constant object at 0x7f9cd48cde80>
Negative interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x7f9cd48cd970>
Positive interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x7f9cd48cd850>
negative particle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x7f9cd48cd790>
positive particle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x7f9cd48cd520>
negative primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x7f9cd48cd7c0>
positive primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x7f9cd48cd460>
electrolyte transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman_transport_efficiency.Bruggeman object at 0x7f9cd48cdb50>
electrode transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman_transport_efficiency.Bruggeman object at 0x7f9cd491dc70>
transverse convection <pybamm.models.submodels.convection.transverse.no_convection.NoConvection object at 0x7f9cd48b0ac0>
through-cell convection <pybamm.models.submodels.convection.through_cell.no_convection.NoConvection object at 0x7f9cd48c48e0>
negative primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x7f9cd481fd30>
positive primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x7f9cd481f8e0>
negative interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x7f9cd481fd00>
negative interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.CurrentForInverseButlerVolmer object at 0x7f9cd481f9d0>
positive interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x7f9cd481fa60>
positive interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.CurrentForInverseButlerVolmer object at 0x7f9cd481f850>
negative primary particle <pybamm.models.submodels.particle.x_averaged_polynomial_profile.XAveragedPolynomialProfile object at 0x7f9cd4952160>
negative primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x7f9cd481f6a0>
positive primary particle <pybamm.models.submodels.particle.fickian_diffusion.FickianDiffusion object at 0x7f9cd481f7c0>
positive primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x7f9cd481f4c0>
negative electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x7f9cd481f550>
positive electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x7f9cd481f1c0>
electrolyte diffusion <pybamm.models.submodels.electrolyte_diffusion.constant_concentration.ConstantConcentration object at 0x7f9cd48c4a60>
leading-order electrolyte conductivity <pybamm.models.submodels.electrolyte_conductivity.leading_order_conductivity.LeadingOrder object at 0x7f9cd481fd90>
negative surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x7f9cd481fe20>
positive surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x7f9cd481fe80>
thermal <pybamm.models.submodels.thermal.isothermal.Isothermal object at 0x7f9cd481fdf0>
current collector <pybamm.models.submodels.current_collector.homogeneous_current_collector.Uniform object at 0x7f9cd481ff10>
negative primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd48330d0>
positive primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd4833100>
negative primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd48331f0>
positive primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x7f9cd48331c0>
negative lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x7f9cd4833280>
positive lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x7f9cd48332e0>
total interface <pybamm.models.submodels.interface.total_interfacial_current.TotalInterfacialCurrent object at 0x7f9cd481ff70>

Building the model also sets up the equations, boundary and initial conditions for the model. For example, if we look at model.rhs before building we see that it is empty

[7]:

model.rhs

[7]:

{}

If we try to use this empty model, PyBaMM will give an error. So, before proceeding we must build the model

[8]:

model.build_model()

Now if we look at model.rhs we see that it contains an entry relating to the concentration in each particle, as expected for the SPM

[9]:

model.rhs

[9]:

{Variable(0x3825da4a5fc4eb0b, Discharge capacity [A.h], children=[], domains={}): Multiplication(0x7678edd47e530eec, *, children=['0.0002777777777777778', 'Current function [A]'], domains={}),
 Variable(-0x7fb8d0e6e9632372, Throughput capacity [A.h], children=[], domains={}): Multiplication(-0x7c65e8600b424661, *, children=['0.0002777777777777778', 'abs(Current function [A])'], domains={}),
 Variable(0x69f725db1a464db8, Average negative particle concentration [mol.m-3], children=[], domains={'primary': ['current collector']}): MatrixMultiplication(0xf98a766c86b2483, @, children=['mass(Average negative particle concentration [mol.m-3])', '-3.0 * Current function [A] / (Number of electrodes connected in parallel to make a cell * Electrode width [m] * Electrode height [m]) / Negative electrode thickness [m] / x-average(3.0 * Negative electrode active material volume fraction / Negative particle radius [m]) / Faraday constant [C.mol-1] / x-average(Negative particle radius [m])'], domains={'primary': ['current collector']}),
 Variable(0x48143b39c7603013, X-averaged positive particle concentration [mol.m-3], children=[], domains={'primary': ['positive particle'], 'secondary': ['current collector']}): Divergence(0x17c75a81711ad510, div, children=['Positive particle diffusivity [m2.s-1] * grad(X-averaged positive particle concentration [mol.m-3])'], domains={'primary': ['positive particle'], 'secondary': ['current collector']})}

Now the model can be used in a simulation and solved in the usual way, and we still have access to model defaults such as the default geometry and default spatial methods which are used in the simulation

[10]:

simulation = pybamm.Simulation(model)
simulation.solve([0, 3600])
simulation.plot()

[10]:

<pybamm.plotting.quick_plot.QuickPlot at 0x7f9cd496abe0>

Building a custom model from submodels#

Instead of editing a pre-existing model, you may wish to build your own model from scratch by combining existing submodels of you choice. In this section, we build a Single Particle Model in which the diffusion is assumed infinitely fast in both particles.

To begin, we load a base lithium-ion model. This sets up the basic model structure behind the scenes, and also sets the default parameters to be those corresponding to a lithium-ion battery. Note that the base model does not select any default submodels, so there is no need to pass build=False.

[11]:

model = pybamm.lithium_ion.BaseModel()

Submodels can be added to the model.submodels dictionary in the same way that we changed the submodels earlier.

We use the simplest model for the external circuit, which is the explicit “current control” submodel

[12]:

model.submodels["external circuit"] = pybamm.external_circuit.ExplicitCurrentControl(
    model.param, model.options
)

We want to build a 1D model, so select the Uniform current collector model (if the current collectors are behaving uniformly, then a 1D model is appropriate). We also want the model to be isothermal, so select the thermal model accordingly. Further, we assume that the porosity and active material are constant in space and time.

[13]:

model.submodels["current collector"] = pybamm.current_collector.Uniform(model.param)
model.submodels["thermal"] = pybamm.thermal.isothermal.Isothermal(model.param)
model.submodels["porosity"] = pybamm.porosity.Constant(model.param, model.options)
model.submodels["negative active material"] = pybamm.active_material.Constant(
    model.param, "negative", model.options
)
model.submodels["positive active material"] = pybamm.active_material.Constant(
    model.param, "positive", model.options
)

We assume that the current density varies linearly in the electrodes. This corresponds the the leading-order terms in Ohm’s law in the limit in which the SPM is derived in [3]

[14]:

model.submodels["negative electrode potentials"] = pybamm.electrode.ohm.LeadingOrder(
    model.param, "negative"
)
model.submodels["positive electrode potentials"] = pybamm.electrode.ohm.LeadingOrder(
    model.param, "positive"
)

We assume uniform concentration in both the negative and positive particles. We also have to separately specify a model for the total concentration in each electrode, which is calculated from the concentration in the particles (not a separate ODE)

[15]:

options = {**model.options, "particle": "uniform profile"}
model.submodels["negative primary particle"] = (
    pybamm.particle.XAveragedPolynomialProfile(model.param, "negative", options)
)
model.submodels["positive primary particle"] = (
    pybamm.particle.XAveragedPolynomialProfile(model.param, "positive", options)
)

model.submodels["negative total particle concentration"] = (
    pybamm.particle.TotalConcentration(model.param, "negative", options)
)
model.submodels["positive total particle concentration"] = (
    pybamm.particle.TotalConcentration(model.param, "positive", options)
)

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.

[16]:

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
)
model.submodels["positive interface"] = pybamm.kinetics.InverseButlerVolmer(
    model.param, "positive", "lithium-ion main", options=model.options
)
model.submodels["negative interface current"] = (
    pybamm.kinetics.CurrentForInverseButlerVolmer(
        model.param, "negative", "lithium-ion main"
    )
)
model.submodels["positive interface current"] = (
    pybamm.kinetics.CurrentForInverseButlerVolmer(
        model.param, "positive", "lithium-ion main"
    )
)
model.submodels["negative interface utilisation"] = pybamm.interface_utilisation.Full(
    model.param, "negative", model.options
)
model.submodels["positive interface utilisation"] = pybamm.interface_utilisation.Full(
    model.param, "positive", model.options
)

We don’t want any particle mechanics, SEI formation or lithium plating in this model

[17]:

model.submodels["Negative particle mechanics"] = pybamm.particle_mechanics.NoMechanics(
    model.param, "negative", model.options
)
model.submodels["Positive particle mechanics"] = pybamm.particle_mechanics.NoMechanics(
    model.param, "positive", model.options
)
model.submodels["Negative sei"] = pybamm.sei.NoSEI(
    model.param, "negative", model.options
)
model.submodels["Positive sei"] = pybamm.sei.NoSEI(
    model.param, "positive", model.options
)
model.submodels["Negative sei on cracks"] = pybamm.sei.NoSEI(
    model.param, "negative", model.options, cracks=True
)
model.submodels["Positive sei on cracks"] = pybamm.sei.NoSEI(
    model.param, "positive", model.options, cracks=True
)
model.submodels["Negative lithium plating"] = pybamm.lithium_plating.NoPlating(
    model.param, "Negative"
)
model.submodels["Positive lithium plating"] = pybamm.lithium_plating.NoPlating(
    model.param, "Positive"
)

Finally, for the electrolyte we assume that diffusion is infinitely fast so that the concentration is uniform, and also use the leading-order model for charge conservation, which leads to a linear variation in ionic current in the electrodes

[18]:

model.submodels["electrolyte diffusion"] = (
    pybamm.electrolyte_diffusion.ConstantConcentration(model.param)
)
model.submodels["electrolyte conductivity"] = (
    pybamm.electrolyte_conductivity.LeadingOrder(model.param)
)

Now that we have set all of the submodels we can build the model

[19]:

model.build_model()

We can then use the model in a simulation in the usual way

[20]:

simulation = pybamm.Simulation(model)
simulation.solve([0, 3600])
simulation.plot()

[20]:

<pybamm.plotting.quick_plot.QuickPlot at 0x7f9ccfd7ba30>

References#

The relevant papers for this notebook are:

[21]:

pybamm.print_citations()

[1] Joel A. E. Andersson, Joris Gillis, Greg Horn, James B. Rawlings, and Moritz Diehl. CasADi – A software framework for nonlinear optimization and optimal control. Mathematical Programming Computation, 11(1):1–36, 2019. doi:10.1007/s12532-018-0139-4.
[2] Charles R. Harris, K. Jarrod Millman, Stéfan J. van der Walt, Ralf Gommers, Pauli Virtanen, David Cournapeau, Eric Wieser, Julian Taylor, Sebastian Berg, Nathaniel J. Smith, and others. Array programming with NumPy. Nature, 585(7825):357–362, 2020. doi:10.1038/s41586-020-2649-2.
[3] Scott G. Marquis, Valentin Sulzer, Robert Timms, Colin P. Please, and S. Jon Chapman. An asymptotic derivation of a single particle model with electrolyte. Journal of The Electrochemical Society, 166(15):A3693–A3706, 2019. doi:10.1149/2.0341915jes.
[4] Venkat R. Subramanian, Vinten D. Diwakar, and Deepak Tapriyal. Efficient macro-micro scale coupled modeling of batteries. Journal of The Electrochemical Society, 152(10):A2002, 2005. doi:10.1149/1.2032427.
[5] Valentin Sulzer, Scott G. Marquis, Robert Timms, Martin Robinson, and S. Jon Chapman. Python Battery Mathematical Modelling (PyBaMM). Journal of Open Research Software, 9(1):14, 2021. doi:10.5334/jors.309.