{
"cells": [
{
"cell_type": "markdown",
"id": "0",
"metadata": {},
"source": [
"# Estimating ECM Parameters from Multi-Pulse HPPC Data\n",
"\n",
"This notebook provides example usage for estimating stationary parameters for a two RC branch Thevenin model using multi-pulse HPPC data.\n",
"\n",
"### Setting up the Environment\n",
"\n",
"Before we begin, we need to ensure that we have all the necessary tools. We will install PyBOP from its development branch and upgrade some dependencies:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "1",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"/home/engs2510/Documents/Git/PyBOP/.nox/notebooks-overwrite/bin/python3: No module named pip\r\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Note: you may need to restart the kernel to use updated packages.\n",
"/home/engs2510/Documents/Git/PyBOP/.nox/notebooks-overwrite/bin/python3: No module named pip\r\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Note: you may need to restart the kernel to use updated packages.\n"
]
}
],
"source": [
"%pip install --upgrade pip ipywidgets openpyxl pandas -q\n",
"%pip install pybop -q\n",
"\n",
"import pandas as pd\n",
"import plotly.graph_objects as go\n",
"import pybamm\n",
"\n",
"import pybop\n",
"\n",
"pybop.PlotlyManager().pio.renderers.default = \"notebook_connected\""
]
},
{
"cell_type": "markdown",
"id": "2",
"metadata": {},
"source": [
"In this example, we use the default parameter value for the \"Open-circuit voltage [V] as provided by the original PyBaMM class. The other relevant parameters for the ECM model implementation are updated as per the cell specification."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "3",
"metadata": {},
"outputs": [],
"source": [
"# Load the parameters\n",
"parameter_set = pybop.ParameterSet.pybamm(\"ECM_Example\")\n",
"parameter_set.update(\n",
" {\n",
" \"Cell capacity [A.h]\": 3,\n",
" \"Nominal cell capacity [A.h]\": 3,\n",
" \"Element-1 initial overpotential [V]\": 0,\n",
" \"Upper voltage cut-off [V]\": 4.2,\n",
" \"Lower voltage cut-off [V]\": 2.5,\n",
" \"R0 [Ohm]\": 1e-3,\n",
" \"R1 [Ohm]\": 3e-3,\n",
" \"C1 [F]\": 5e2,\n",
" \"Open-circuit voltage [V]\": pybop.empirical.Thevenin().default_parameter_values[\n",
" \"Open-circuit voltage [V]\"\n",
" ],\n",
" }\n",
")\n",
"# Optional arguments - only needed for two RC pairs\n",
"parameter_set.update(\n",
" {\n",
" \"R2 [Ohm]\": 2e-3,\n",
" \"C2 [F]\": 3e4,\n",
" \"Element-2 initial overpotential [V]\": 0,\n",
" },\n",
" check_already_exists=False,\n",
")"
]
},
{
"cell_type": "markdown",
"id": "4",
"metadata": {},
"source": [
"## Identifying the Parameters\n",
"\n",
"Now that the initial parameter set is constructed, we can start the PyBOP fitting process. First, we define the model class with two RC elements. One important thing here to note is \"maximum solver timestep\" (dt_max) needs to be set correctly to have a good fit."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "5",
"metadata": {},
"outputs": [],
"source": [
"model = pybop.empirical.Thevenin(\n",
" parameter_set=parameter_set,\n",
" options={\"number of rc elements\": 2},\n",
" solver=pybamm.CasadiSolver(mode=\"safe\", dt_max=40),\n",
")"
]
},
{
"cell_type": "markdown",
"id": "6",
"metadata": {},
"source": [
"We use multiple HPPC pulses from the dataset: Kollmeyer, Phillip; Skells, Michael (2020), “Samsung INR21700 30T 3Ah Li-ion Battery Data”, Mendeley Data, V1, doi: 10.17632/9xyvy2njj3.1 "
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "7",
"metadata": {},
"outputs": [],
"source": [
"file_loc = r\"../data/Samsung_INR21700/multipulse_hppc.xlsx\"\n",
"df = pd.read_excel(file_loc, index_col=None, na_values=[\"NA\"])\n",
"df = df.drop_duplicates(subset=[\"Time\"], keep=\"first\")\n",
"\n",
"dataset = pybop.Dataset(\n",
" {\n",
" \"Time [s]\": df[\"Time\"].to_numpy(),\n",
" \"Current function [A]\": df[\"Current\"].to_numpy(),\n",
" \"Voltage [V]\": df[\"Voltage\"].to_numpy(),\n",
" }\n",
")"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "8",
"metadata": {},
"outputs": [],
"source": [
"r0_guess = 0.005\n",
"parameters = pybop.Parameters(\n",
" pybop.Parameter(\n",
" \"R0 [Ohm]\",\n",
" prior=pybop.Gaussian(r0_guess, r0_guess / 10),\n",
" bounds=[0, 0.1],\n",
" ),\n",
" pybop.Parameter(\n",
" \"R1 [Ohm]\",\n",
" prior=pybop.Gaussian(r0_guess, r0_guess / 10),\n",
" bounds=[0, 0.1],\n",
" ),\n",
" pybop.Parameter(\n",
" \"R2 [Ohm]\",\n",
" prior=pybop.Gaussian(r0_guess, r0_guess / 10),\n",
" bounds=[0, 0.1],\n",
" ),\n",
" pybop.Parameter(\n",
" \"C1 [F]\",\n",
" prior=pybop.Gaussian(500, 100),\n",
" bounds=[100, 1000],\n",
" ),\n",
" pybop.Parameter(\n",
" \"C2 [F]\",\n",
" prior=pybop.Gaussian(2000, 500),\n",
" bounds=[1000, 10000],\n",
" ),\n",
")"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "9",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
" \n",
" "
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# To see current vs time profile.\n",
"fig1 = go.Figure()\n",
"# Add a line trace for current vs. time\n",
"fig1.add_trace(\n",
" go.Scatter(\n",
" x=df[\"Time\"].to_numpy(),\n",
" y=df[\"Current\"].to_numpy(),\n",
" mode=\"lines\", # 'lines', 'markers', or 'lines+markers'\n",
" name=\"Current vs Time\",\n",
" )\n",
")\n",
"\n",
"# Customize layout\n",
"fig1.update_layout(\n",
" title=\"Current vs Time\",\n",
" xaxis_title=\"Time (s)\",\n",
" yaxis_title=\"Current (A)\",\n",
" template=\"plotly\", # Use a Plotly template (optional)\n",
")\n",
"\n",
"# Show the plot\n",
"fig1.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "10",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# To see voltage vs time profile.\n",
"fig2 = go.Figure()\n",
"# Add a line trace for current vs. time\n",
"fig2.add_trace(\n",
" go.Scatter(\n",
" x=df[\"Time\"].to_numpy(),\n",
" y=df[\"Voltage\"].to_numpy(),\n",
" mode=\"lines\", # 'lines', 'markers', or 'lines+markers'\n",
" name=\"Voltage vs Time\",\n",
" )\n",
")\n",
"\n",
"# Customize layout\n",
"fig2.update_layout(\n",
" title=\"Voltage vs Time\",\n",
" xaxis_title=\"Time (s)\",\n",
" yaxis_title=\"Voltage (V)\",\n",
" template=\"plotly\", # Use a Plotly template (optional)\n",
")\n",
"\n",
"# Show the plot\n",
"fig2.show()"
]
},
{
"cell_type": "markdown",
"id": "11",
"metadata": {},
"source": [
"The `FittingProblem` class provides us with a single class that holds all of the objects we need to evaluate our selected `SumSquaredError` cost function.\n",
"\n",
"Initial state can be either \"Initial SoC\" or \"Initial open-circuit voltage [V]\". In this example, we get the initial OCV by accessing the voltage data. However, user can simply use a value instead, e.g., {\"Initial open-circuit voltage [V]\": 4.1}. Similarly, if SOC input is required, {\"Initial SoC\": 0.95}."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "12",
"metadata": {},
"outputs": [],
"source": [
"model.build(\n",
" initial_state={\"Initial open-circuit voltage [V]\": df[\"Voltage\"].to_numpy()[0]}\n",
")\n",
"problem = pybop.FittingProblem(\n",
" model,\n",
" parameters,\n",
" dataset,\n",
")\n",
"\n",
"cost = pybop.SumSquaredError(problem)"
]
},
{
"cell_type": "markdown",
"id": "13",
"metadata": {},
"source": [
"Next, we construct the optimisation class with our algorithm of choice and run it. In this case, we select the CMA-ES method as it provides global optimisation capability. For the sake of reducing the runtime of this example, we limit the maximum iterations to 100; however, feel free to update this value."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "14",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/home/engs2510/Documents/Git/PyBOP/.nox/notebooks-overwrite/lib/python3.12/site-packages/pybamm/solvers/base_solver.py:762: SolverWarning:\n",
"\n",
"Explicit interpolation times not implemented for CasADi solver with 'safe' mode\n",
"\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Initial parameters: [6.23380007e-03 4.72246102e-03 4.76953108e-03 4.20430114e+02\n",
" 2.15280196e+03]\n",
"Estimated parameters: [1.00042463e-02 3.37755166e-03 1.89253796e-02 3.74961536e+02\n",
" 1.65220599e+03]\n"
]
}
],
"source": [
"optim = pybop.XNES(\n",
" cost,\n",
" sigma0=[1e-3, 1e-3, 1e-3, 10, 10],\n",
" max_unchanged_iterations=20,\n",
" max_iterations=100,\n",
")\n",
"x, final_cost = optim.run()\n",
"print(\"Initial parameters:\", optim.x0)\n",
"print(\"Estimated parameters:\", x)"
]
},
{
"cell_type": "markdown",
"id": "15",
"metadata": {},
"source": [
"## Plotting and Visualisation\n",
"\n",
"PyBOP provides various plotting utilities to visualize the results of the optimisation."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "16",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"pybop.quick_plot(problem, problem_inputs=x, title=\"Optimised Comparison\");"
]
},
{
"cell_type": "markdown",
"id": "17",
"metadata": {},
"source": [
"### Convergence and Parameter Trajectories\n",
"\n",
"To assess the optimisation process, we can plot the convergence of the cost function and the trajectories of the parameters:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "18",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
""
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"pybop.plot_convergence(optim)\n",
"pybop.plot_parameters(optim);"
]
},
{
"cell_type": "markdown",
"id": "19",
"metadata": {},
"source": [
"### Conclusion\n",
"\n",
"This notebook illustrates how to perform parameter estimation for multi-pulse HPPC data, providing insights into the optimisation process through various visualisations."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"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.12.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}