{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## A NMC/Gr parameterisation example using PyBOP\n",
    "\n",
    "This notebook introduces a synthetic re-parameterisation of the single-particle model with corrupted observations. To start, we import the PyBOP package for parameterisation and the PyBaMM package to generate the initial synethic data,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%pip install git+https://github.com/pybop-team/PyBOP.git@develop -q \n",
    "%pip install pybamm -q"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we import the added packages plus any additional dependencies,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pybop\n",
    "import pybamm\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generate Synthetic Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We need to generate the synthetic data required for later reparameterisation. To do this we will run the PyBaMM forward model and store the generated data. This will be integrated into PyBOP in a future release for fast synthetic generation. For now, we define the PyBaMM model with a default parameter set,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "synthetic_model = pybamm.lithium_ion.SPM()\n",
    "params = synthetic_model.default_parameter_values"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now modify individual parameters with the bespoke values and run the simulation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "params.update(\n",
    "    {\n",
    "        \"Negative electrode active material volume fraction\": 0.52,\n",
    "        \"Positive electrode active material volume fraction\": 0.63,\n",
    "    }\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define the experiment and run the forward model to capture the synthetic data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "experiment = pybamm.Experiment(\n",
    "    [\n",
    "        (\n",
    "            \"Discharge at 2C for 5 minutes (1 second period)\",\n",
    "            \"Rest for 2 minutes (1 second period)\",\n",
    "            \"Charge at 1C for 5 minutes (1 second period)\",\n",
    "            \"Rest for 2 minutes (1 second period)\",\n",
    "        ),\n",
    "    ]\n",
    "    * 2\n",
    ")\n",
    "sim = pybamm.Simulation(synthetic_model, experiment=experiment, parameter_values=params)\n",
    "synthetic_sol = sim.solve()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the synthetic data,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sim.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, let's corrupt the synthetic data with 5mV of gaussian noise centered around zero,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corrupt_V = synthetic_sol[\"Terminal voltage [V]\"].data\n",
    "corrupt_V += np.random.normal(0,0.005,len(corrupt_V))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Identify the Parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, to blind fit the synthetic parameters we need to define the observation variables as well as update the forward model to be of PyBOP type (This composes PyBaMM's model class). For the observed voltage variable, we used the newly corrupted voltage array, "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = pybop.lithium_ion.SPM()\n",
    "observations = [\n",
    "            pybop.Observed(\"Time [s]\", synthetic_sol[\"Time [s]\"].data),\n",
    "            pybop.Observed(\"Current function [A]\", synthetic_sol[\"Current [A]\"].data),\n",
    "            pybop.Observed(\"Voltage [V]\", corrupt_V),\n",
    "        ]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we define the targetted forward model parameters for estimation. Furthermore, PyBOP provides functionality to define a prior for the parameters. The initial parameters values used in the estimiation will be randomly drawn from the prior distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fit_params = [\n",
    "    pybop.Parameter(\n",
    "        \"Negative electrode active material volume fraction\",\n",
    "        prior=pybop.Gaussian(0.5, 0.02),\n",
    "        bounds=[0.375, 0.625],\n",
    "    ),\n",
    "    pybop.Parameter(\n",
    "        \"Positive electrode active material volume fraction\",\n",
    "        prior=pybop.Gaussian(0.65, 0.02),\n",
    "        bounds=[0.525, 0.75],\n",
    "    ),\n",
    "]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now construct PyBOP's parameterisation class. This class provides the parameterisation methods needed to fit the forward model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "parameterisation = pybop.Parameterisation(\n",
    "    model, observations=observations, fit_parameters=fit_params\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we run the estimation algorithm. For this example, we use a root-mean square cost function with the BOBYQA algorithm implemented in NLOpt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results, last_optim, num_evals = parameterisation.rmse(\n",
    "    signal=\"Voltage [V]\", method=\"nlopt\"\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plotting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First, run SPM forward model with the estimated parameters,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "params.update(\n",
    "        {\"Negative electrode active material volume fraction\": results[0], \n",
    "        \"Positive electrode active material volume fraction\": results[1]}\n",
    "        )\n",
    "optsol = sim.solve()[\"Terminal voltage [V]\"].data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, we plot the estimated forward model against the corrupted synthetic observation,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(corrupt_V, label='Groundtruth')\n",
    "plt.plot(optsol, label='Estimated')\n",
    "plt.xlabel('Time (s)')\n",
    "plt.ylabel('Voltage (V)')\n",
    "plt.legend()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "python",
   "version": "3.10.12"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}