{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# OCV Fitting\n",
    "\n",
    "Open-circuit voltage fitting is the key step for conducting Degradation Mode Analysis (DMA). PyProBE has a number of built-in methods for this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture\n",
    "%pip install matplotlib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import polars as pl\n",
    "\n",
    "import pyprobe\n",
    "from pyprobe.analysis import degradation_mode_analysis as dma\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this example, we are going to generate a synthetic aged OCV curve. We will use the half cell OCV fits from Chen 2020 [1]."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def graphite_LGM50_ocp_Chen2020(sto):\n",
    "    \"\"\"Chen2020 graphite ocp fit.\"\"\"\n",
    "    u_eq = (\n",
    "        1.9793 * np.exp(-39.3631 * sto)\n",
    "        + 0.2482\n",
    "        - 0.0909 * np.tanh(29.8538 * (sto - 0.1234))\n",
    "        - 0.04478 * np.tanh(14.9159 * (sto - 0.2769))\n",
    "        - 0.0205 * np.tanh(30.4444 * (sto - 0.6103))\n",
    "    )\n",
    "\n",
    "    return u_eq\n",
    "\n",
    "\n",
    "def nmc_LGM50_ocp_Chen2020(sto):\n",
    "    \"\"\"Chen2020 nmc ocp fit.\"\"\"\n",
    "    u_eq = (\n",
    "        -0.8090 * sto\n",
    "        + 4.4875\n",
    "        - 0.0428 * np.tanh(18.5138 * (sto - 0.5542))\n",
    "        - 17.7326 * np.tanh(15.7890 * (sto - 0.3117))\n",
    "        + 17.5842 * np.tanh(15.9308 * (sto - 0.3120))\n",
    "    )\n",
    "\n",
    "    return u_eq\n",
    "\n",
    "\n",
    "z = np.linspace(0, 1, 1000)  # stoichiometry vector\n",
    "\n",
    "# generate complete ocp curves\n",
    "ocp_pe = nmc_LGM50_ocp_Chen2020(z)\n",
    "ocp_ne = graphite_LGM50_ocp_Chen2020(z)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will now define a set of stoichiometry limits to generate our synthetic OCV curve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_pts = 1000\n",
    "# positive electrode\n",
    "x_pe_lo = 0.85  # stoichiometry at low cell SOC\n",
    "x_pe_hi = 0.27  # stoichiometry at high cell SOC\n",
    "x_pe = np.linspace(x_pe_lo, x_pe_hi, n_pts)  # stoichiometry range\n",
    "\n",
    "# negative electrode\n",
    "x_ne_lo = 0.03  # stoichiometry at low cell SOC\n",
    "x_ne_hi = 0.9  # stoichiometry at high cell SOC\n",
    "x_ne = np.linspace(x_ne_lo, x_ne_hi, n_pts)  # stoichiometry range\n",
    "\n",
    "# full cell voltage and capacity\n",
    "voltage = nmc_LGM50_ocp_Chen2020(x_pe) - graphite_LGM50_ocp_Chen2020(x_ne)\n",
    "capacity = np.linspace(0, 5, n_pts)  # capacity range in Ah\n",
    "\n",
    "plt.figure()\n",
    "plt.plot(capacity, voltage)\n",
    "plt.xlabel(\"Capacity (Ah)\")\n",
    "plt.ylabel(\"Voltage (V)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will now generate a fit to this voltage curve.\n",
    "\n",
    "First, we must create instances of the OCP class to hold our electrode OCP information:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "neg_ocp = dma.OCP(graphite_LGM50_ocp_Chen2020)\n",
    "pos_ocp = dma.OCP(nmc_LGM50_ocp_Chen2020)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this example, we have directly declared the OCPs since we already have callable functions for them. Alternatively you can use the `from_data` or `from_expression` to define your OCP from data points or a Sympy expression.\n",
    "\n",
    "Now we can fit to the ocv curve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# put the voltage and capacity data into a Result object (not necessary in normal use)\n",
    "OCV_result = pyprobe.Result(\n",
    "    lf=pl.DataFrame({\"Capacity [Ah]\": capacity, \"Voltage [V]\": voltage}),\n",
    "    info={},\n",
    ")\n",
    "\n",
    "stoichiometry_limits, fitted_curve = dma.run_ocv_curve_fit(\n",
    "    input_data=OCV_result,\n",
    "    ocp_ne=neg_ocp,\n",
    "    ocp_pe=pos_ocp,\n",
    "    fitting_target=\"OCV\",\n",
    "    optimizer=\"minimize\",\n",
    "    optimizer_options={\n",
    "        \"x0\": np.array([0.9, 0.1, 0.1, 0.9]),\n",
    "        \"bounds\": [(0, 1), (0, 1), (0, 1), (0, 1)],\n",
    "    },\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This produces two result objects. The first is the fitted stoichiometry limits:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(stoichiometry_limits.data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And the second is a result object containing the fitted OCP curve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots()\n",
    "fitted_curve.plot(x=\"SOC\", y=\"Input Voltage [V]\", ax=ax, label=\"Input\")\n",
    "fitted_curve.plot(\n",
    "    x=\"SOC\",\n",
    "    y=\"Fitted Voltage [V]\",\n",
    "    ax=ax,\n",
    "    color=\"red\",\n",
    "    label=\"Fit\",\n",
    "    linestyle=\"--\",\n",
    ")\n",
    "ax.set_ylabel(\"Voltage (V)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also fit to differentiated voltage data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "stoichiometry_limits, fitted_curve = dma.run_ocv_curve_fit(\n",
    "    input_data=OCV_result,\n",
    "    ocp_ne=neg_ocp,\n",
    "    ocp_pe=pos_ocp,\n",
    "    fitting_target=\"dQdV\",\n",
    "    optimizer=\"differential_evolution\",\n",
    "    optimizer_options={\"bounds\": [(0.8, 1), (0.2, 0.3), (0, 0.1), (0.86, 1)]},\n",
    ")\n",
    "print(stoichiometry_limits.data)\n",
    "\n",
    "fig, ax = plt.subplots()\n",
    "fitted_curve.plot(x=\"SOC\", y=\"Input dSOCdV [1/V]\", ax=ax, label=\"Input\")\n",
    "fitted_curve.plot(\n",
    "    x=\"SOC\",\n",
    "    y=\"Fitted dSOCdV [1/V]\",\n",
    "    ax=ax,\n",
    "    color=\"red\",\n",
    "    label=\"Fit\",\n",
    "    linestyle=\"--\",\n",
    ")\n",
    "ax.set_ylabel(\"dSOCdV (1/V)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1.  Chen CH, Planella FB, O’Regan K, Gastol D, Widanage WD, Kendrick E. Development of Experimental Techniques for Parameterization of Multi-scale Lithium-ion Battery Models. Journal of The Electrochemical Society. 2020;167(8): 080534. https://doi.org/10.1149/1945-7111/AB9050."
   ]
  }
 ],
 "metadata": {
  "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}