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.
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
[notice] A new release of pip is available: 25.1.1 -> 25.3
[notice] To update, run: pip install --upgrade pip
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 0x31f612840>
discharge and throughput variables <pybamm.models.submodels.external_circuit.discharge_throughput.DischargeThroughput object at 0x31f6128a0>
porosity <pybamm.models.submodels.porosity.constant_porosity.Constant object at 0x31c710620>
Negative interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x31e326990>
Positive interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x31e754260>
negative primaryparticle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x31f6127e0>
positive primaryparticle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x31cc8a300>
negative primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x31f612a20>
positive primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x31dc9a0c0>
electrolyte transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman.Bruggeman object at 0x31ea3dd90>
electrode transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman.Bruggeman object at 0x31e9dca70>
transverse convection <pybamm.models.submodels.convection.transverse.no_convection.NoConvection object at 0x31e3b2e70>
through-cell convection <pybamm.models.submodels.convection.through_cell.no_convection.NoConvection object at 0x31da78620>
negative primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x31b161b20>
positive primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x31b161c10>
negative interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x16a1f8770>
negative interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.base_inverse.CurrentForInverseKinetics object at 0x31d4cb500>
positive interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x31b190bc0>
positive interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.base_inverse.CurrentForInverseKinetics object at 0x31e3b2a80>
negative primary particle <pybamm.models.submodels.particle.fickian_diffusion.FickianDiffusion object at 0x10b6c10a0>
negative primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x31f6129c0>
positive primary particle <pybamm.models.submodels.particle.fickian_diffusion.FickianDiffusion object at 0x31d4f1730>
positive primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x31b1633e0>
negative electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x31e3260c0>
positive electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x31e80c2f0>
electrolyte diffusion <pybamm.models.submodels.electrolyte_diffusion.constant_concentration.ConstantConcentration object at 0x31b190b30>
leading-order electrolyte conductivity <pybamm.models.submodels.electrolyte_conductivity.leading_order_conductivity.LeadingOrder object at 0x31e756540>
negative surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x10b6c2510>
positive surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x31c21e270>
thermal <pybamm.models.submodels.thermal.isothermal.Isothermal object at 0x10b6c1040>
surface temperature <pybamm.models.submodels.thermal.surface.ambient.Ambient object at 0x31e39d0d0>
current collector <pybamm.models.submodels.current_collector.homogeneous_current_collector.Uniform object at 0x31b190890>
negative primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f5f6c30>
negative primary sei thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31e05f980>
positive primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f5f6180>
positive primary sei thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f6128d0>
negative primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f6127b0>
negative primary sei on cracks thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f612870>
positive primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f612810>
positive primary sei on cracks thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f612780>
negative primary lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x31f612690>
positive primary lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x31dbd4080>
total interface <pybamm.models.submodels.interface.total_interfacial_current.TotalInterfacialCurrent object at 0x31e7dcc80>
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 0x31f746000>
discharge and throughput variables <pybamm.models.submodels.external_circuit.discharge_throughput.DischargeThroughput object at 0x31f745fa0>
porosity <pybamm.models.submodels.porosity.constant_porosity.Constant object at 0x31f746030>
Negative interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x10b49dd90>
Positive interface utilisation <pybamm.models.submodels.interface.interface_utilisation.full_utilisation.Full object at 0x31f6f8320>
negative primaryparticle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x31f745e20>
positive primaryparticle mechanics <pybamm.models.submodels.particle_mechanics.no_mechanics.NoMechanics object at 0x31f745c40>
negative primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x31e373d40>
positive primary active material <pybamm.models.submodels.active_material.constant_active_material.Constant object at 0x31f745b20>
electrolyte transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman.Bruggeman object at 0x10b6c17f0>
electrode transport efficiency <pybamm.models.submodels.transport_efficiency.bruggeman.Bruggeman object at 0x31b160fb0>
transverse convection <pybamm.models.submodels.convection.transverse.no_convection.NoConvection object at 0x10b6c24b0>
through-cell convection <pybamm.models.submodels.convection.through_cell.no_convection.NoConvection object at 0x31f745a90>
negative primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x31f7457f0>
positive primary open-circuit potential <pybamm.models.submodels.interface.open_circuit_potential.single_ocp.SingleOpenCircuitPotential object at 0x31f7458b0>
negative interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x31f7459a0>
negative interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.base_inverse.CurrentForInverseKinetics object at 0x31f745760>
positive interface <pybamm.models.submodels.interface.kinetics.inverse_kinetics.inverse_butler_volmer.InverseButlerVolmer object at 0x12dc79a90>
positive interface current <pybamm.models.submodels.interface.kinetics.inverse_kinetics.base_inverse.CurrentForInverseKinetics object at 0x31f716180>
negative primary particle <pybamm.models.submodels.particle.x_averaged_polynomial_profile.XAveragedPolynomialProfile object at 0x31f6a1490>
negative primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x31f745430>
positive primary particle <pybamm.models.submodels.particle.fickian_diffusion.FickianDiffusion object at 0x31f7452b0>
positive primary total particle concentration <pybamm.models.submodels.particle.total_particle_concentration.TotalConcentration object at 0x31f7453d0>
negative electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x31f7456a0>
positive electrode potential <pybamm.models.submodels.electrode.ohm.leading_ohm.LeadingOrder object at 0x31f744fe0>
electrolyte diffusion <pybamm.models.submodels.electrolyte_diffusion.constant_concentration.ConstantConcentration object at 0x31f6c8d70>
leading-order electrolyte conductivity <pybamm.models.submodels.electrolyte_conductivity.leading_order_conductivity.LeadingOrder object at 0x31f745370>
negative surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x31f7451f0>
positive surface potential difference <pybamm.models.submodels.electrolyte_conductivity.surface_potential_form.explicit_surface_form_conductivity.Explicit object at 0x31e755130>
thermal <pybamm.models.submodels.thermal.isothermal.Isothermal object at 0x31f745280>
surface temperature <pybamm.models.submodels.thermal.surface.ambient.Ambient object at 0x31f745070>
current collector <pybamm.models.submodels.current_collector.homogeneous_current_collector.Uniform object at 0x31f745160>
negative primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f744470>
negative primary sei thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f744320>
positive primary sei <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f744230>
positive primary sei thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f744140>
negative primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f746060>
negative primary sei on cracks thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f6ca900>
positive primary sei on cracks <pybamm.models.submodels.interface.sei.no_sei.NoSEI object at 0x31f7460f0>
positive primary sei on cracks thickness <pybamm.models.submodels.interface.sei.sei_thickness.SEIThickness object at 0x31f746120>
negative primary lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x31f733e90>
positive primary lithium plating <pybamm.models.submodels.interface.lithium_plating.no_plating.NoPlating object at 0x31f733920>
total interface <pybamm.models.submodels.interface.total_interfacial_current.TotalInterfacialCurrent object at 0x31f716690>
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(-0x66e7b96f57f14d08, Discharge capacity [A.h], children=[], domains={}): Multiplication(0x48c07e1071c1540, *, children=['0.0002777777777777778', 'Current function [A]'], domains={}),
Variable(-0x773248a54645943c, Throughput capacity [A.h], children=[], domains={}): Multiplication(-0x5800ae8c72a1768c, *, children=['0.0002777777777777778', 'abs(Current function [A])'], domains={}),
Variable(-0x7027188ff58ae776, Average negative particle concentration [mol.m-3], children=[], domains={'primary': ['current collector']}): MatrixMultiplication(-0x228b804158130f33, @, children=['mass(Average negative particle concentration [mol.m-3])', '-3.0 * Current function [A] / (Electrode width [m] * Electrode height [m] * Number of electrodes connected in parallel to make a cell) / 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(-0x3f6460d6557d6997, X-averaged positive particle concentration [mol.m-3], children=[], domains={'primary': ['positive particle'], 'secondary': ['current collector']}): Divergence(-0x1006b5f1f45d4611, 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 0x3280195b0>
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.CurrentForInverseKinetics(
model.param, "negative", "lithium-ion main"
)
)
model.submodels["positive interface current"] = (
pybamm.kinetics.CurrentForInverseKinetics(
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 0x329a126f0>
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] Von DAG Bruggeman. Berechnung verschiedener physikalischer konstanten von heterogenen substanzen. i. dielektrizitätskonstanten und leitfähigkeiten der mischkörper aus isotropen substanzen. Annalen der physik, 416(7):636–664, 1935.
[3] 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.
[4] Alan C. Hindmarsh. The PVODE and IDA algorithms. Technical Report, Lawrence Livermore National Lab., CA (US), 2000. doi:10.2172/802599.
[5] Alan C. Hindmarsh, Peter N. Brown, Keith E. Grant, Steven L. Lee, Radu Serban, Dan E. Shumaker, and Carol S. Woodward. SUNDIALS: Suite of nonlinear and differential/algebraic equation solvers. ACM Transactions on Mathematical Software (TOMS), 31(3):363–396, 2005. doi:10.1145/1089014.1089020.
[6] 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.
[7] 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.
[8] 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.