化学反応速度論の基礎

そのとおり。化学反応速度論は、燃料と酸化剤がどのくらいの速さで反応するかを数学的に記述する学問だ。燃焼CFDでは、この反応速度を正しくモデル化しないと、火炎温度も排出ガス組成もまるで的外れな結果になる。

中核はArrhenius型反応速度式だ。素反応ごとに前指数因子 $A$、温度指数 $n$、活性化エネルギー $E_a$ の3パラメータで速度定数を表す。

反応速度定数 $k$ は次のように書ける。

ここで $A$ は前指数因子（頻度因子）、$n$ は温度指数、$E_a$ は活性化エネルギー [J/mol]、$R$ は一般ガス定数 8.314 J/(mol K) だ。たとえば水素の基本反応 H$_2$ + O$_2$ 系では $E_a$ が数十 kJ/mol のオーダーになる。

そうだ。$E_a$ が大きいほど指数項 $\exp(-E_a/RT)$ が小さくなり、低温では反応がほとんど進まない。逆に温度が上がると指数的に速度が増大する。これが燃焼の着火遅れや火炎の安定性に直結するんだ。

化学種 $i$ の質量分率 $Y_i$ の時間変化は、全素反応の寄与を合算して次のように書ける。

そうだ。$[X_j]$ はモル濃度、$M_{w,i}$ は分子量。複数の素反応が絡み合って、化学種の生成・消滅速度 $\dot{\omega}_i$ を構成する。

まさにそれがStiffness（硬さ）の問題だ。たとえばメタン/空気の詳細反応機構GRI-Mech 3.0は53化学種・325素反応を含むが、ラジカル反応の時定数は $10^{-9}$ 秒オーダー、一方で主反応は $10^{-3}$ 秒オーダーと6桁以上の差がある。

そのとおり。陽的Euler法やRunge-Kutta法では微小な時間刻みが必要になり、実用的な計算時間に収まらない。そこで陰的解法（BDF法、Rosenbrock法など）や、次回の記事で詳しく扱うISAT（In-Situ Adaptive Tabulation）のようなテーブル化手法が使われる。

反応機構には以下の階層がある。

計算コストが膨大になるからね。100化学種を含む詳細機構を100万セルの3Dメッシュで解こうとすると、各セルで100元の連立ODEを毎ステップ解く必要がある。だから実務では縮約機構を使うか、ISATやFGM（Flamelet Generated Manifold）で次元を削減する。

うむ。数値手法の詳細は次の記事で扱う。Arrhenius式のパラメータ設定を誤ると着火温度が数百K狂うこともあるから、実験データとの照合は必ず行おう。

燃焼CFDにおける化学反応の数値積分は大きく3つのアプローチに分かれる。(1) 直接積分（DI）、(2) テーブル化手法、(3) オペレータ分割法だ。

各CFDセルで化学種のODE系を毎タイムステップ解く方法だ。Stiff系に対応するため、陰的多段階法を使う。

代表的な陰的ソルバーを比較しよう。

そうだ。Ansys Fluentでは「Stiff Chemistry Solver」としてCVODEベースの積分器が組み込まれている。STAR-CCM+でもSUNDIALSライブラリを内部で利用している。ヤコビアン行列の評価がボトルネックになるから、解析的ヤコビアンの自動生成が重要だ。

Popeが1997年に提案した手法で、一度計算した化学反応の入出力を二分木構造でテーブルに記録する。新しい組成点が来たとき、既存の記録から線形近似で答えを返せるなら直接積分を省略する。

許容誤差 $\epsilon_{\text{tol}}$ はユーザーが設定する。Fluentのデフォルトは $10^{-4}$ 程度だ。ISATにより直接積分の呼び出し回数を90%以上削減できることも多い。ただしテーブルサイズがメモリを圧迫する場合があり、大規模並列計算では各プロセスが独立にテーブルを持つためメモリ効率が課題になる。

流体輸送と化学反応を別々のステップで解く手法だ。Strang分割が代表的で、CFDの輸送ステップと化学反応のODE積分を交互に実行する。

1ステップの流れはこうなる。

OpenFOAMの reactingFoam ソルバーではオペレータ分割が標準実装されている。化学反応ソルバーとして ode、EulerImplicit、noChemistrySolver が選択可能で、chemistryProperties 辞書で設定する。CONVERGEでもSAGE（Stochastic Adaptive Grid Engine）ソルバーが化学反応の直接積分にオペレータ分割を採用している。

代表的な縮約手法を整理しておこう。

そのとおり。0D/1Dの着火遅れ検証にはCVODEで詳細機構を直接解き、3D LESではDRGEPで縮約した30化学種程度の機構をISATと組み合わせる、というのが実務的な落としどころだ。

実務的な解析フローを順を追って解説しよう。

1. 燃料の特定と反応機構の選定 -- 対象燃料（メタン、n-ヘプタン、水素など）に適した反応機構をCantera等で0D検証する

CanteraやCHEMKIN-PROを使って、定容断熱着火のシミュレーションを行う。初期温度・圧力・当量比を変えて着火遅れ時間 $\tau_{\text{ign}}$ を計算し、実験のショック管データと比較する。ここで反応機構の妥当性を確認してからCFDに進む。

主要なフォーマットを整理しよう。

そうだ。FluentではCHEMKIN形式の反応機構ファイル（.inp）と熱力学データファイル（therm.dat）、輸送特性ファイル（tran.dat）の3つをインポートする。STAR-CCM+ではDARSライブラリ経由でCHEMKIN形式を読み込む。

メタン/空気、当量比1.0、圧力20 atmの場合の着火遅れ時間の目安はこうなる。

DRM-19は19化学種・84反応でGRI-Mech 3.0の主要パスを再現できる。3DのRANS計算にはちょうどいい規模だ。ただしNOx予測が必要な場合は、Zeldovich機構（thermal NOx）の3反応を別途追加する必要がある。

メタンの断熱火炎温度は約2230 Kだが、実際のバーナーでは輻射で200-400 K低下する。これがNOx生成にも大きく影響するから、燃焼解析では輻射モデルは必須と考えたほうがいい。

主要CFDツールの化学反応計算機能を比較してみよう。

CONVERGEは内燃機関向けに開発された経緯があり、ガソリンサロゲート燃料のような数百化学種の詳細機構を3Dで直接扱えるよう最適化されている。AMR（Adaptive Mesh Refinement）と化学反応の自動負荷分散を組み合わせて、HPCクラスタで効率的に計算できる設計だ。

1. Models > Species > Species Transport を有効化

1. chemkinToFoam コマンドでCHEMKIN形式を変換

そうだ。PythonベースなのでJupyter Notebookで着火遅れのパラメトリックスタディを手軽に回せる。3D CFDに投入する前の反応機構スクリーニングにはCanteraが最適だ。

大きく3つの方向性がある。(1) 機械学習による化学反応の高速化、(2) 非平衡効果の考慮、(3) 代替燃料の反応機構開発だ。

DNN（Deep Neural Network）やGPR（Gaussian Process Regression）を使って、化学反応のODE積分結果をサロゲートモデルとして学習させる研究が進んでいる。ISATの代替として、ニューラルネットワークが組成空間の写像を学習する。

最近の成果として注目されているのが以下の手法だ。

訓練データの範囲外への外挿精度が最大の課題だ。燃焼では局所消炎や再着火のような非定常現象で組成空間が大きく変動するから、学習データのカバレッジが重要になる。

Non-equilibrium plasma assisted combustion（PAC）は、ナノ秒パルス放電でO原子やOH・N2の振動励起状態を生成し、着火促進や火炎安定化を実現する技術だ。従来のArrhenius式では記述できない電子衝突解離や振動-並進緩和を反応機構に組み込む必要がある。

水素、アンモニア（NH3）、e-fuel（合成燃料）の反応機構開発が活発だ。

そうだ。アンモニアの層流燃焼速度は水素の約1/5しかない。H2を10-30%混合することで着火特性を改善するが、NOx（特にfuel NOx）の抑制が大きな課題だ。この分野の反応機構はまだ成熟していないから、最新の文献を常にフォローする必要がある。

エネルギー転換期にあって、燃焼の化学反応速度論は再び活発な研究分野になっている。CFDエンジニアも基礎化学の動向に目を配ることが大切だよ。

燃焼CFDの化学反応関連のトラブルは、大きく分けて発散問題、精度問題、性能問題の3カテゴリーに分かれる。

症状: 温度や化学種質量分率が非物理的な値（負の質量分率、10000 K超の温度）になり、ソルバーがクラッシュする。

対策チェックリスト:

数値拡散と化学反応ソース項の競合が原因だ。特にマイナー化学種（ラジカル）で起きやすい。対策としてSpecies Boundingに加え、2次精度の差分スキーム（Fluent: Second Order Upwind、OpenFOAM: limitedLinear）を使うことが重要だ。

考えられる原因と対策:

原因と対策を表にまとめよう。

そうだ。3Dで問題が起きたときも、まず該当条件で0D計算を回して反応機構自体の問題かCFD側の問題かを切り分ける。これが最も効率的なデバッグ手順だ。