"""Implicit 1-diode I–V solver (generation branch, Voc bracketing, IV sweep)."""
from __future__ import annotations
import math
from collections.abc import Callable
import numpy as np
from psc_sim.iv import residual_iv_implicit
from psc_sim.parameters import CellParameters
from psc_sim.simulation.types import IVResult, SweepSpec
IV_TRUNCATION_HINT = (
"高電圧の一部では、DC の 1 ダイオード＋Rs＋Rsh の陰式に厳密な解がなく I–V を打ち切っています"
"（回路図の C・Rion はこの I–V 計算に含まれません）。"
"Rsh を大きくする・I₀ を下げる・電圧レンジを狭めると改善することがあります。"
)
def _bracket_root(
f: Callable[[float], float],
lo: float,
hi: float,
*,
expand_hi: Callable[[float], float] | None = None,
expand_lo: Callable[[float], float] | None = None,
max_expand: int = 16,
max_bisect: int = 80,
tol_f: float = 1e-12,
tol_x: float = 1e-13,
) -> float | None:
"""Bisection on [lo, hi] after optional bracket expansion. Returns None if no sign change."""
if expand_lo is not None:
for _ in range(max_expand):
if f(lo) >= 0.0:
break
lo = expand_lo(lo)
else:
return None
flo = f(lo)
fhi = f(hi)
if expand_hi is not None:
for _ in range(max_expand):
if flo == 0.0:
return lo
if flo * fhi <= 0.0:
break
hi = expand_hi(hi)
fhi = f(hi)
else:
return None
elif flo * fhi > 0.0:
return None
flo = f(lo)
for _ in range(max_bisect):
mid = 0.5 * (lo + hi)
fm = f(mid)
if abs(fm) < tol_f or abs(hi - lo) < tol_x:
return mid
if flo * fm <= 0.0:
hi = mid
else:
lo = mid
flo = fm
return 0.5 * (lo + hi)
[docs]
def terminal_voc(p: CellParameters, rs: float | None = None, *, v_hi: float = 2.0) -> float | None:
"""Return the I=0 terminal voltage for the diode + shunt branch."""
rs_use = p.Rs if rs is None else rs
def f0(v: float) -> float:
return residual_iv_implicit(0.0, v, p, rs_use)
result = _bracket_root(
f0,
0.0,
max(v_hi, 0.1),
expand_hi=lambda h: h * 1.5,
tol_f=1e-13,
tol_x=1e-10,
)
return result
[docs]
def solve_generation_current(
V: float,
p: CellParameters,
rs: float,
i_hint: float | None = None,
) -> float:
"""Solve the physical photovoltaic branch between generated current and Voc."""
f_hi = residual_iv_implicit(0.0, V, p, rs)
if abs(f_hi) < 1e-13:
return 0.0
if f_hi > 0.0:
return float("nan")
lo = -max(2.0 * abs(p.Iph), abs(i_hint or 0.0) * 1.4, 1e-9)
def f_i(i: float) -> float:
return residual_iv_implicit(i, V, p, rs)
result = _bracket_root(
f_i,
lo,
0.0,
expand_lo=lambda x: x * 1.8,
max_expand=24,
tol_f=1e-12,
tol_x=1e-13,
)
return float("nan") if result is None else result
[docs]
def solve_current(
V: float,
p: CellParameters,
rs: float | None = None,
*,
i_hint: float | None = None,
) -> float:
"""Single public API for terminal current at bias V (generation branch)."""
rs_use = p.Rs if rs is None else rs
return solve_generation_current(V, p, rs_use, i_hint)
[docs]
def sweep_iv_curve(
p: CellParameters,
v_min: float = -0.3,
v_max: float = 1.25,
step: float = 0.01,
rs: float | None = None,
*,
spec: SweepSpec | None = None,
) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
"""Voltage sweep with Voc truncation beyond the physical open-circuit point."""
if spec is not None:
v_min, v_max, step = spec.v_min, spec.v_max, spec.step
result = sweep_iv_result(p, SweepSpec(v_min, v_max, step), rs)
return result.V, result.I, result.P
[docs]
def sweep_iv_result(
p: CellParameters,
spec: SweepSpec,
rs: float | None = None,
) -> IVResult:
"""Voltage sweep returning IVResult with truncation flag."""
rs_use = p.Rs if rs is None else rs
v_min, v_max, step = spec.v_min, spec.v_max, spec.step
voc = terminal_voc(p, rs_use, v_hi=max(v_max, 0.2)) if rs is None or rs_use == p.Rs else None
V = np.arange(v_min, v_max + step * 0.5, step, dtype=float)
if (
voc is not None
and v_min <= voc <= v_max
and not np.any(np.isclose(V, voc, atol=1e-9, rtol=0.0))
):
V = np.sort(np.append(V, voc))
I = np.empty_like(V)
prev: float | None = None
for k in range(V.size):
v_raw = float(V[k])
vq = round(v_raw, 6)
voc_tol = max(1e-9, abs(step) * 1e-5)
if voc is not None and v_raw > voc + voc_tol:
I[k] = float("nan")
continue
if voc is not None and abs(v_raw - voc) <= voc_tol:
i = 0.0
else:
i = solve_generation_current(vq, p, rs_use, prev)
if math.isfinite(i):
I[k] = i
prev = i
else:
I[k] = float("nan")
P = V * I
return IVResult.from_arrays(V, I, P)
def has_iv_truncation(I: np.ndarray) -> bool:
return bool(np.any(~np.isfinite(I)))