化学反応速度論の基礎 — 理論と支配方程式

カテゴリ: 流体解析(CFD) | 2026-01-15
chemical-kinetics-theory
理論と物理の世界へ

概要

🧑‍🎓

先生! 化学反応速度論って、燃焼シミュレーションの根っこにある理論ですよね?


🎓

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


🧑‍🎓

具体的にはどんな式が中心になるんですか?


🎓

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


支配方程式

🧑‍🎓

まずはArrhenius式から教えてください。


🎓

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


$$ k(T) = A\,T^n \exp\!\left(-\frac{E_a}{R\,T}\right) $$

🎓

ここで $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$ の時間変化は、全素反応の寄与を合算して次のように書ける。


$$ \frac{dY_i}{dt} = \frac{M_{w,i}}{\rho} \sum_{r=1}^{N_r} (\nu''_{i,r} - \nu'_{i,r})\, k_{f,r} \prod_{j} [X_j]^{\nu'_{j,r}} $$

🧑‍🎓

$\nu'$ と $\nu''$ は反応物と生成物の化学量論係数ですね。


🎓

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


Stiff化学反応系

🧑‍🎓

燃焼の化学反応って、時間スケールが極端に違うと聞きました。


🎓

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


🧑‍🎓

そんなに差があると、普通の陽解法では解けないですよね。


🎓

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


反応機構の階層

🧑‍🎓

反応機構にもいろいろなレベルがあるんですか?


🎓

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


区分化学種数反応数代表例用途
グローバル1段2-51-2Westbrook-Dryer概算・初期検討
縮約機構10-3020-100DRM-19, Lu-Law3D RANS/LES
スケルタル機構30-100100-500skeletal-iso-octane詳細RANS
詳細機構50-300+300-3000+GRI-Mech 3.0, USC Mech II0D/1D, DNS
🧑‍🎓

3Dの燃焼LESで詳細機構を直接解くのは難しいということですか?


🎓

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


🧑‍🎓

化学反応速度論の基礎理論、しっかり理解できました。Stiffnessの問題が核心だということがよく分かりました。


🎓

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


Coffee Break よもやま話

ライト兄弟は最初の「CFDエンジニア」だった?

ライト兄弟は1901年に自作の風洞で200以上の翼型を試験しました。当時のコンピュータは? もちろん存在しません。彼らは手作業で揚力と抗力を測定し、最適な翼型を見つけ出した。現代のCFDエンジニアがFluent1発で計算する揚力係数を、ライト兄弟は何百回もの風洞実験で手に入れたのです。

各項の物理的意味
  • 時間項 $\partial(\rho\phi)/\partial t$:蛇口をひねった瞬間を思い浮かべてください。最初は水がバタバタと不安定に出て、しばらくすると安定した流れになりますよね? この「変化している最中」を記述するのが時間項です。心臓の拍動で血流が脈打つのも、エンジンのバルブが開閉するたびに流れが変動するのも、すべて非定常現象。では定常解析とは? 「十分時間が経って流れが落ち着いた後」だけを見る——つまりこの項をゼロにする。計算コストが大幅に下がるため、まず定常で解いてみるのがCFDの基本戦略です。
  • 対流項 $\nabla \cdot (\rho \mathbf{u} \phi)$:川に落ち葉を落としたらどうなりますか? 流れに乗って下流に運ばれますよね。これが「対流」——流体の動きが物を運ぶ効果です。暖房の温風が部屋の端まで届くのも、空気という「運び屋」が熱を対流で輸送しているから。ここが面白いところ——この項は「速度×速度」を含むため非線形です。つまり、流れが速くなるとこの項が急激に強くなり、制御が難しくなる。これが乱流の根本原因です。よくある勘違い:「対流と伝導は同じようなもの」→ 全然違います! 対流は流れが運ぶ、伝導は分子が伝える。桁違いの効率差があります。
  • 拡散項 $\nabla \cdot (\Gamma \nabla \phi)$:コーヒーにミルクを入れて放置したことはありますか? かき混ぜなくても、しばらく経つと自然に混ざりますよね。あれが分子拡散です。では次の質問——ハチミツとお水、どちらが流しやすいですか? 当然お水ですよね。ハチミツは粘性($\mu$)が高いから流れにくい。粘性が大きいと拡散項が強くなり、流体は「もったりした」動きになります。レイノルズ数が小さい流れ(ゆっくり、ドロドロ)では拡散が支配的。逆にRe数が大きい流れでは対流が圧倒し、拡散は脇役になります。
  • 圧力項 $-\nabla p$:注射器のピストンを押すと、液体が針先から勢いよく出ますよね? なぜでしょう? ピストン側が高圧、針先が低圧——この圧力差が流体を押す力になるからです。ダムの放水も同じ原理。天気図で等圧線がギュッと密になっている場所では? そう、強風が吹きます。「圧力差があるところに流れが生まれる」——これがナビエ-ストークス方程式の圧力項の物理的意味。ここでの勘違いポイント:CFDの「圧力」は絶対圧ではなくゲージ圧のことが多い。圧縮性解析に切り替えたとたんに結果がおかしくなる場合、絶対圧/ゲージ圧の混同が原因かもしれません。
  • ソース項 $S_\phi$:暖められた空気が上に昇る——なぜでしょう? 周囲より軽く(密度が低く)なったから、浮力で押し上げられるのです。この浮力はソース項として方程式に追加されます。他にも、ガスコンロの炎で化学反応熱が発生する、工場の電磁ポンプで金属溶湯にローレンツ力がかかる…これらはすべて「外部から流体にエネルギーや力を注入する」作用であり、ソース項で表現します。ソース項を忘れるとどうなるか? 自然対流の解析で浮力を入れ忘れると、流体は一切動かない——冬の部屋で暖房をつけたのに暖かい空気が上に行かない、という物理的にありえない結果になります。
仮定条件と適用限界
  • 連続体仮定:クヌッセン数 Kn < 0.01(分子平均自由行程 ≪ 代表長さ)で成立
  • ニュートン流体仮定:せん断応力と歪み速度が線形関係(非ニュートン流体では粘度モデルが必要)
  • 非圧縮性仮定(Ma < 0.3の場合):密度を一定として扱う。マッハ数0.3以上では圧縮性効果を考慮
  • ブシネスク近似(自然対流):密度変化を浮力項のみで考慮し、他の項では一定密度を使用
  • 適用外ケース:希薄気体(Kn > 0.1)、超音速・極超音速流れ(衝撃波捕捉が必要)、自由表面流れ(VOF/Level Set等が必要)
次元解析と単位系
変数SI単位注意点・換算メモ
速度 $u$m/s入口条件で体積流量から換算する際、断面積の単位に注意
圧力 $p$Paゲージ圧と絶対圧の区別。圧縮性解析では絶対圧を使用
密度 $\rho$kg/m³空気: 約1.225 kg/m³@20°C、水: 約998 kg/m³@20°C
粘性係数 $\mu$Pa·s動粘性係数 $\nu = \mu/\rho$ [m²/s] との混同に注意
レイノルズ数 $Re$無次元$Re = \rho u L / \mu$。層流/乱流遷移の判定指標
CFL数無次元$CFL = u \Delta t / \Delta x$。時間刻みの安定性に直結

数値例:円管内層流(d=10mm, L=1m, 水μ=0.001Pa·s, Q=0.1L/min)

Re = ρuD/μ = 998×0.021×0.01/0.001 ≈ 212層流) 圧力損失 ΔP = 128μLQ/(πd⁴) ≈ 68.2 Pa

乱流モデル別の精度比較(後向きステップ、再付着長さ):

k-ε標準5.8h(実験6.1h)
-4.9%
k-ω SST6.0h
-1.6%
RSM6.05h
-0.8%
LES6.12h
+0.3%
実験値6.1h
基準

k-ω SSTは精度とコストのバランスが良く、多くの実務で最初の選択肢になります。

簡易計算ツール:流体力学基礎

レイノルズ数 Re = ρuL/μ を計算し、層流/乱流の判定を行います。

CFL数 = u·Δt/Δx を計算し、時間刻みの安定性を確認します。

円管内の層流ハーゲン-ポアズイユ流れの圧力損失 ΔP = 128μLQ/(πd⁴)

CFDメッシュの品質管理や乱流モデルの選定に悩む時間を、もっと創造的な設計作業に使えたら。 — Project NovaSolverはそんな実務者の声から生まれました。

Project NovaSolver — CAE実務の課題に向き合う研究開発

「化学反応速度論の基礎をもっと効率的に解析できないか?」——私たちは実務者の声に耳を傾け、既存ワークフローの改善を目指す次世代CAEプロジェクトに取り組んでいます。具体的な機能はまだ公開前ですが、開発の進捗をお届けします。

進捗通知を受け取る →