[Bug]: Uncertain compatibility between pybamm.lithium_ion.SPMe() and Marquis 2019 paper

[ejfdickinson]

PyBaMM Version

24.1

Python Version

3.9.16

Describe the bug

pybamm.lithium_ion.SPMe claims equivalence with the Marquis 2019 research paper.

However, the voltage output from PyBaMM's model does not seem to equate with the expression of cell voltage in the Marquis SPMe model (Marquis 2019, eq. [38]). Rather, cell voltage is still reported as the difference of the current collector potentials (as defined by Marquis 2019 for the DFN model, eq. [5]). This excludes the integration between eqs. [9] and [10], so it gives a different result from the one implied by the Marquis 2019 model.

• Is there a variable that gives access to the cell voltage as given in the Marquis 2019 model?
• Or is, it not in fact the case that the SPMe implementation is compatible with Marquis 2019, in full?

Steps to Reproduce

The below code snippet demonstrates that the cell voltage output of the pybamm.lithium_ion.SPMe model equates to the difference of the current collector potentials, as in the DFN model. This is not equivalent to the specification of the Marquis 2019 model (eq. 38), which incorporates averages over each electrode.

import pybamm
import numpy as np

parameter_values = pybamm.ParameterValues("Chen2020")
model = pybamm.lithium_ion.SPMe()
experiment = pybamm.Experiment(["Discharge at 1C until 3 V"])

sim = pybamm.Simulation(
    parameter_values=parameter_values,
    model=model,
    experiment=experiment
)
sol = sim.solve()

print(
    np.max(
        np.abs(sol["Voltage [V]"].entries - (sol["Positive current collector potential [V]"].entries-sol["Negative current collector potential [V]"].entries))
    )
)

Relevant log output

No response

[kawaMANMI]

Hi @ejfdickinson , I am not sure how familiar you are with this note

https://docs.pybamm.org/en/v23.5_a/source/examples/notebooks/plotting/plot-voltage-components.html

to get the voltage component."

[ejfdickinson]

@kawaMANMI

Thanks for the note.

The issue I raised is that the equation at the very top of that note (physical interpretation: "the voltage of a cell is given by the difference in solid-phase potentials at the two edges of the domain") is applicable only to the DFN case.

It's not how cell voltage is computed in the Marquis 2019 asymptotic SPMe model, because of the way that the SPMe derivation averages the particles of each electrode through-thickness to a single particle, while retaining macroscopic resolution of the electrolyte variables.

If the equation from the voltage components note is how pybamm.lithium_ion.SPMe computes cell voltage - and it is certainly what is placed in the "Voltage [V]" global output - then the implemented model doesn't match the paper. That's a documentation bug, and maybe an implementation bug if the specified intention was to implement SPMe according to Marquis 2019. It really depends what PyBaMM's SPMe is supposed to do. Is that set out in equations anywhere?

Or, if we can get cell voltage defined properly for the SPMe, but it's not in "Voltage [V]", in what solution variable is it accessible?

[ejfdickinson]

Perhaps another way of thinking of it is:

• the cell voltage is always the difference of current collector potentials (from physical common sense)
• but in an asymptotically derived SPMe, the positive current collector potential is not evaluated as $\left.\phi_\mathrm{s,p}\right|_{x=L}$, because we are averaging over the positive electrode particles (contacting the electric conductor) rather than treating them as spatially resolved.

[valentinsulzer]

This may not be made explicit enough in Marquis2019, but we account for the difference between "average solid-phase electrode potential" (phi_av) and "solid-phase potential at the current collector interface" (phi_cc) by assuming a quadratic solid-phase potential in each electrode (uniform interfacial current density -> linear current density -> quadratic potential). See https://github.com/pybamm-team/PyBaMM/blob/84eb098e3140d774e39ab3d3c842c22462faf4aa/pybamm/models/submodels/electrode/ohm/composite_ohm.py#L52-L66 for the code where this is implemented.

In the negative, phi_cc = 0 and phi_av is calculated based on that. In the positive, phi_av is calculated based on averaging over the electrode and solving / inverting BV, and phi_cc is calculated proportionally to phi_av. The difference between phi_av and phi_cc ("solid phase ohmic losses") comes out as IL/3sigma, which is equations (56)-(57) in Marquis 2019.

[ejfdickinson]

@tinosulzer Could you clarify what model pybamm.lithium_ion.SPMe actually implements? Is it the canonical SPMe in eqs [40] of Marquis 2019?

If so, should I interpret your explanation as being that I should expect the PyBaMM Terminal voltage [V] output to obey Marquis 2019 [40j], and that it should also equate to the difference of current collector potentials?

edit: Or is it the SPMe(S) collated in eqs [39], and hence Terminal voltage [V] should obey [38]?

[valentinsulzer]

The canonical SPMe, with caveats explained in https://ora.ox.ac.uk/objects/uuid:f659e447-d78e-4bcb-9c1a-24d870984380/files/sft848q83t

But the only difference between the canonical SPMe and the SPMe(S) is the treatment of the electrolyte concentration, explained in the last paragraph of section 3. Equations (28) are the relevant ones for your question

[ejfdickinson]

Thanks @tinosulzer - I'll try to find time to look closely next week, and get back to you if I still feel there are unresolved issues.

As a practical follow-up, it would be great to update the SPMe docs to say exactly what model and equations from the referenced paper are implemented.

The wider context here is that we've had great difficulty reproducing the PyBaMM results for SPMe in independent modelling environments. However, it's been quite ambiguous what is meant by "from Marquis 2019" in terms of the actual equation set solved, and I haven't been confident reverse-engineering this from the PyBaMM source. With your input I think we can be sure we're comparing like-with-like.

[valentinsulzer]

Depending on what you can do in the other modeling environments I would suggest comparing other variables (concentrations, potentials, etc) to pinpoint the source of the discrepancy. And changing parameters to absurd values to eliminate or accentuate gradients and overpotentials to see differences

[valentinsulzer]

Closing due to inactivity, feel free to reopen if new issues come up