"""Fit measured EIS to the lumped Rs + (Rsh||Cj) + (Rion + Cion) model using scipy."""
from __future__ import annotations
import math
from typing import Any
import numpy as np
from scipy.optimize import least_squares
from psc_sim.eis_lumped import eis_impedance
from psc_sim.parameters import CellParameters
[docs]
def fit_lumped_eis(
freqs_hz: np.ndarray,
zre_meas: np.ndarray,
zim_meas: np.ndarray,
p0: CellParameters,
*,
max_nfev: int = 500,
) -> tuple[CellParameters, Any]:
"""
Adjust Rs, Rsh, Cj, Rion, Cion to minimize (Zre,Zim) RMSE.
I0,n,Iph,T are held fixed from ``p0``.
"""
mask = np.isfinite(freqs_hz) & (freqs_hz > 0.0) & np.isfinite(zre_meas) & np.isfinite(zim_meas)
f = np.asarray(freqs_hz, dtype=float)[mask]
zr = np.asarray(zre_meas, dtype=float)[mask]
zi = np.asarray(zim_meas, dtype=float)[mask]
if len(f) < 5:
raise ValueError("Need at least 5 finite EIS points")
def pack(x: np.ndarray) -> CellParameters:
Rs, lrsh, lcj, Rion, lcion = (
float(x[0]),
float(x[1]),
float(x[2]),
float(x[3]),
float(x[4]),
)
return CellParameters(
I0=p0.I0,
n=p0.n,
Iph=p0.Iph,
T=p0.T,
Rs=Rs,
Rsh=10.0**lrsh,
Cj=10.0**lcj,
Rion=Rion,
Cion=10.0**lcion,
)
x0 = np.array(
[
p0.Rs,
math.log10(max(p0.Rsh, 1.0)),
math.log10(max(p0.Cj, 1e-15)),
p0.Rion,
math.log10(max(p0.Cion, 1e-15)),
],
dtype=float,
)
lb = np.array([0.0, 1.0, -14.0, 1e-3, -14.0], dtype=float)
ub = np.array([100.0, 8.0, -2.0, 2e6, -2.0], dtype=float)
def residuals(x: np.ndarray) -> np.ndarray:
pm = pack(x)
zre_m, zim_m = eis_impedance(pm, f)
err = np.concatenate([zre_m - zr, zim_m - zi])
s = 1.0 / max(float(np.max(np.abs(zr))), float(np.max(np.abs(zi))), 1.0)
return err * s
out = least_squares(
residuals,
x0,
bounds=(lb, ub),
loss="linear",
max_nfev=max_nfev,
verbose=0,
)
return pack(out.x), out