{
 "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": [
    "import pyprobe\n",
    "import polars as pl\n",
    "from pyprobe.analysis import degradation_mode_analysis as dma\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "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",
    "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(base_dataframe=pl.DataFrame({'Capacity [Ah]': capacity, 'Voltage [V]': voltage}), info= {})\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 = pyprobe.Plot()\n",
    "fig.add_line(fitted_curve, x='SOC', y='Input Voltage [V]',label = \"Input\")\n",
    "fig.add_line(fitted_curve, x='SOC', y='Fitted Voltage [V]', color='red', label = 'Fit', dash='dash')\n",
    "fig.show_image()\n",
    "# fig.show() # This will show the plot interactively, it is commented out for the sake of the documentation"
   ]
  },
  {
   "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={\n",
    "        'bounds': [(0.8, 1), (0.2, 0.3), (0, 0.1), (0.86, 1)]\n",
    "    }\n",
    ")\n",
    "print(stoichiometry_limits.data)\n",
    "fig = pyprobe.Plot()\n",
    "fig.add_line(fitted_curve, x='SOC', y='Input dSOCdV [1/V]',label = \"Input\")\n",
    "fig.add_line(fitted_curve, x='SOC', y='Fitted dSOCdV [1/V]', color='red', label = 'Fit', dash='dash')\n",
    "fig.show_image()\n",
    "# fig.show() # This will show the plot interactively, it is commented out for the sake of the documentation"
   ]
  },
  {
   "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": {
  "kernelspec": {
   "display_name": "pyprobe-dev",
   "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}