個数密度関数法(PBM) — 理論と支配方程式
概要
先生、個数密度関数法(PBM)って何ですか?
PBM(Population Balance Model)は、気泡・液滴・粒子のサイズ分布の時間空間変化を追跡するモデルだ。合体(coalescence)、分裂(breakup)、核生成(nucleation)、成長(growth)によってサイズ分布が変化する過程を記述する。
どんな場面で使うんですか?
気泡塔での気泡径分布、乳化プロセスでの液滴径分布、晶析での結晶サイズ分布、エアロゾルの粒径変化など、分散相のサイズ分布が重要な問題全般で使う。Euler-Euler法と組み合わせて使うのが一般的だ。
支配方程式
PBMの方程式を教えてください。
体積 $V$ を持つ粒子(気泡・液滴)の個数密度 $n(V, \mathbf{x}, t)$ の輸送方程式は次の通りだ。
右辺の各項は以下を表す。
- $B_{coal}$: 合体による生成(2つの小さな粒子が合体してVのサイズになる)
- $D_{coal}$: 合体による消滅(サイズVの粒子が他と合体してより大きくなる)
- $B_{break}$: 分裂による生成(大きな粒子が分裂してVのサイズが生まれる)
- $D_{break}$: 分裂による消滅(サイズVの粒子が分裂して小さくなる)
合体と分裂の速度はどうモデル化するんですか?
代表的なクロージャモデルを示そう。
| プロセス | モデル例 | 駆動メカニズム |
|---|---|---|
| 合体頻度 | Prince & Blanch (1990) | 乱流衝突 + 膜排出 |
| 合体頻度 | Luo (1993) | 乱流エネルギー |
| 分裂頻度 | Luo & Svendsen (1996) | 乱流渦による分裂 |
| 分裂頻度 | Martínez-Bazán (2010) | 慣性-表面張力バランス |
Sauter平均径
サイズ分布から代表径を求めるにはどうするんですか?
Sauter平均径 $d_{32}$ が最もよく使われる。
$m_k = \int_0^\infty V^{k/3} n(V) dV$ は分布の $k$ 次モーメントだ。$d_{32}$ は体積-表面積比に対応し、物質移動や反応速度の評価に適している。
ライト兄弟は最初の「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-ω SSTは精度とコストのバランスが良く、多くの実務で最初の選択肢になります。
CFDメッシュの品質管理や乱流モデルの選定に悩む時間を、もっと創造的な設計作業に使えたら。 — Project NovaSolverはそんな実務者の声から生まれました。
Project NovaSolver — CAE実務の課題に向き合う研究開発
「個数密度関数法(PBM)をもっと効率的に解析できないか?」——私たちは実務者の声に耳を傾け、既存ワークフローの改善を目指す次世代CAEプロジェクトに取り組んでいます。具体的な機能はまだ公開前ですが、開発の進捗をお届けします。
進捗通知を受け取る →