熱特性のばらつきを考慮した
電池の熱暴走CAEモデリング

LIBの安全性（熱暴走）シミュレーションは数十年にわたり多数報告されている^1),2)。これらの報告では、セルの比熱や熱伝導率といった熱物性値に加え、電極の DSC や ARC 測定により取得した発熱プロファイルを用い、環境への放熱を考慮した熱伝導方程式に基づいてセル温度の時間変化を予測している。モジュールやパックでの類焼発生の有無を評価するためには、セルモデルを拡張して複数セルを配置した解析を行うのが一般的である。

また、LIBの高容量化や大規模なモジュールでの運用を考慮した場合、熱暴走時の噴出ガスによる周囲への伝熱は大きくなると予想される。そこで、我々は電解液の気化や正極の熱分解に伴う酸素放出など、熱暴走時のガス流れの計測および検証を含めた高精度シミュレーションモデルの構築を進めてきた^3),4)。

いっぽうで、シミュレーションの入力となる熱パラメータには測定ばらつきが存在することが知られており、特定の電極材料を対象にしたDSC測定では同一条件下での測定でも潜熱が12〜25%程度ばらつくという報告もある⁵⁾。しかしながら、過去多くの報告では代表値を用いた確定論的アプローチが採用されており、セル間での類焼挙動に対するパラメータのばらつきの影響は明確ではない。この課題を踏まえ、実測に基づく熱パラメータのばらつきを従来の熱暴走シミュレーションモデルに対して取り込むことで、物性ばらつきによる影響を評価可能なアプローチを開発した。

図１にモデル構築の流れを示す。まず、セルの熱容量、熱伝導率、電極のDSC測定を行った。ここでは、これら熱特性に関わる4つのパラメータ（比熱、熱伝導率、電極分解反応の活性化エネルギー、潜熱）に着目し、各熱特性が正規分布の測定のばらつきがあると仮定し、これらの確率分布を考慮したLIB熱暴走モデルの解析をすることで、パラメータごとにセルの温度プロファイルを得た。これらの解析結果から、各パラメータでのセルの熱暴走時刻や最大温度のデータセットを整理した。

また、上記のデータセットを用いて機械学習によるサロゲートモデルを構築することで、これら熱パラメータがセルの熱暴走時刻および最高温度に与える影響度の評価を行った。これにより、出力分布の解釈が容易になるとともに、影響度の大きいパラメータを抽出することが可能となる。

ヒータ加熱試験を再現する解析モデルとして、角型リチウムイオン電池3セルを対象に、1次元の熱連鎖モデルを構築した。図２に解析ジオメトリを示す。トリガーセル（Cell#1）に一定の発熱量を与え、熱分解反応による自己発熱、隣接セルへの熱伝搬、環境への放熱を計算する。セル温度T_iに対するエネルギー方程式は以下の通りである(添え字のiはセル番号を示す)。

ここでρ、C、kはそれぞれセル密度、比熱、熱伝導率を表し、q_i,in,q_i,reac,q_i,ht,q_i,radはそれぞれヒータによる入熱量、熱分解反応熱、対流による放熱量、輻射による放熱量を表す。

熱分解反応の未反応度c_jおよび反応熱q_reacは以下で仮定した。

ここで、k_jは反応速度定数、Q_j,totは反応潜熱である(添え字のjは反応ピークの番号を示す)。

対流による熱輸送q_i,htは環境温度T_ambとセル温度T_iの差に比例し、輻射による熱輸送q_i,radはStefan-Boltzmann式に従うと仮定する。

ここで、hは対流熱伝達率、aは比表面積、εは輻射率、σはStefan-Boltzmann定数を表す。

なお、試験系に依存する対流熱伝達率, 輻射率、セル間接触熱抵抗は、既報の実測類焼試験結果に基づきパラメータ同定および妥当性検証を行った^4),6)。

B-2章で構築したモデルから熱パラメータのパラメトリック計算を実施し、測定値のばらつきが各セルの温度プロファイルに与える影響を評価した。図４にトリガーセル（#1）、中央セル（#2）、端部セル（#3）の温度プロファイルを示す。実測は太線、解析結果400ケースは薄い線で示す。

トリガーセルでは温度上昇がほぼ同時刻に生じ、ばらつきも小さい。これは、トリガーセルの熱暴走時刻がヒータ加熱量に支配され、熱特性パラメータのばらつきの影響を受けにくいためと考えられる。

中央セルでは、実測で約300sで温度が急上昇するが、シミュレーションではトリガーセルと比較して大きなばらつきが見られる。端部セルでは、シミュレーションにおける熱暴走の発生時刻のばらつきがさらに大きくなっており、トリガーセルから遠方のセルほどばらつきが拡大していることが分かる。

図５に各セルの熱暴走時刻（t_TR）と最大温度（T_max）のヒストグラムを示す。t_TRは、トリガーセル、中央セル、端部セルの順に遅れており、熱暴走が伝播していることが分かる。

最大温度のヒストグラムには、解析を行った全ケースの結果をCase1:トリガーセルのみ熱暴走したケース、Case2：トリガーセルと中央セルが熱暴走したケース、Case3：全てのセルで熱暴走したケースに分類し、それぞれのケースの最大温度を整理した。この図からは、熱暴走するセルの数が増えるほど最大温度も高温側で分布することが分かる。また、中央セルの温度を見ると約3%のケースで最大温度が180℃程度となっており熱暴走を起こす温度まで上昇していないことが分かった。さらに端部セルの場合では約25%が熱暴走を起こす温度まで上昇しなかったことが分かる。これは、入力分布の中で発熱量や伝熱量が小さいパラメータセット(例えば、潜熱や熱伝導率が小さい組み合わせなど)では、中央セルや端部セルが熱暴走に至る温度まで上昇せず、セル間や外気への熱移動によって冷却されていることを示していると考えられる。

以上のことから、測定値のばらつきがセル間で伝播し、熱暴走の有無に大きく影響することが示唆される。

本章では上記の温度プロファイルのデータベースを活用し機械学習を用いた予測モデルを構築することで、熱パラメータの感度を評価する。図６に、熱暴走モデルの計算値とサロゲートモデル予測値の散布図を示す。各プロットはグラフ上で概ね45°線上にあり、サロゲートモデルが各セルのT_maxおよびt_TRを良好に予測していることが分かる。トリガーセルのt_TRの決定係数R²はやや低い（R²=0.555）が、t_TRの変動幅が小さいためであり、平均二乗誤差RMSEは他セルと同程度であることが分かる。

図７は、それぞれT_maxとt_TRに対する各入力のSobolインデックスを示している(インデックスの値が大きいほど、感度が大きい)。最大温度T_maxと熱暴走時刻t_TRのいずれにおいても、熱伝導率kと潜熱Q_totの影響が大きく、さらにトリガーセルからの距離に応じて連続的に変化していることが分かる。これは、伝熱・発熱特性が温度-時間（位置）プロファイルに影響することを示している。いっぽう、比熱Cと活性化エネルギーEのSobolインデックスは低く、これらの測定値のばらつきが T_maxや t_TRの変動に与える影響は比較的小さい。これは、仮定した測定ばらつきの範囲が小さい(実測精度が高い)ためであり、C は平均値比で約1.2%、E は約0.16% の標準偏差と非常に小さい設定であった(ほかのパラメータは10%程度)。このため、これらの範囲では T_maxと t_TRに対する影響も小さくなったと考えられる。