フェーズ3電解液を加えた SPMe への拡張

フェーズ2で作った単一粒子モデル(SPM)は、電解液の中身を一切見ていませんでした。 この章では「電解液に何が起きているか」を DFN の液相物質保存式で記述し、 その結果を2つの電圧補正項として SPM に足します。 こうして得られるモデルが SPMe(Single Particle Model with electrolyte)です。 追加コストはわずかなのに、高 C レートでの電圧降下と容量低下をかなり正しく捉えられるようになります。

3.1 SPM はいつ破れるか — 電解液に勾配が立つ

フェーズ2の SPM は、DFN の5本の支配方程式(フェーズ1)のうち固相拡散と界面反応だけを残し、 電解液については次の2つを仮定していました。

• 電解液塩濃度 $c_e(x,t)$ [mol/m³] はどこでも初期値 $c_{e0}=1000$ mol/m³ のまま変わらない。
• 液相電位 $\phi_e(x,t)$ [V] は $x$ 方向に一様(電解液での電圧損失はゼロ)。

低レート(0.1C 程度)ならこれで十分ですが、レートを上げるとこの2つの仮定が両方とも崩れます。 理由を放電時の Li$^+$ の動きから追いかけてみましょう。 座標はフェーズ1と同じく、負極集電体を $x=0$、正極集電体を $x=L=L_n+L_s+L_p$ にとり、 放電電流密度 $I>0$ [A/m²]、酸化方向の界面電流密度 $j>0$ [A/m²] とします。

放電中、電解液の中では何が起きているか

放電時、負極($0 \le x \le L_n$)では酸化反応($j>0$)によって粒子から Li$^+$ が 電解液に放出され、正極($L_n+L_s \le x \le L$)では還元反応($j<0$)によって Li$^+$ が電解液から吸収されます。つまり電解液から見ると、 負極側が「湧き出し」、正極側が「吸い込み」です。Li$^+$ はセルを 負極 → セパレータ → 正極と横断しなければなりません。

このとき Li$^+$ を運ぶ手段は2つあります。電場による泳動(migration)と、 濃度勾配による拡散(diffusion)です。ここで効いてくるのが 輸率 $t_+^0$ [−](フェーズ1で導入。電解液を流れる電流のうち Li$^+$ の泳動が運ぶ割合)です。 本教材の電解液(LiPF$_6$/EC:EMC)では $t_+^0 = 0.2594$、 つまり泳動が運べるのは電流の約 26% だけで、 残りの約 74% は陰イオン(PF$_6^-$)が逆向きに動くことで運ばれます。 しかし陰イオンは電極反応に参加できないので、負極側に置き去りにされて塩(Li$^+$ と PF$_6^-$ のペア)が溜まり、 正極側では塩が減っていきます。この過不足を解消できるのは拡散だけです。 拡散が濃度差を必要とする以上、電流を流し続ける限り濃度勾配は必ず立ちます。

時間スケールで見る:いつ勾配は「間に合ってしまう」か

濃度勾配が育つのにかかる時間は、拡散時定数 $\tau_e \sim \varepsilon_e L^2 / D_e^{\mathrm{eff}}$ で見積もれます。 ここで $L = L_n + L_s + L_p = 172.8$ μm はセル厚さ、 $\varepsilon_e$ [−] は空隙率(電解液体積分率)、 $D_e^{\mathrm{eff}} = D_e \varepsilon_e^{b}$ は実効電解液拡散係数 [m²/s] (フェーズ1、$b=1.5$ は Bruggeman 指数)です。 本教材のパラメータでは最も遅い負極領域で $D_e^{\mathrm{eff}} \approx 1.77\times10^{-10} \times 0.25^{1.5} \approx 2.2\times10^{-11}$ m²/s なので、 $\tau_e \approx 0.25 \times (172.8\times10^{-6})^2 / (2.2\times10^{-11}) \approx 340$ 秒。 一方、放電にかかる時間は $3600/C$ 秒($C$ は C レート)です。

• 0.1C(10時間放電):$\tau_e \ll$ 放電時間。勾配は立つが浅く、影響は小さい。
• 1C(1時間放電):放電の前半で勾配が定常に達し、数十 mV の損失を生む。
• 3C 以上:定常勾配そのものが深く、正極側の枯渇が現実になる。SPM の電圧予測は 数百 mV 単位で外れ、放電容量も過大評価になる。

まとめると、SPM の誤差は C レートの増加とともに単調に拡大します。 おおまかには電解液のオーム損が $I$ に比例して(このセルでは 1C で約 20 mV)、 濃度起因の損失はそれより速く(勾配が深くなるほど対数関数を通して増幅されて)効いてきます。 3.6節の図3.2でこの傾向を自分の手で確かめます。

3.2 SPMe とは — SPM に何を足し、何をやめないか

SPMe の設計思想はとても割り切っています。DFN の5式(フェーズ1の最終形)のうち、 液相物質保存(式2)だけを PDE として本当に解き、 液相電荷保存(式4)は「解く」のではなく「$x$ で積分して電圧式への補正項に変換」します。 固相電荷保存(式3)は引き続き無視し(電極内で $\phi_s$ 一様)、 反応分布も引き続き一様($j$ が電極内で $x$ によらない)と仮定します。

足すもの

1. 液相物質保存 PDE を解く(SPEC の式2。フェーズ1で導出済み): $$\varepsilon_e \frac{\partial c_e}{\partial t} = \frac{\partial}{\partial x}\!\left(D_e^{\mathrm{eff}}\frac{\partial c_e}{\partial x}\right) + (1-t_+^0)\frac{a j}{F}$$ これで $c_e(x,t)$ の時間発展(3.1節の「塩の偏り」)が手に入ります。
2. 電圧式に電解液起因の補正項を2つ加える(3.5節で導出): 濃度差に由来する $\Delta\phi_{e,\mathrm{conc}}$ と、イオン輸送の抵抗に由来する $\Delta\phi_{e,\mathrm{ohm}}$ です。

やめないもの(SPM から引き継ぐ簡略化)

• 反応分布は依然一様:負極内・正極内で $j$ は $x$ によらない定数。 したがって固相は引き続き「各電極の代表粒子1個ずつ」で足ります。
• $\phi_e(x)$ の厳密解は解かない:液相電荷保存を場の方程式として解くのではなく、 積分した「端子電圧への正味の寄与」だけを取り出します。$\phi_e(x)$ の形そのものは DFN(フェーズ5)まで持ち越しです。
• 固相のオーム損は無視:$\sigma^{\mathrm{eff}}$ [S/m] が十分大きいとして $\phi_s$ は各電極内で一様のまま。

3.3 電解液拡散 PDE の設定(3領域・境界条件・ソース項)

まず解くべき方程式を領域ごとに書き下ろします。使う記号を再掲します (定義と単位はフェーズ1の記号表と同じです)。

3領域での方程式

液相物質保存(式2)の反応項 $aj$ は活物質がある電極内でのみ現れます。セパレータには 活物質がない($aj=0$)ので、次の3本になります。

ここで $j_n$、$j_p$ は負極・正極の界面電流密度 [A/m²] です(SPMe では各電極内で定数。すぐ下で値を求めます)。

領域境界の連続条件

負極/セパレータ境界($x=L_n$)とセパレータ/正極境界($x=L_n+L_s$)では、 次の2つの連続条件を課します。

両端(集電体)のノイマン条件

物理的根拠:$x=0$ は銅箔、$x=L$ はアルミ箔の集電体で、どちらも電子は通しますが イオンも溶媒も通しません。塩の流束はゼロ、つまり $-D_e^{\mathrm{eff}} \partial c_e/\partial x = 0$ です。$D_e^{\mathrm{eff}} \ne 0$ なので 勾配ゼロのノイマン条件(Neumann condition)になります。 塩はセルの外に出られない — この事実が3.4節の「総塩量保存」に直結します。

一様反応の仮定でソース項が定数になる

仮定を認めると、$j_n$、$j_p$ の値は液相電荷保存(式4)の電流バランスだけから決まります。 液相の見かけ電流密度 $i_e(x)$ [A/m²] を使って1ステップずつ確認します。

3.4 定常濃度プロファイルの解析解

定電流を流し続けると、湧き出し・吸い込みと拡散が釣り合って濃度分布は時間変化しなくなります ($\partial c_e/\partial t = 0$)。このときの分布は手で厳密に解けます。 結果は数値解の検証(図3.1の点線)に使えるうえ、3.5節の電圧補正の大きさを見積もる道具にもなります。 表記を軽くするため、この節では領域ごとの実効拡散係数を $D_n \equiv D_e^{\mathrm{eff},-}$、$D_s \equiv D_e^{\mathrm{eff,sep}}$、$D_p \equiv D_e^{\mathrm{eff},+}$ と書き、負極のソース密度を $q_n \equiv (1-t_+^0)I/(F L_n)$ [mol/(m³·s)]、 正極のそれを $q_p \equiv (1-t_+^0)I/(F L_p)$(符号は $-q_p$ で入る)とおきます。

負極:2次関数

セパレータ:1次関数

正極:2次関数(逆向き)

最後の積分定数:総塩量保存

3.5 電圧補正項の導出(この章の核)

$c_e(x,t)$ が手に入ったので、いよいよ「電解液が端子電圧をどれだけ下げるか」を導きます。 出発点は液相電荷保存=修正オームの法則(SPEC 式4、フェーズ1で導出)です。

電圧の分解:どこに $\phi_e$ が現れるか

端子電圧は $V = \phi_s(L) - \phi_s(0)$ です。過電圧の定義 $\eta = \phi_s - \phi_e - U(\theta)$(フェーズ1)を電極内の反応点で使うと、 各電極で $\phi_s = U + \eta + \phi_e$ と書けます。 ここで仮定 3・4 より $U$ と $\eta$ は電極内で一様ですが、$\phi_e(x)$ は位置に依存します。 電極内のどの粒子も同じ $j$ で反応している(仮定3)ので、代表粒子が感じる液相電位としては 電極内の平均値 $\langle \phi_e \rangle_n \equiv \frac{1}{L_n}\int_0^{L_n} \phi_e\,\mathrm{d}x$、 $\langle \phi_e \rangle_p \equiv \frac{1}{L_p}\int_{L_n+L_s}^{L} \phi_e\,\mathrm{d}x$ を使うのが SPMe の流儀です。すると:

SPM では $\phi_e$ 一様の仮定により $\Delta\phi_e = 0$ でした。SPMe はここを埋めます。

修正オームの法則を勾配について解く

(i) 濃度項の積分

(ii) オーム項の準備:$i_e(x)$ は区分的に線形

オーム項 $-i_e/\kappa^{\mathrm{eff}}$ を積分するには $i_e(x)$ の形が要ります。 一様反応なら、これは3.3節の積分を途中で止めるだけで手に入ります。

(iii) オーム項の積分:$1/3$ はどこから来るか

SPMe の電圧式(最終形)

以上を電圧の分解式に戻すと、SPMe の電圧式が完成します。

ここで $\theta_n$、$\theta_p$ [−] は各代表粒子の表面化学量論、 $\eta_{n,p} = \frac{2RT}{F}\mathrm{asinh}\!\left(\frac{j_{n,p}}{2 i_{0,n,p}}\right)$ [V] は Butler–Volmer($\alpha_a=\alpha_c=0.5$)を逆に解いた反応過電圧です(フェーズ2)。 $\kappa^{\mathrm{eff}}$ を $c_{e0}$ で評価する近似(仮定2)により、 $\Delta\phi_{e,\mathrm{ohm}}$ は定電流放電の間じゅう一定値、 時間変化するのは $\Delta\phi_{e,\mathrm{conc}}$ だけです。

数値も見ておきましょう。$\kappa(c_{e0}) = 0.949$ S/m から $\kappa_n^{\mathrm{eff}} = 0.949 \times 0.25^{1.5} = 0.119$ S/m、 $\kappa_s^{\mathrm{eff}} = 0.306$ S/m、$\kappa_p^{\mathrm{eff}} = 0.184$ S/m。 角括弧の中身(電解液の面積比抵抗)は $4.16 \times 10^{-4}$ Ω·m² で、 1C($I = 48.7$ A/m²)なら $\Delta\phi_{e,\mathrm{ohm}} \approx -20$ mV、5C で約 $-101$ mV です。 一方 $\Delta\phi_{e,\mathrm{conc}}$ は 1C の定常状態($c_e(0) \approx 1700$、$c_e(L) \approx 490$ mol/m³、3.4節の解析解)で 約 $-47$ mV。このセルでは濃度項の方が大きいことがわかります。

3.6 実装 — SPMe を動かす

実装は SPM(フェーズ2)に電解液 PDE を1本足すだけです。全体の流れ:

1. 固相:代表粒子2個の球拡散をフェーズ2と同じ FVM で解く(表面流束は $j_n/F$、$j_p/F$ で一定)。
2. 液相:3.3節の PDE を3領域 FVM で解く。領域境界の流束連続は、界面の $D_e^{\mathrm{eff}}$ を隣接セルの調和平均で評価すれば自動的に満たされる(セル中心の 有限体積法の定石。詳細はフェーズ4)。ソース項は 3.3節の定数。
3. 電圧:各時刻で $\theta_{n,p}$ と $c_e(0)$、$c_e(L)$ を取り出し、 3.5節の式で $V = V_{\mathrm{SPM}} + \Delta\phi_{e,\mathrm{conc}} + \Delta\phi_{e,\mathrm{ohm}}$ を組み立てる。

下のインタラクティブ図は、この手順をブラウザ内の JavaScript (時間積分は Crank–Nicolson、フェーズ4で学ぶ方法の1つ)で実行したものです。

実行できる Python:SPMe 完全実装(solve_ivp 版)

仕上げに、この章の内容をすべて含む SPMe を1本の Python にまとめます。 状態ベクトルは $y = [\underbrace{c_s^n}_{15}, \underbrace{c_s^p}_{15}, \underbrace{c_e}_{30}]$ の 60 成分。固相2粒子(各 $N_r = 15$ セルの FVM)と電解液(3領域 $N_x = 30$ セルの FVM)を 連結した ODE 系を solve_ivp(method="BDF") に渡します (拡散方程式系は「硬い(stiff)」ので陰的ソルバが必要 — 理由はフェーズ4で)。 1C と 3C の放電曲線を SPM 電圧(点線)と重ねて描き、各レートの放電容量を表示します。 実行には初回のみ Pyodide のダウンロードで数十秒かかります(計算自体は数秒です)。 コードはそのまま Jupyter にコピーしても動きます(matplotlib 版は末尾のコメント)。

# ============================================================
# SPMe(単一粒子モデル + 電解液)の完全実装
#   固相: 2粒子 x 球対称 FVM (Nr=15) / 電解液: 3領域 FVM (Nx=30)
#   を連結した ODE 系を solve_ivp(method="BDF") で解く。
#   1C と 3C の放電曲線を SPM 電圧(点線)と重ねて描く。
#   Jupyter でもそのまま動く(matplotlib 版は末尾のコメント参照)。
# ============================================================
import numpy as np
from scipy.integrate import solve_ivp

# ---- パラメータ(LG M50, Chen et al. 2020。js/params.js と同一値)----
F, Rg, T = 96485.33212, 8.314462618, 298.15   # C/mol, J/(mol K), K
L_n, L_s, L_p = 85.2e-6, 12.0e-6, 75.6e-6     # 電極厚さ [m]
A_cell, cap_Ah = 0.1027, 5.0                  # 面積 [m^2], 公称容量 [Ah]
R_n, R_p = 5.86e-6, 5.22e-6                   # 粒子半径 [m]
eps_n, eps_s, eps_p = 0.25, 0.47, 0.335       # 空隙率 [-]
epss_n, epss_p, brug = 0.75, 0.665, 1.5       # 活物質分率, Bruggeman
cs_max_n, cs_max_p = 33133.0, 63104.0         # 最大 Li 濃度 [mol/m^3]
Ds_n, Ds_p = 3.3e-14, 4.0e-15                 # 固相拡散係数 [m^2/s]
th_n0, th_p0 = 0.9014, 0.2700                 # 100% SOC の化学量論
ce0, t_plus = 1000.0, 0.2594                  # 初期塩濃度 [mol/m^3], 輸率
k_n, k_p = 6.716e-12, 3.545e-11               # 反応速度定数
a_n, a_p = 3*epss_n/R_n, 3*epss_p/R_p         # 比界面積 [1/m]
I_1C = cap_Ah / A_cell                        # 1C 電流密度 [A/m^2]
V_CUT = 2.5                                   # 下限電圧 [V]

def ocpN(t):   # 黒鉛 OCP [V]
    return (1.9793*np.exp(-39.3631*t) + 0.2482
            - 0.0909*np.tanh(29.8538*(t-0.1234))
            - 0.04478*np.tanh(14.9159*(t-0.2769))
            - 0.0205*np.tanh(30.4444*(t-0.6103)))

def ocpP(t):   # NMC811 OCP [V]
    return (-0.8090*t + 4.4875 - 0.0428*np.tanh(18.5138*(t-0.5542))
            - 17.7326*np.tanh(15.7890*(t-0.3117))
            + 17.5842*np.tanh(15.9308*(t-0.3120)))

# 電解液物性は初期濃度 ce0 で評価して固定(3.5節の近似)
c1 = ce0/1000.0
De    = 8.794e-11*c1**2 - 3.972e-10*c1 + 4.862e-10   # [m^2/s]
kappa = 0.1297*c1**3 - 2.51*c1**1.5 + 3.329*c1       # [S/m]
R_ell = (L_n/(3*kappa*eps_n**brug) + L_s/(kappa*eps_s**brug)
         + L_p/(3*kappa*eps_p**brug))                # 電解液オーム抵抗 [Ohm m^2]

# ---- 固相 FVM(球対称、Nr セル)----
Nr = 15
rf_n = np.linspace(0.0, R_n, Nr+1); rf_p = np.linspace(0.0, R_p, Nr+1)  # セル面
Vsh_n = (rf_n[1:]**3 - rf_n[:-1]**3)/3.0   # 殻体積(共通因子 4π は省略)
Vsh_p = (rf_p[1:]**3 - rf_p[:-1]**3)/3.0
dr_n, dr_p = R_n/Nr, R_p/Nr

def particle_rhs(cs, Ds, rf, Vsh, dr, jmol):
    # Vsh dc/dt = 面フラックス差。表面境界: -Ds dc/dr = jmol(流出が正)
    q = Ds*rf[1:-1]**2*np.diff(cs)/dr      # 内部面の r^2 D dc/dr
    d = np.zeros_like(cs)
    d[:-1] += q; d[1:] -= q
    d[-1] -= rf[-1]**2*jmol
    return d/Vsh

# ---- 電解液 FVM(3領域、Nx=30)----
nn, ns, npos = 15, 3, 12
dx   = np.r_[np.full(nn, L_n/nn), np.full(ns, L_s/ns), np.full(npos, L_p/npos)]
epsx = np.r_[np.full(nn, eps_n),  np.full(ns, eps_s),  np.full(npos, eps_p)]
Deff = De*epsx**brug
G = 1.0/(dx[:-1]/(2*Deff[:-1]) + dx[1:]/(2*Deff[1:]))   # 界面(調和平均)

def electrolyte_rhs(ce, I):
    src = np.r_[np.full(nn,  (1-t_plus)*I/(F*L_n)),     # 3.3節の一様ソース
                np.zeros(ns),
                np.full(npos, -(1-t_plus)*I/(F*L_p))]
    q = G*np.diff(ce)                                    # 両端はノイマン(流束 0)
    d = np.zeros_like(ce)
    d[:-1] += q; d[1:] -= q
    return (d + src*dx)/(epsx*dx)

def rhs(t, y, I):
    return np.r_[particle_rhs(y[:Nr],      Ds_n, rf_n, Vsh_n, dr_n,  I/(a_n*L_n*F)),
                 particle_rhs(y[Nr:2*Nr],  Ds_p, rf_p, Vsh_p, dr_p, -I/(a_p*L_p*F)),
                 electrolyte_rhs(y[2*Nr:], I)]

def voltage(y, I):   # (V_SPM, V_SPMe) を返す
    csn, csp, ce = y[:Nr], y[Nr:2*Nr], y[2*Nr:]
    csn_s = np.clip(csn[-1] - ( I/(a_n*L_n*F))*dr_n/(2*Ds_n), 1.0, cs_max_n-1.0)
    csp_s = np.clip(csp[-1] - (-I/(a_p*L_p*F))*dr_p/(2*Ds_p), 1.0, cs_max_p-1.0)
    i0n = k_n*F*np.sqrt(ce0*csn_s*(cs_max_n-csn_s))      # i0 は ce0 で評価(SPM と共通)
    i0p = k_p*F*np.sqrt(ce0*csp_s*(cs_max_p-csp_s))
    eta_n = 2*Rg*T/F*np.arcsinh( I/(a_n*L_n)/(2*i0n))
    eta_p = 2*Rg*T/F*np.arcsinh(-I/(a_p*L_p)/(2*i0p))
    V_spm = ocpP(csp_s/cs_max_p) - ocpN(csn_s/cs_max_n) + eta_p - eta_n
    dphi_conc = 2*Rg*T/F*(1-t_plus)*np.log(max(ce[-1], 1.0)/max(ce[0], 1.0))
    return V_spm, V_spm + dphi_conc - I*R_ell

def cap_at_cut(cap, V):  # 2.5 V 到達容量 [Ah](線形補間)
    i = np.argmax(V <= V_CUT)
    if V[i] > V_CUT: return cap[-1]
    w = (V[i-1]-V_CUT)/(V[i-1]-V[i])
    return cap[i-1] + w*(cap[i]-cap[i-1])

y0 = np.r_[np.full(Nr, th_n0*cs_max_n), np.full(Nr, th_p0*cs_max_p), np.full(30, ce0)]
traces = []
for crate, col in [(1.0, "#2563eb"), (3.0, "#dc2626")]:
    I = crate*I_1C
    ev = lambda t, y, I: voltage(y, I)[0] - V_CUT        # SPM 電圧(遅い方)で停止
    ev.terminal, ev.direction = True, -1
    t_max = 1.3*3600/crate
    sol = solve_ivp(rhs, (0.0, t_max), y0, args=(I,), method="BDF",
                    t_eval=np.linspace(0.0, t_max, 240), events=ev,
                    rtol=1e-6, atol=1e-2)
    Vs = np.array([voltage(sol.y[:, i], I) for i in range(sol.y.shape[1])])
    V_spm, V_spme = Vs[:, 0], Vs[:, 1]
    cap = I*A_cell*sol.t/3600.0
    # --- SPMe の終了点: 電圧カットオフ、または正極端の電解液枯渇(初期値の1%)---
    # c_e→0 ではその電流を物理的に維持できない。一様反応の SPMe は枯渇を
    # 過大評価しがちなことに注意(DFN との答え合わせはフェーズ5・図5.2)。
    ce_hist = sol.y[-1, :]                     # 正極端の c_e(t)
    idep = int(np.argmax(ce_hist <= 0.01*ce0))
    depleted = bool(ce_hist[idep] <= 0.01*ce0)
    c1_ = cap_at_cut(cap, V_spm)
    c2_ = cap_at_cut(cap, V_spme)
    if depleted:
        c2_ = min(c2_, cap[idep])
    note = "(正極端の電解液枯渇で終了)" if depleted else ""
    print(f"{crate:.0f}C: 放電容量 SPM = {c1_:.3f} Ah / SPMe = {c2_:.3f} Ah{note}")
    iend = max(idep - 1, 0) if depleted else len(sol.t) - 1
    ce_e = sol.y[2*Nr:, iend]
    print(f"    SPMe 終了時の c_e: 負極端 {ce_e[0]:7.1f} / 正極端 {ce_e[-1]:7.1f} mol/m^3")
    m = V_spme >= V_CUT - 1e-9                            # SPMe はカットオフで切って描く
    if depleted:
        m &= np.arange(len(V_spme)) <= idep               # 枯渇以降は描かない
    traces.append({"x": cap[m].tolist(), "y": V_spme[m].tolist(),
                   "name": f"SPMe {crate:.0f}C", "mode": "lines",
                   "line": {"color": col, "width": 2.5}})
    traces.append({"x": cap.tolist(), "y": V_spm.tolist(),
                   "name": f"SPM {crate:.0f}C", "mode": "lines",
                   "line": {"color": col, "dash": "dot", "width": 1.5}})

plot_spec = {"traces": traces,
             "layout": {"title": {"text": "SPMe vs SPM 放電曲線(1C / 3C)"},
                        "xaxis": {"title": {"text": "放電容量 [Ah]"}},
                        "yaxis": {"title": {"text": "端子電圧 V [V]"},
                                  "range": [2.4, 4.3]}}}

# --- Jupyter で matplotlib を使う場合は以下を有効化 ---
# import matplotlib.pyplot as plt
# for tr in plot_spec["traces"]:
#     ls = ":" if tr["line"].get("dash") == "dot" else "-"
#     plt.plot(tr["x"], tr["y"], ls, color=tr["line"]["color"], label=tr["name"])
# plt.xlabel("放電容量 [Ah]"); plt.ylabel("端子電圧 [V]")
# plt.ylim(2.4, 4.3); plt.legend(); plt.grid(alpha=0.3); plt.show()

実行すると、1C では SPM と SPMe の容量差はわずか(約 0.02 Ah)ですが、3C では 0.2 Ah 超に開くはずです。また 3C の出力にある「終了時の $c_e$:正極端」が 負の値になっている点に注目してください。$c_e < 0$ は物理的にあり得ません。 これは $D_e$ と $\kappa$ を初期濃度で固定した線形化モデルが枯渇領域で破綻しているサインで、 コード中では対数の引数を 1 mol/m³ でクリップして電圧の急落として扱っています。 現実のセル(と濃度依存物性を入れた DFN)では、枯渇が近づくと局所の輸送・反応が変化して 完全な負濃度にはなりません — この違いもフェーズ5で見ます。

3.7 SPMe の限界

SPMe は「電解液の物質輸送」は取り込みましたが、SPM から引き継いだ 一様反応の仮定はそのままです。これが次に破れるボトルネックになります。

反応分布の不均一 — 高レートで反応はセパレータ側に集中する

現実の電極では、$j(x)$ は5本の方程式の連立で自己無撞着に決まる分布です。 電流がイオンとして電極の奥(集電体側)まで届くには、電解液のオーム損 $\propto x/\kappa^{\mathrm{eff}}$ を余分に払わなければなりません。 つまりセパレータに近い粒子ほど「電気的に近くて安い」ので、 レートが上がるほど反応はセパレータ側へ偏ります。 さらに高レートでは正極の奥で $c_e$ が枯れて $i_0 \propto c_e^{1/2}$ が潰れ、 反応できる場所がますますセパレータ側へ押し出される — この正のフィードバックが、 電極の奥が実質的に「使えない」状態を作ります。

SPMe はこれをまったく表現できません。$j$ を一様と決め打ちしているため:

• 局所的に大きい $j$ が生む余分な過電圧(と電圧の早落ち)を捉えられない。
• 負極表面の局所電位が Li 析出(lithium plating)の条件に達するかどうかという 劣化解析には使えない($\phi_s - \phi_e$ の分布を持っていないため)。
• 電極厚さを増やした「厚膜電極」の設計検討では、まさに反応分布が主役になるため SPMe では答えが出ない。

また本章では $\kappa^{\mathrm{eff}}$ と $D_e$ を初期濃度 $c_{e0}$ で固定しました。 濃度が半分になれば $\kappa$ も大きく下がるので、勾配が深い高レートではオーム項も 本当は時間・位置依存です。この近似の緩和(濃度依存物性)も DFN で行います。

どれくらい効くのか、数字で見ておきます(いずれも本教材と同じ Chen2020 セル)。 本章の定数物性 SPMe は 2C で早くも正極端の枯渇を予測します(図3.2)。 一方、濃度依存物性まで入れた PyBaMM の SPMe は 2C なら 4.76 Ah 走り ($D_e$ は薄い電解液ほど大きく、枯渇を自己抑制するため)、3C では枯渇して 0.25 Ah で 止まります。そして反応分布まで解く DFN は、同じ 3C でも反応を組み替えて 2.30 Ah まで放電できます。「電解液を入れるか」だけでなく 「物性の濃度依存」「反応分布」がそれぞれ1段ずつ答えを変える — この階段をフェーズ5の図5.2で一望します。

章末:理解度チェックと次章への橋渡し