{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fitting electrolyte conductivity (Landesfeind 2019 form)\n",
    "\n",
    "[Landesfeind & Gasteiger (2019)](https://doi.org/10.1149/2.0571912jes) report measurements of the four binary-electrolyte transport properties — ionic conductivity $\\kappa$, salt diffusion coefficient $D_e$, thermodynamic factor $\\chi$ and cation transference number $t_+^0$ — for three solvent systems (`EC:DMC (1:1)`, `EC:EMC (3:7)`, `EMC:FEC (19:1)`) over a wide $(c_e, T)$ range, and provide a **closed-form fit** for each property. For the ionic conductivity (in mS/cm, with $c$ in mol/L and $T$ in K) the closed form is\n",
    "\n",
    "$$\\kappa(c, T) = \\underbrace{p_1\\,(1 + (T - p_2))}_{A(T)}\\;c\\;\n",
    "\\frac{\\overbrace{1 + p_3\\,\\sqrt{c} + p_4\\,(1 + p_5\\,e^{1000/T})\\,c}^{B(c,T)}}\n",
    "{\\underbrace{1 + p_6\\,e^{1000/T}\\,c^{4}}_{C(c,T)}}\\,,$$\n",
    "\n",
    "with six fitted coefficients $p_1, \\dots, p_6$ per solvent system. (The diffusivity, TDF and transference number have analogous closed forms with their own coefficients — 4, 9 and 9 respectively.)\n",
    "\n",
    "`iwp.direct_entries.landesfeind_electrolyte(c_e, system)` packages all four functions plus the published coefficients for the chosen system as a single `DirectEntry`. Critically, the coefficients are stored as `pybamm.Parameter`s rather than hard-coded numbers, which means we can **override any of them with `iwp.Parameter` in a pipeline** to refit them to our own data — without reimplementing the closed form. This notebook demonstrates that workflow on the conductivity coefficients.\n",
    "\n",
    "We work at a single isotherm ($T = 293.15$ K, the closest measurement temperature in Landesfeind 2019 to $25$ °C). We compare three things:\n",
    "\n",
    "1. **Digitized paper data** at $T = 293.15$ K, plotted against the **published equation** (no fit — the published $p_1\\dots p_6$ are already in the direct entry).\n",
    "2. **Mock data** generated from the published equation at $T = 293.15$ K with added Gaussian noise.\n",
    "3. A **fit to the mock data**, recovering a subset of the published coefficients.\n",
    "\n",
    "At a single $T$ the six coefficients are not all separately identifiable: $p_2$ enters only via the prefactor $A(T) = p_1\\,(1 + (T - p_2))$ and is perfectly correlated with $p_1$ when $T$ is fixed; likewise the temperature-coupled terms $p_5\\,e^{1000/T}$ and $p_6\\,e^{1000/T}$ collapse to constants. We therefore hold $p_2$, $p_5$, $p_6$ at their published values and fit the remaining three ($p_1$, $p_3$, $p_4$)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ionworkspipeline as iwp\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import pybamm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. The Landesfeind direct entry\n",
    "\n",
    "`landesfeind_electrolyte` returns a `DirectEntry` whose `parameters` dict contains the conductivity, diffusivity, TDF and transference number functions, plus the published numeric values for every coefficient those functions reference. We grab the conductivity function — which still has $p_1, \\dots, p_6$ as `pybamm.Parameter`s — so we can evaluate it later with the published coefficients (for reference curves) and inside the fitting objective (with $p_1, p_3, p_4$ replaced by `InputParameter`s the optimizer sets)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "electrolyte_entry = iwp.direct_entries.landesfeind_electrolyte(\n",
    "    c_e=1000, system=\"EC:EMC (3:7)\"\n",
    ")\n",
    "kappa_func = electrolyte_entry.parameters[\"Electrolyte conductivity [S.m-1]\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Digitized paper data at $T = 293.15$ K\n",
    "\n",
    "The CSV `landesfeind_2019_ec_emc_3_7.csv` was sourced from the [CALiSol-23 dataset](https://github.com/Pele0599/CALiSol-23), which digitized the conductivity points from the figures of Landesfeind & Gasteiger 2019. We filter to the $T = 293.15$ K isotherm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "all_paper_data = pd.read_csv(\"landesfeind_2019_ec_emc_3_7.csv\")\n",
    "T_iso = 293.15  # K\n",
    "paper_data = all_paper_data[all_paper_data[\"Temperature [K]\"] == T_iso].reset_index(\n",
    "    drop=True\n",
    ")\n",
    "paper_data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Generate mock data\n",
    "\n",
    "We sample the published equation on a denser concentration grid at the same temperature and add 5% Gaussian noise. This is what a clean experimental campaign at this isotherm might produce."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rng = np.random.default_rng(0)\n",
    "c_mock = np.linspace(100, 3000, 12)  # mol/m^3\n",
    "T_mock = np.full_like(c_mock, T_iso)\n",
    "\n",
    "c_vec = pybamm.Vector(c_mock[:, np.newaxis])\n",
    "T_vec = pybamm.Vector(T_mock[:, np.newaxis])\n",
    "kappa_published = (\n",
    "    iwp.ParameterValues(electrolyte_entry.parameters)\n",
    "    .process_symbol(kappa_func(c_vec, T_vec))\n",
    "    .evaluate()\n",
    "    .flatten()\n",
    ")\n",
    "kappa_mock = kappa_published * (1 + 0.05 * rng.standard_normal(c_mock.size))\n",
    "mock_data = pd.DataFrame(\n",
    "    {\n",
    "        \"Electrolyte concentration [mol.m-3]\": c_mock,\n",
    "        \"Temperature [K]\": T_mock,\n",
    "        \"Electrolyte conductivity [S.m-1]\": kappa_mock,\n",
    "    }\n",
    ")\n",
    "mock_data.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. A custom fitting objective\n",
    "\n",
    "Because we are fitting an algebraic relationship rather than a PyBaMM simulation, we subclass `FittingObjective` and skip the simulation in `build`. The conductivity function is taken straight from `all_parameter_values` (which the pipeline populates from the direct entry) — we don't reimplement the equation. We bake the data points $(c_e, T)$ into the symbol as `pybamm.Vector`s once in `build`, so each optimizer step is a single vectorized `evaluate` rather than a Python loop. The fitted coefficients $p_1, p_3, p_4$ remain `InputParameter`s that the optimizer feeds in via the `inputs` dict."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class LandesfeindConductivityObjective(iwp.objectives.FittingObjective):\n",
    "    def process_data(self):\n",
    "        data = self.data[\"data\"]\n",
    "        self._processed_data = {\n",
    "            \"Electrolyte conductivity [S.m-1]\": data[\n",
    "                \"Electrolyte conductivity [S.m-1]\"\n",
    "            ].values\n",
    "        }\n",
    "\n",
    "    def build(self, all_parameter_values):\n",
    "        data = self.data[\"data\"]\n",
    "        c_e = data[\"Electrolyte concentration [mol.m-3]\"].values[:, np.newaxis]\n",
    "        T = data[\"Temperature [K]\"].values[:, np.newaxis]\n",
    "        kappa = all_parameter_values[\"Electrolyte conductivity [S.m-1]\"]\n",
    "        self._processed_kappa = all_parameter_values.process_symbol(\n",
    "            kappa(pybamm.Vector(c_e), pybamm.Vector(T))\n",
    "        )\n",
    "\n",
    "    def run(self, inputs, full_output=False):\n",
    "        return {\n",
    "            \"Electrolyte conductivity [S.m-1]\": self._processed_kappa.evaluate(\n",
    "                inputs=inputs\n",
    "            ).flatten()\n",
    "        }"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Fit $p_1$, $p_3$, $p_4$ to the mock data\n",
    "\n",
    "The fit parameters use the same names as the coefficients defined inside `landesfeind_electrolyte`, so the pipeline replaces them with `InputParameter`s while leaving $p_2$, $p_5$, $p_6$ at the published values from the direct entry."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fit_parameters = {\n",
    "    name: iwp.Parameter(name, initial_value=v0, bounds=bnds)\n",
    "    for name, v0, bnds in [\n",
    "        (\"Landesfeind electrolyte conductivity p1\", 0.3, (1e-3, 5)),\n",
    "        (\"Landesfeind electrolyte conductivity p3\", -0.5, (-5, 5)),\n",
    "        (\"Landesfeind electrolyte conductivity p4\", 0.5, (-2, 2)),\n",
    "    ]\n",
    "}\n",
    "\n",
    "objective = LandesfeindConductivityObjective(mock_data)\n",
    "datafit = iwp.DataFit(objective, parameters=fit_parameters)\n",
    "\n",
    "pipeline = iwp.Pipeline(\n",
    "    {\n",
    "        \"electrolyte\": electrolyte_entry,\n",
    "        \"fit\": datafit,\n",
    "    }\n",
    ")\n",
    "fitted_values = pipeline.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compare fitted vs published values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "comparison = pd.DataFrame(\n",
    "    {\n",
    "        \"published\": [electrolyte_entry.parameters[name] for name in fit_parameters],\n",
    "        \"fit (mock data)\": [fitted_values[name] for name in fit_parameters],\n",
    "    },\n",
    "    index=list(fit_parameters),\n",
    ")\n",
    "comparison"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Plot at $T = 293.15$ K\n",
    "\n",
    "We plot, on the same axes:\n",
    "\n",
    "- The 6 digitized paper data points and the published equation curve (no fit — taken straight from the direct entry).\n",
    "- The 12 mock data points and the curve produced by the fit ($p_1$, $p_3$, $p_4$ from the optimizer; $p_2$, $p_5$, $p_6$ from the direct entry's published values)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c_plot = np.linspace(50, 3000, 200)\n",
    "c_vec = pybamm.Vector(c_plot[:, np.newaxis])\n",
    "T_vec = pybamm.Vector(np.full_like(c_plot, T_iso)[:, np.newaxis])\n",
    "sym = kappa_func(c_vec, T_vec)\n",
    "published_curve = (\n",
    "    iwp.ParameterValues(electrolyte_entry.parameters)\n",
    "    .process_symbol(sym)\n",
    "    .evaluate()\n",
    "    .flatten()\n",
    ")\n",
    "fitted_curve = fitted_values.process_symbol(sym).evaluate().flatten()\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(7, 5))\n",
    "ax.plot(\n",
    "    paper_data[\"Electrolyte concentration [mol.m-3]\"] / 1000,\n",
    "    paper_data[\"Electrolyte conductivity [S.m-1]\"],\n",
    "    \"o\",\n",
    "    color=\"C0\",\n",
    "    label=\"Paper data\",\n",
    ")\n",
    "ax.plot(c_plot / 1000, published_curve, \"-\", color=\"C0\", label=\"Published equation\")\n",
    "ax.plot(\n",
    "    mock_data[\"Electrolyte concentration [mol.m-3]\"] / 1000,\n",
    "    mock_data[\"Electrolyte conductivity [S.m-1]\"],\n",
    "    \"s\",\n",
    "    color=\"C1\",\n",
    "    label=\"Mock data\",\n",
    ")\n",
    "ax.plot(c_plot / 1000, fitted_curve, \"--\", color=\"C1\", label=\"Fit to mock data\")\n",
    "ax.set_xlabel(\"LiPF$_6$ concentration [mol/L]\")\n",
    "ax.set_ylabel(\"Electrolyte conductivity [S/m]\")\n",
    "ax.set_title(f\"Landesfeind 2019 EC:EMC (3:7) at T = {T_iso} K\")\n",
    "ax.legend()\n",
    "fig.tight_layout()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}