Plett OCP Hysteresis Model

CHANGELOG.md

@@ -14,6 +14,7 @@
 - Renamed "electrode diffusivity" to "particle diffusivity" as a non-breaking change with a deprecation warning ([#3624](https://github.com/pybamm-team/PyBaMM/pull/3624))
 - Add support for BPX version 0.4.0 which allows for blended electrodes and user-defined parameters in BPX([#3414](https://github.com/pybamm-team/PyBaMM/pull/3414))
 - Added `by_submodel` feature in `print_parameter_info` method to allow users to print parameters and types of submodels in a tabular and readable format ([#3628](https://github.com/pybamm-team/PyBaMM/pull/3628))
+- Added `WyciskOpenCircuitPotential` for differential capacity hysteresis state open-circuit potential submodel ([#3593](https://github.com/pybamm-team/PyBaMM/pull/3593))
 
 ## Bug Fixes
 

docs/source/api/models/submodels/interface/open_circuit_potential/index.rst

@@ -7,3 +7,4 @@ Open-circuit potential models
   current_sigmoid_ocp
   single_ocp
   msmr_ocp
+  wycisk_ocp

docs/source/api/models/submodels/interface/open_circuit_potential/wycisk_ocp.rst

@@ -0,0 +1,5 @@
+Wycisk Open Circuit Potential
+=============================
+
+.. autoclass:: pybamm.open_circuit_potential.WyciskOpenCircuitPotential
+    :members:

docs/source/examples/index.rst

@@ -50,6 +50,7 @@ The notebooks are organised into subfolders, and can be viewed in the galleries
     notebooks/models/compare-particle-diffusion-models.ipynb
     notebooks/models/composite_particle.ipynb
     notebooks/models/coupled-degradation.ipynb
+    notebooks/models/differential-capacity-hysteresis-state.ipynb
     notebooks/models/DFN-with-particle-size-distributions.ipynb
     notebooks/models/DFN.ipynb
     notebooks/models/electrode-state-of-health.ipynb

docs/source/examples/notebooks/models/differential-capacity-hysteresis-state.ipynb

@@ -0,0 +1,269 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Differential Capacity Hysteresis State model"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Note: you may need to restart the kernel to use updated packages.\n"
+     ]
+    }
+   ],
+   "source": [
+    "%pip install \"pybamm[plot,cite]\" -q    # install PyBaMM if it is not installed\n",
+    "import pybamm"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Model Equations\n",
+    "\n",
+    "Herein the model equations for the Differential Capacity Hysteresis State open-circuit potential model are outlined, as described in Wycisk (2022).\n",
+    "\n",
+    "### Hysteresis State Variable\n",
+    "\n",
+    "This approach utilizes a state variable to represent the degree of hysteresis at a given time and stoichiometry, $h(z,t)$. The hysteresis is treated separately from the open-circuit potential, where the potential of the electrode is written as\n",
+    "\n",
+    "$$ U = U_{avg}^0(z) + H(z) \\cdot h(z,t) - \\eta $$\n",
+    "\n",
+    "Where $H(z)$ is a function representing the hysteresis as a function of stoichiometry, $z$, and where $\\eta$ represents the sum of the overpotentials. $U_{avg}^0(z)$ is simply the average of the delithiation and lithiation open-circuit potential branches. $H(z)$ can be determined by finding the half-difference value between the lithiation and delithiation branches across the entire stoichiometry range. The state variable $h(z,t)$ is both stoichiometry and time-dependant, and spans between the range of -1 and 1. The hysteresis state variable $h(z,t)$ can be expressed in differential form with respect to time as\n",
+    "\n",
+    "$$ \\frac{dh(z,t)}{dt} = \\left(\\frac{k(z) \\cdot I(t)}{Q_{cell}}\\right)\\left(1-\\text{sgn}\\left(\\frac{dz(t)}{dt}\\right) h(z,t)\\right) $$\n",
+    "\n",
+    "where $ k(z) $ is expressed as \n",
+    "\n",
+    "$$ k(z) = K \\cdot \\frac{1}{\\left(C_{diff}\\left(z\\right)\\right)^{x}} $$\n",
+    "\n",
+    "And where $C_{diff}(z)$ is the differential capacity with respect to potential, expressed as \n",
+    "\n",
+    "$$ C_{diff}(z) = \\frac{dQ}{dU_{avg}^0(z)} $$\n",
+    "\n",
+    "Here, $Q$ is the capacity of the phase or active material experiencing the voltage hysteresis. The remaining parameters are $K$ and $x$ which are both fitting parameters that affect the response of the hysteresis state decay when passing charge in either direction.\n",
+    "\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Comparing the DCHS and Current-Sigmoid model approaches\n",
+    "\n",
+    "The behavior of the DCHS model is different than the current-sigmoid model approach for open-circuit potential in systems with hysteresis. Where the current-sigmoid model switches between hysteresis states simply based on the instantaneous current, the DCHS model switches based on the amount of charge passed through the active material phase while also relying on the previous hysteresis state. To assess this differentiated performance, we will compare it to the current-sigmoid model by adapting the Chen2020_composite parameter set.\n",
+    "\n",
+    "First we generate the model, and specify the open-circuit potential methods for the negative and positive electrodes. To maintain consistency with the parameter set, two phases for the negative electrode will be defined."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "model_DCHS = pybamm.lithium_ion.DFN(\n",
+    "    {\n",
+    "        \"open-circuit potential\": ((\"single\", \"Wycisk\"), \"single\"),\n",
+    "        \"particle phases\": (\"2\", \"1\"),\n",
+    "    }\n",
+    ")\n",
+    "\n",
+    "model_current_sigmoid = pybamm.lithium_ion.DFN(\n",
+    "    {\n",
+    "        \"open-circuit potential\": ((\"single\", \"current sigmoid\"), \"single\"),\n",
+    "        \"particle phases\": (\"2\", \"1\"),\n",
+    "    }\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Next, lets define the modifications to the parameter set"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "parameters_DCHS = pybamm.ParameterValues(\"Chen2020_composite\")\n",
+    "parameters_current_sigmoid = pybamm.ParameterValues(\"Chen2020_composite\")\n",
+    "\n",
+    "\n",
+    "# get the lithiation and delithiation functions\n",
+    "lithiation_ocp = parameters_DCHS[\"Secondary: Negative electrode lithiation OCP [V]\"]\n",
+    "delithiation_ocp = parameters_DCHS[\"Secondary: Negative electrode delithiation OCP [V]\"]\n",
+    "\n",
+    "\n",
+    "# define an additional OCP function\n",
+    "def ocp_avg(sto):\n",
+    "    return (lithiation_ocp(sto) + delithiation_ocp(sto)) / 2\n",
+    "\n",
+    "\n",
+    "# add additional parameters\n",
+    "parameters_DCHS.update(\n",
+    "    {\n",
+    "        \"Secondary: Negative electrode OCP [V]\": ocp_avg,\n",
+    "    },\n",
+    "    check_already_exists=False,\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Next, we need to add the additional parameters required by the model"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "parameters_DCHS.update(\n",
+    "    {\n",
+    "        \"Secondary: Negative particle hysteresis decay rate\": 0.005,\n",
+    "        \"Secondary: Negative particle hysteresis switching factor\": 10,\n",
+    "    },\n",
+    "    check_already_exists=False,\n",
+    ")\n",
+    "\n",
+    "\n",
+    "experiment = pybamm.Experiment(\n",
+    "    [\n",
+    "        (\"Discharge at 1 C for 1 hour or until 2.5 V\", \"Rest for 15 minutes\"),\n",
+    "        (\n",
+    "            \"Charge at 1C until 4.2 V\",\n",
+    "            \"Hold at 4.2 V until 0.05 C\",\n",
+    "            \"Rest for 15 minutes\",\n",
+    "        ),\n",
+    "    ]\n",
+    ")\n",
+    "\n",
+    "\n",
+    "simulation_dchs = pybamm.Simulation(\n",
+    "    model_DCHS, experiment=experiment, parameter_values=parameters_DCHS\n",
+    ")\n",
+    "solution_dchs = simulation_dchs.solve(calc_esoh=False)\n",
+    "\n",
+    "simulation_current_sigmoid = pybamm.Simulation(\n",
+    "    model_current_sigmoid,\n",
+    "    experiment=experiment,\n",
+    "    parameter_values=parameters_current_sigmoid,\n",
+    ")\n",
+    "\n",
+    "solution_current_sigmoid = simulation_current_sigmoid.solve(calc_esoh=False)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Now plotting the results and the hysteresis state "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "application/vnd.jupyter.widget-view+json": {
+       "model_id": "7fa481d981ff4f9da3b688120bdbf930",
+       "version_major": 2,
+       "version_minor": 0
+      },
+      "text/plain": [
+       "interactive(children=(FloatSlider(value=0.0, description='t', max=3.1492654802910014, step=0.03149265480291001…"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "output_variables = [\n",
+    "    \"X-averaged negative electrode secondary hysteresis state\",\n",
+    "    \"Negative electrode secondary open-circuit potential [V]\",\n",
+    "    \"Negative electrode secondary stoichiometry\",\n",
+    "    \"Terminal voltage [V]\",\n",
+    "    \"X-averaged negative electrode secondary open-circuit potential [V]\",\n",
+    "]\n",
+    "\n",
+    "pybamm.QuickPlot(solution_dchs, output_variables=output_variables).dynamic_plot()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "application/vnd.jupyter.widget-view+json": {
+       "model_id": "6745f5d4c7e64bbc8a57011b282bdf2e",
+       "version_major": 2,
+       "version_minor": 0
+      },
+      "text/plain": [
+       "interactive(children=(FloatSlider(value=0.0, description='t', max=3.1492654802910014, step=0.03149265480291001…"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "output_variables = [\n",
+    "    \"Terminal voltage [V]\",\n",
+    "    \"Current [A]\",\n",
+    "    \"Negative electrode secondary open-circuit potential [V]\",\n",
+    "]\n",
+    "pybamm.QuickPlot(\n",
+    "    [solution_current_sigmoid, solution_dchs],\n",
+    "    labels=[\"Current sigmoid\", \"DCHS\"],\n",
+    "    output_variables=output_variables,\n",
+    ").dynamic_plot()"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "venv",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.11.0"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

pybamm/CITATIONS.bib

@@ -691,3 +691,16 @@ @article{landesfeind2019temperature
   year={2019},
   publisher={The Electrochemical Society}
 }
+
+@article{Wycisk2022,
+title = {Modified Plett-model for modeling voltage hysteresis in lithium-ion cells},
+journal = {Journal of Energy Storage},
+volume = {52},
+pages = {105016},
+year = {2022},
+issn = {2352-152X},
+doi = {https://doi.org/10.1016/j.est.2022.105016},
+url = {https://www.sciencedirect.com/science/article/pii/S2352152X22010192},
+author = {Dominik Wycisk and Marc Oldenburger and Marc Gerry Stoye and Toni Mrkonjic and Arnulf Latz},
+keywords = {Lithium-ion battery, Voltage hysteresis, Plett-model, Silicon–graphite anode},
+}

pybamm/models/full_battery_models/base_battery_model.py

@@ -102,7 +102,7 @@ class BatteryModelOptions(pybamm.FuzzyDict):
                 reactions.
             * "open-circuit potential" : str
                 Sets the model for the open circuit potential. Can be "single"
-                (default), "current sigmoid", or "MSMR". If "MSMR" then the "particle"
+                (default), "current sigmoid", "Wycisk", or "MSMR". If "MSMR" then the "particle"
                 option must also be "MSMR". A 2-tuple can be provided for different
                 behaviour in negative and positive electrodes.
             * "operating mode" : str
@@ -263,7 +263,7 @@ def __init__(self, extra_options):
                 "stress and reaction-driven",
             ],
             "number of MSMR reactions": ["none"],
-            "open-circuit potential": ["single", "current sigmoid", "MSMR"],
+            "open-circuit potential": ["single", "current sigmoid", "MSMR", "Wycisk"],
             "operating mode": [
                 "current",
                 "voltage",

pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py

@@ -248,6 +248,9 @@ def set_open_circuit_potential_submodel(self):
                     ocp_model = ocp_submodels.SingleOpenCircuitPotential
                 elif ocp_option == "current sigmoid":
                     ocp_model = ocp_submodels.CurrentSigmoidOpenCircuitPotential
+                elif ocp_option == "Wycisk":
+                    pybamm.citations.register("Wycisk2022")
+                    ocp_model = ocp_submodels.WyciskOpenCircuitPotential
                 elif ocp_option == "MSMR":
                     ocp_model = ocp_submodels.MSMROpenCircuitPotential
                 self.submodels[f"{domain} {phase} open-circuit potential"] = ocp_model(

pybamm/models/submodels/interface/open_circuit_potential/__init__.py

@@ -2,3 +2,4 @@
 from .single_ocp import SingleOpenCircuitPotential
 from .current_sigmoid_ocp import CurrentSigmoidOpenCircuitPotential
 from .msmr_ocp import MSMROpenCircuitPotential
+from .wycisk_ocp import WyciskOpenCircuitPotential