{ "cells": [ { "cell_type": "markdown", "metadata": { "id": "expmkveO04pw" }, "source": [ "## Parameter identification for various models\n", "\n", "To investigate the performance of parameter identification for different electrochemical models we will start with synthetic data from the highest order model in PyBOP (Many-particle DFN) and try to identify the correct parameter values on the reduced order models.\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 and upgrade dependencies:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "X87NUGPW04py", "outputId": "0d785b07-7cff-4aeb-e60a-4ff5a669afbf" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "/Users/engs2510/Documents/Git/Second_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" ] }, { "name": "stdout", "output_type": "stream", "text": [ "/Users/engs2510/Documents/Git/Second_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 pybamm -q\n", "%pip install pybop -q" ] }, { "cell_type": "markdown", "metadata": { "id": "jAvD5fk104p0" }, "source": [ "Next, we import the added packages plus any additional dependencies," ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "SQdt4brD04p1" }, "outputs": [], "source": [ "import numpy as np\n", "import pybamm\n", "\n", "import pybop\n", "\n", "go = pybop.plot.PlotlyManager().go\n", "pybop.plot.PlotlyManager().pio.renderers.default = \"notebook_connected\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's fix the random seed in order to generate consistent output during development, although this does not need to be done in practice." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "np.random.seed(8)" ] }, { "cell_type": "markdown", "metadata": { "id": "X8-tubYY04p_" }, "source": [ "## Optimising the Parameters" ] }, { "cell_type": "markdown", "metadata": { "id": "PQqhvSZN04p_" }, "source": [ "First, we define the model to be used for the parameter optimisation," ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "zuvGHWID04p_" }, "outputs": [], "source": [ "parameter_set = pybop.ParameterSet.pybamm(\"Chen2020\")\n", "parameter_set = pybamm.get_size_distribution_parameters(parameter_set)\n", "synth_model = pybop.lithium_ion.DFN(\n", " parameter_set=parameter_set, options={\"particle size\": \"distribution\"}\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Simulating Forward Model\n", "\n", "We can then simulate the model using the `predict` method, with a default constant current to generate voltage data." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "n_points = 450\n", "t_eval = np.linspace(0, 1600 + 1000, n_points)\n", "current = np.concatenate(\n", " [np.ones(200) * parameter_set[\"Nominal cell capacity [A.h]\"], np.zeros(250)]\n", ")\n", "initial_state = {\"Initial SoC\": 0.5}" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "dataset = pybop.Dataset(\n", " {\n", " \"Time [s]\": t_eval,\n", " \"Current function [A]\": current,\n", " }\n", ")" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "synth_model.build(dataset=dataset, initial_state=initial_state)\n", "synth_model.signal = [\"Voltage [V]\"]\n", "values = synth_model.simulate(t_eval=t_eval, inputs={})" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Adding Noise to Voltage Data\n", "\n", "To make the parameter estimation more realistic, we add Gaussian noise to the data." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [ { "data": { "text/html": [ " \n", " " ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "sigma = 0.001\n", "corrupt_values = values[\"Voltage [V]\"].data + np.random.normal(\n", " 0, sigma, len(values[\"Voltage [V]\"].data)\n", ")\n", "go.Figure(\n", " data=go.Scatter(x=t_eval, y=corrupt_values, mode=\"lines\"),\n", " layout=go.Layout(title=\"Corrupted Voltage\", width=800, height=600),\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Identifying the Parameters" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will now set up the parameter estimation process by defining the datasets for optimisation and selecting the model parameters we wish to estimate." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Creating a Dataset\n", "\n", "The dataset for optimisation is composed of time, current, and the noisy voltage data:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "dataset = pybop.Dataset(\n", " {\n", " \"Time [s]\": t_eval,\n", " \"Current function [A]\": current,\n", " \"Voltage [V]\": corrupt_values,\n", " }\n", ")" ] }, { "cell_type": "markdown", "metadata": { "id": "ffS3CF_704qA" }, "source": [ "Next, we define the model parameters for optimisation. Furthermore, PyBOP provides functionality to define a prior for the parameters. The initial parameter values used in the optimisation will be randomly drawn from the prior distribution." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "WPCybXIJ04qA" }, "outputs": [], "source": [ "parameters = pybop.Parameters(\n", " pybop.Parameter(\n", " \"Positive electrode thickness [m]\",\n", " prior=pybop.Gaussian(7.56e-05, 0.05e-05),\n", " bounds=[65e-06, 85e-06],\n", " true_value=parameter_set[\"Positive electrode thickness [m]\"],\n", " ),\n", " pybop.Parameter(\n", " \"Negative electrode thickness [m]\",\n", " prior=pybop.Gaussian(8.52e-05, 0.05e-05),\n", " bounds=[75e-06, 95e-06],\n", " true_value=parameter_set[\"Negative electrode thickness [m]\"],\n", " ),\n", ")" ] }, { "cell_type": "markdown", "metadata": { "id": "n4OHa-aF04qA" }, "source": [ "We can now define the output signal, the problem (which combines the model with the dataset) and construct a cost function which in this example is the `GravimetricEnergyDensity()` used to maximise the gravimetric energy density of the cell." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "models = [\n", " pybop.lithium_ion.SPM(parameter_set=parameter_set),\n", " pybop.lithium_ion.SPMe(parameter_set=parameter_set),\n", "]" ] }, { "cell_type": "markdown", "metadata": { "id": "eQiGurUV04qB" }, "source": [ "Let's construct PyBOP's optimisation class for each model. This class provides the methods needed to fit the forward model. For this example, we use an evolution strategy (XNES) as the optimiser." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "etMzRtx404qA" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Running Single Particle Model\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "Halt: No significant change for 15 iterations.\n", "OptimisationResult:\n", " Initial parameters: [7.61508933e-05 8.59364316e-05]\n", " Optimised parameters: [6.50355739e-05 8.10453665e-05]\n", " Final cost: 0.33026298961068246\n", " Optimisation time: 22.5766921043396 seconds\n", " Number of iterations: 40\n", " SciPy result available: No\n", "Running Single Particle Model with Electrolyte\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "Halt: No significant change for 15 iterations.\n", "OptimisationResult:\n", " Initial parameters: [7.61508933e-05 8.59364316e-05]\n", " Optimised parameters: [7.08088980e-05 8.53871336e-05]\n", " Final cost: 0.005807228768217647\n", " Optimisation time: 22.70540690422058 seconds\n", " Number of iterations: 23\n", " SciPy result available: No\n" ] } ], "source": [ "optims = []\n", "xs = []\n", "for model in models:\n", " print(f\"Running {model.name}\")\n", " model.set_initial_state(initial_state)\n", " problem = pybop.FittingProblem(model, parameters, dataset)\n", " cost = pybop.SumSquaredError(problem)\n", " optim = pybop.XNES(\n", " cost, verbose=True, max_iterations=60, max_unchanged_iterations=15\n", " )\n", " results = optim.run()\n", " optims.append(optim)\n", " xs.append(results.x)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "N3FtAhrT04qB" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "| Model: Single Particle Model | Results: [6.50355739e-05 8.10453665e-05] |\n", "| Model: Single Particle Model with Electrolyte | Results: [7.08088980e-05 8.53871336e-05] |\n" ] } ], "source": [ "for optim, x in zip(optims, xs):\n", " print(f\"| Model: {optim.cost.problem.model.name} | Results: {x} |\")" ] }, { "cell_type": "markdown", "metadata": { "id": "KxKURtH704qC" }, "source": [ "## Plotting and Visualisation\n", "\n", "PyBOP provides various plotting utilities to visualise the results of the optimisation." ] }, { "cell_type": "markdown", "metadata": { "id": "-cWCOiqR04qC" }, "source": [ "### Comparing System Response\n", "\n", "We can quickly plot the system's response using the estimated parameters compared to the initial parameters:\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "ZVfozY0A04qC" }, "outputs": [ { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "for optim, x in zip(optims, xs):\n", " pybop.plot.quick(\n", " optim.cost.problem, problem_inputs=x, title=optim.cost.problem.model.name\n", " )" ] }, { "cell_type": "markdown", "metadata": { "id": "ntIvAJmA04qD" }, "source": [ "### Cost Landscape\n", "\n", "Finally, we can visualise the cost landscape and the path taken by the optimiser:\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 467 }, "id": "tJUJ80Ve04qD", "outputId": "855fbaa2-1e09-4935-eb1a-8caf7f99eb75" }, "outputs": [ { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "for optim in optims:\n", " pybop.plot.surface(optim, title=optim.cost.problem.model.name)" ] } ], "metadata": { "colab": { "provenance": [] }, "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": 0 }