diff --git a/.all-contributorsrc b/.all-contributorsrc index 7f9a01b02..c3fb03f35 100644 --- a/.all-contributorsrc +++ b/.all-contributorsrc @@ -130,6 +130,15 @@ "contributions": [ "example" ] + }, + { + "login": "Dibyendu-IITKGP", + "name": "Dibyendu-IITKGP", + "avatar_url": "https://avatars.githubusercontent.com/u/32595915?v=4", + "profile": "https://github.com/Dibyendu-IITKGP", + "contributions": [ + "example" + ] } ], "contributorsPerLine": 7, diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index f063b3d4d..8d6b58c4b 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -4,7 +4,7 @@ ci: repos: - repo: https://github.com/astral-sh/ruff-pre-commit - rev: "v0.6.4" + rev: "v0.6.5" hooks: - id: ruff args: [--fix, --show-fixes] diff --git a/CHANGELOG.md b/CHANGELOG.md index 2a8834327..2e36303e2 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,6 +4,8 @@ ## Bug Fixes +- [#505](https://github.com/pybop-team/PyBOP/pull/505) - Bug fixes for `LogPosterior` with transformed `GaussianLogLikelihood` likelihood. + ## Breaking Changes # [v24.9.1](https://github.com/pybop-team/PyBOP/tree/v24.9.0) - 2024-09-16 diff --git a/README.md b/README.md index c57d06423..e410640c8 100644 --- a/README.md +++ b/README.md @@ -140,6 +140,7 @@ Thanks goes to these wonderful people ([emoji key](https://allcontributors.org/d Muhammed Nedim Sogut
Muhammed Nedim Sogut

💻 MarkBlyth
MarkBlyth

💻 f-g-r-i-m-m
f-g-r-i-m-m

💡 + Dibyendu-IITKGP
Dibyendu-IITKGP

💡 diff --git a/examples/notebooks/LG_M50_ECM/data/LGM50_5Ah_OCV.mat b/examples/data/LG_M50_ECM/data/LGM50_5Ah_OCV.mat similarity index 100% rename from examples/notebooks/LG_M50_ECM/data/LGM50_5Ah_OCV.mat rename to examples/data/LG_M50_ECM/data/LGM50_5Ah_OCV.mat diff --git a/examples/notebooks/LG_M50_ECM/data/LGM50_5Ah_Pulse.mat b/examples/data/LG_M50_ECM/data/LGM50_5Ah_Pulse.mat similarity index 100% rename from examples/notebooks/LG_M50_ECM/data/LGM50_5Ah_Pulse.mat rename to examples/data/LG_M50_ECM/data/LGM50_5Ah_Pulse.mat diff --git a/examples/notebooks/LG_M50_ECM/data/LGM50_5Ah_RateTest.mat b/examples/data/LG_M50_ECM/data/LGM50_5Ah_RateTest.mat similarity index 100% rename from examples/notebooks/LG_M50_ECM/data/LGM50_5Ah_RateTest.mat rename to examples/data/LG_M50_ECM/data/LGM50_5Ah_RateTest.mat diff --git a/examples/data/Samsung_INR21700/multipulse_hppc.xlsx b/examples/data/Samsung_INR21700/multipulse_hppc.xlsx new file mode 100644 index 000000000..aa7fd79c7 Binary files /dev/null and b/examples/data/Samsung_INR21700/multipulse_hppc.xlsx differ diff --git a/examples/data/Samsung_INR21700/sample_drive_cycle.xlsx b/examples/data/Samsung_INR21700/sample_drive_cycle.xlsx new file mode 100644 index 000000000..5e3555c28 Binary files /dev/null and b/examples/data/Samsung_INR21700/sample_drive_cycle.xlsx differ diff --git a/examples/data/Samsung_INR21700/sample_hppc_pulse.xlsx b/examples/data/Samsung_INR21700/sample_hppc_pulse.xlsx new file mode 100644 index 000000000..d5118b61f Binary files /dev/null and b/examples/data/Samsung_INR21700/sample_hppc_pulse.xlsx differ diff --git a/examples/notebooks/LG_M50_ECM/1-single-pulse-circuit-model.ipynb b/examples/notebooks/1-single-pulse-circuit-model.ipynb similarity index 76% rename from examples/notebooks/LG_M50_ECM/1-single-pulse-circuit-model.ipynb rename to examples/notebooks/1-single-pulse-circuit-model.ipynb index 91b44943f..2478ce755 100644 --- a/examples/notebooks/LG_M50_ECM/1-single-pulse-circuit-model.ipynb +++ b/examples/notebooks/1-single-pulse-circuit-model.ipynb @@ -25,20 +25,14 @@ "name": "stdout", "output_type": "stream", "text": [ - "/Users/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", "zsh:1: no matches found: pybop[plot]\r\n" ] }, @@ -46,14 +40,15 @@ "name": "stdout", "output_type": "stream", "text": [ - "Note: you may need to restart the kernel to use updated packages.\n" + "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" ] }, { "name": "stdout", "output_type": "stream", "text": [ - "/Users/engs2510/Documents/Git/PyBOP/.nox/notebooks-overwrite/bin/python3: No module named pip\r\n" + "\r\n" ] }, { @@ -140,9 +135,17 @@ "metadata": {}, "outputs": [], "source": [ - "ocp = loadmat(\"data/LGM50_5Ah_OCV.mat\", simplify_cells=True, mat_dtype=False)\n", - "pulse_data = loadmat(\"data/LGM50_5Ah_Pulse.mat\", simplify_cells=True, mat_dtype=False)\n", - "rate_data = loadmat(\"data/LGM50_5Ah_RateTest.mat\", simplify_cells=True, mat_dtype=False)" + "ocp = loadmat(\n", + " \"../data/LG_M50_ECM/data/LGM50_5Ah_OCV.mat\", simplify_cells=True, mat_dtype=False\n", + ")\n", + "pulse_data = loadmat(\n", + " \"../data/LG_M50_ECM/data/LGM50_5Ah_Pulse.mat\", simplify_cells=True, mat_dtype=False\n", + ")\n", + "rate_data = loadmat(\n", + " \"../data/LG_M50_ECM/data/LGM50_5Ah_RateTest.mat\",\n", + " simplify_cells=True,\n", + " mat_dtype=False,\n", + ")" ] }, { @@ -193,7 +196,7 @@ " require.undef(\"plotly\");\n", " requirejs.config({\n", " paths: {\n", - " 'plotly': ['https://cdn.plot.ly/plotly-2.34.0.min']\n", + " 'plotly': ['https://cdn.plot.ly/plotly-2.35.2.min']\n", " }\n", " });\n", " require(['plotly'], function(Plotly) {\n", @@ -210,9 +213,9 @@ { "data": { "text/html": [ - "
\n", + " " + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "text/html": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "pybop.quick_plot(problem, problem_inputs=x, title=\"Optimised Comparison\");" + ] + }, + { + "cell_type": "markdown", + "id": "16", + "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": "17", + "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": "18", + "metadata": {}, + "source": [ + "## Using the estimated parameter for thermal predictions\n", + "With the estimated RC parameters, the temperature distribution for a given drive cycle can be calculated using the identified Thevenin model. We now we use a `.xlsx` file containing time-series current data as a pybamm experiment. A sample file is used here, but user's may choose to upload customized drive cycle." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "19", + "metadata": {}, + "outputs": [], + "source": [ + "file_loc = r\"../data/Samsung_INR21700/sample_drive_cycle.xlsx\"\n", + "df = pd.read_excel(file_loc, sheet_name=\"Sheet3\", index_col=None, na_values=[\"NA\"])\n", + "\n", + "# Remove duplicate rows, keeping the first occurrence\n", + "df = df.drop_duplicates(subset=[\"Time\"], keep=\"first\")\n", + "\n", + "# Create the pybamm experiment\n", + "experiment = pybamm.Experiment([pybamm.step.current(df.to_numpy())])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "20", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "[Step([[ 0 3]\n", + " [ 1 3]\n", + " [ 2 3]\n", + " ...\n", + " [3598 3]\n", + " [3599 3]\n", + " [3600 3]], duration=3600, period=1, direction=Discharge)]" + ] + }, + "execution_count": null, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "experiment.steps" + ] + }, + { + "cell_type": "markdown", + "id": "21", + "metadata": {}, + "source": [ + "Update the estimated RC values. These values will be used to calculate heat generation and corresponding temperature distribution in the thermal submodel. Given `model.predict` is a light wrapper on the `PyBaMM.Simulation` class, we interact with it the same way. Visualisation of voltage response and cell temperature is plotted below using the PyBaMM solution." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "22", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "2024-09-17 13:35:23.641 - [WARNING] callbacks.on_experiment_infeasible_time(240): \n", + "\n", + "\tExperiment is infeasible: default duration (3600 seconds) was reached during 'Step([[ 0 3]\n", + " [ 1 3]\n", + " [ 2 3]\n", + " ...\n", + " [3598 3]\n", + " [3599 3]\n", + " [3600 3]], duration=3600, period=1, direction=Discharge)'. The returned solution only contains up to step 1 of cycle 1. Please specify a duration in the step instructions.\n" + ] + } + ], + "source": [ + "sol = model.predict(\n", + " inputs=x,\n", + " experiment=experiment,\n", + " parameter_set=parameter_set,\n", + " initial_state={\"Initial SoC\": 0.95},\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "23", + "metadata": {}, + "source": [ + "## Conclusion\n", + "\n", + "This notebook illustrates how to extract EC parameters from an HPPC pulse using XNES in PyBOP, providing insights into the optimisation process through various visualisations. The estimated parameters are then used to run a thermal submodel." + ] + } + ], + "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 +} diff --git a/examples/notebooks/equivalent_circuit_identification_multipulse.ipynb b/examples/notebooks/equivalent_circuit_identification_multipulse.ipynb new file mode 100644 index 000000000..1cd5413e3 --- /dev/null +++ b/examples/notebooks/equivalent_circuit_identification_multipulse.ipynb @@ -0,0 +1,595 @@ +{ + "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 +} diff --git a/examples/scripts/spm_MAP.py b/examples/scripts/spm_MAP.py index dd8216674..07b3eb970 100644 --- a/examples/scripts/spm_MAP.py +++ b/examples/scripts/spm_MAP.py @@ -25,6 +25,7 @@ bounds=[0.3, 0.8], initial_value=0.653, true_value=parameter_set["Negative electrode active material volume fraction"], + transformation=pybop.LogTransformation(), ), pybop.Parameter( "Positive electrode active material volume fraction", @@ -32,6 +33,7 @@ bounds=[0.4, 0.7], initial_value=0.657, true_value=parameter_set["Positive electrode active material volume fraction"], + transformation=pybop.LogTransformation(), ), ) @@ -44,7 +46,7 @@ ), ] ) -values = model.predict(initial_state={"Initial SoC": 0.7}, experiment=experiment) +values = model.predict(initial_state={"Initial SoC": 0.5}, experiment=experiment) corrupt_values = values["Voltage [V]"].data + np.random.normal( 0, sigma, len(values["Voltage [V]"].data) ) @@ -60,7 +62,7 @@ # Generate problem, cost function, and optimisation class problem = pybop.FittingProblem(model, parameters, dataset) -cost = pybop.LogPosterior(pybop.GaussianLogLikelihoodKnownSigma(problem, sigma0=sigma)) +cost = pybop.LogPosterior(pybop.GaussianLogLikelihood(problem)) optim = pybop.IRPropMin( cost, sigma0=0.05, diff --git a/noxfile.py b/noxfile.py index d270c0fda..2032476f2 100644 --- a/noxfile.py +++ b/noxfile.py @@ -65,6 +65,7 @@ def examples(session): @nox.session def notebooks(session): """Run the Jupyter notebooks.""" + session.install("openpyxl", "ipywidgets") session.install("-e", ".[all,dev]", silent=False) if PYBOP_SCHEDULED: session.run("pip", "install", f"pybamm=={PYBAMM_VERSION}", silent=False) @@ -80,6 +81,7 @@ def notebooks(session): @nox.session(name="notebooks-overwrite") def notebooks_overwrite(session): """Run the Jupyter notebooks.""" + session.install("openpyxl", "ipywidgets") session.install("-e", ".[all,dev]", silent=False) if PYBOP_SCHEDULED: session.run("pip", "install", f"pybamm=={PYBAMM_VERSION}", silent=False) diff --git a/pybop/costs/_likelihoods.py b/pybop/costs/_likelihoods.py index 4a17a3e05..7b4fe3a49 100644 --- a/pybop/costs/_likelihoods.py +++ b/pybop/costs/_likelihoods.py @@ -141,7 +141,7 @@ def _add_single_sigma(self, index, value): Parameter( f"Sigma for output {index+1}", initial_value=value, - prior=Uniform(0.5 * value, 1.5 * value), + prior=Uniform(1e-3 * value, 1e3 * value), bounds=[1e-8, 3 * value], ) ) @@ -235,10 +235,12 @@ def __init__( super().__init__(problem=log_likelihood.problem) self.gradient_step = gradient_step - # Store the likelihood and prior + # Store the likelihood, prior, update parameters and transformation + self.join_parameters(log_likelihood.parameters) self._log_likelihood = log_likelihood - self._parameters = self._log_likelihood.parameters - self._has_separable_problem = self._log_likelihood.has_separable_problem + + for attr in ["transformation", "_has_separable_problem"]: + setattr(self, attr, getattr(log_likelihood, attr)) if log_prior is None: self._prior = JointLogPrior(*self._parameters.priors()) diff --git a/tests/integration/test_optimisation_options.py b/tests/integration/test_optimisation_options.py index 258bdc17b..d59523fcb 100644 --- a/tests/integration/test_optimisation_options.py +++ b/tests/integration/test_optimisation_options.py @@ -13,8 +13,10 @@ class TestOptimisation: @pytest.fixture(autouse=True) def setup(self): - self.ground_truth = np.asarray([0.55, 0.55]) + np.random.normal( - loc=0.0, scale=0.05, size=2 + self.ground_truth = np.clip( + np.asarray([0.55, 0.55]) + np.random.normal(loc=0.0, scale=0.05, size=2), + a_min=0.4, + a_max=0.75, ) @pytest.fixture diff --git a/tests/integration/test_transformation.py b/tests/integration/test_transformation.py index 5c616a5b7..26ce2a7ac 100644 --- a/tests/integration/test_transformation.py +++ b/tests/integration/test_transformation.py @@ -68,6 +68,7 @@ def noise(self, sigma, values): pybop.SumofPower, pybop.Minkowski, pybop.LogPosterior, + pybop.LogPosterior, # Second for GaussianLogLikelihood ] ) def cost_cls(self, request): @@ -90,13 +91,19 @@ def cost(self, model, parameters, init_soc, cost_cls): problem = pybop.FittingProblem(model, parameters, dataset) # Construct the cost + first_map = True if cost_cls is pybop.GaussianLogLikelihoodKnownSigma: return cost_cls(problem, sigma0=self.sigma0) elif cost_cls is pybop.GaussianLogLikelihood: return cost_cls(problem) + elif cost_cls is pybop.LogPosterior and first_map: + first_map = False + return cost_cls(log_likelihood=pybop.GaussianLogLikelihood(problem)) elif cost_cls is pybop.LogPosterior: return cost_cls( - pybop.GaussianLogLikelihoodKnownSigma(problem, sigma0=self.sigma0) + log_likelihood=pybop.GaussianLogLikelihoodKnownSigma( + problem, sigma0=self.sigma0 + ) ) else: return cost_cls(problem) @@ -114,7 +121,7 @@ def test_thevenin_transformation(self, optimiser, cost): optim = optimiser( cost=cost, sigma0=[0.03, 0.03, 1e-3] - if isinstance(cost, pybop.GaussianLogLikelihood) + if isinstance(cost, (pybop.GaussianLogLikelihood, pybop.LogPosterior)) else [0.03, 0.03], max_unchanged_iterations=35, absolute_tolerance=1e-6, @@ -125,7 +132,7 @@ def test_thevenin_transformation(self, optimiser, cost): x, final_cost = optim.run() # Add sigma0 to ground truth for GaussianLogLikelihood - if isinstance(optim.cost, pybop.GaussianLogLikelihood): + if isinstance(optim.cost, (pybop.GaussianLogLikelihood, pybop.LogPosterior)): self.ground_truth = np.concatenate( (self.ground_truth, np.asarray([self.sigma0])) ) diff --git a/tests/unit/test_posterior.py b/tests/unit/test_posterior.py index b43e21ecd..8d45ddc28 100644 --- a/tests/unit/test_posterior.py +++ b/tests/unit/test_posterior.py @@ -68,15 +68,13 @@ def prior(self): def test_log_posterior_construction(self, likelihood, prior): # Test log posterior construction posterior = pybop.LogPosterior(likelihood, prior) + keys = likelihood.parameters.keys() assert posterior._log_likelihood == likelihood assert posterior._prior == prior - - # Test log posterior construction without parameters - likelihood.parameters.priors = None - - with pytest.raises(TypeError, match="'NoneType' object is not callable"): - pybop.LogPosterior(likelihood, log_prior=None) + assert posterior.parameters[keys[0]] == likelihood.parameters[keys[0]] + assert posterior.has_separable_problem == likelihood.has_separable_problem + assert posterior.transformation == likelihood.transformation @pytest.mark.unit def test_log_posterior_construction_no_prior(self, likelihood):