層流管内流れ（Hagen-Poiseuille）
層流管内流れ（Hagen-Poiseuille）の理論基礎
概要
先生、Hagen-Poiseuille 流れってCFDの最初のベンチマークに使われますよね？
その通り。円管内の完全発達層流は、Navier-Stokes方程式の厳密解が存在する数少ない問題の一つだ。CFDコードの検証に最適で、離散化誤差を理論解と直接比較できる。
支配方程式と厳密解
導出を教えてください。
定常・軸対称・完全発達の条件下で、円柱座標の軸方向NS方程式は次のように簡略化される。
壁面 $r = R$ で $u = 0$（no-slip）、中心 $r = 0$ で $du/dr = 0$（対称性）の境界条件で解くと、
これが Hagen-Poiseuille の放物線型速度分布だ。最大速度は $U_{max} = R^2 (-dp/dx) / (4\mu)$ で、管中心に現れる。
体積流量と平均速度
体積流量はどうなりますか？
速度分布を積分すると、
これが Hagen-Poiseuille の法則だ。流量は半径の4乗に比例する。直径を2倍にすると流量は16倍になる。
平均速度は $\bar{U} = Q / (\pi R^2) = U_{max} / 2$ であり、最大速度の半分だ。
摩擦係数
管摩擦係数も理論的に求まるんですよね？
Darcy-Weisbach の摩擦係数 $f$ は、
Fanning 摩擦係数 $C_f$ を使う流儀もある。$C_f = f/4 = 16/Re_D$ だ。混同しやすいので注意が必要だ。
$f = 64/Re$ は有名ですよね。Re = 2300 を超えると乱流に遷移するんでしたっけ。
一般的にはそう言われるが、正確にはパイプ入口の擾乱レベルに依存する。擾乱が極めて小さい場合は Re = $10^5$ 程度まで層流が維持された実験もある。逆に、入口に大きな擾乱があると Re = 2000 以下でも遷移が起こりうる。
助走区間
完全発達するまでの助走区間ってどのくらいですか？
層流の助走区間長さは次の式で見積もれる。
Re = 100 なら $L_e \approx 6D$、Re = 2000 なら $L_e \approx 120D$ だ。CFDで完全発達流を得たいなら、この長さ以上のパイプを計算するか、周期境界条件を使って「無限に長いパイプ」を模擬する。
周期境界条件を使う方が計算コストは低そうですね。
そうだ。流れ方向に周期境界を設定し、圧力勾配をソース項として与える。OpenFOAM では cyclicAMI 境界と fvOptions の meanVelocityForce で実装できる。
毛細血管とHagen-Poiseuille式——医師も使う流体力学
Hagen-Poiseuille式は「流量は管径の4乗に比例する」という強烈な関係を示しています。直径が半分になると流量は1/16になる。これ、実は循環器科の医師が動脈硬化の深刻さを患者に説明するときに使う論理と同じです。血管が20%狭窄しただけで流量は約59%に落ちる計算になる。CFDでも血管狭窄の解析にHP式は登場しますが、「生体血管は非ニュートン流体で脈動がある」という点で厳密解との乖離が生じます。現場では「HP式で粗見積もり→CFDで精密評価」という二段階アプローチが定番です。
層流管内流れ（Hagen-Poiseuille）の数値計算手法
数値手法
理論解があるのにCFDで解く意味は何ですか？
主に2つの目的がある。
1. コード検証（Code Verification）: 理論解と比較して離散化誤差のオーダーを確認
2. 乱流遷移の研究: 層流解を基底状態として擾乱を加え、遷移過程を追跡
離散化誤差の評価
理論解と比較するとき、具体的にどう評価するんですか？
メッシュを3水準以上で系統的に細かくし、誤差の収束オーダーを確認する。
例として、20セル、40セル、80セルの径方向分割で計算し、誤差を両対数プロットする。2次精度スキームなら傾き $-2$、つまり $\epsilon \propto h^2$ で誤差が減少するはずだ。
| 径方向セル数 | $\Delta r / R$ | $L_2$ 誤差（速度） | 収束オーダー |
|---|---|---|---|
| 10 | 0.1 | $2.5 \times 10^{-3}$ | — |
| 20 | 0.05 | $6.3 \times 10^{-4}$ | 1.99 |
| 40 | 0.025 | $1.6 \times 10^{-4}$ | 1.98 |
| 80 | 0.0125 | $4.0 \times 10^{-5}$ | 2.00 |
収束オーダーが理論値（2次精度なら2）と一致すれば、コードが正しく実装されていると言えるんですね。
その通り。これが Method of Manufactured Solutions (MMS) と並ぶ Code Verification の基本手法だ。
軸対称モデル vs フル3D
円管だから軸対称で計算できますよね？
ただし、助走区間を含む場合や乱流遷移を調べる場合はフル3D計算が必要だ。助走区間の流れは軸対称だが、遷移は3D的な擾乱で駆動される。
OpenFOAM での周期的パイプ流
周期境界条件を使った無限長パイプの設定を教えてください。
OpenFOAM での設定例だ。
```
boundary:
inlet: { type cyclic; neighbourPatch outlet; }
outlet: { type cyclic; neighbourPatch inlet; }
wall: { type wall; }
system/fvOptions:
momentumSource:
{
type meanVelocityForce;
selectionMode all;
fieldName U;
Ubar (1 0 0); // 目標平均速度
}
```
この設定では、入口と出口が周期的に接続され、圧力勾配が自動的に調整されて指定した平均速度 $\bar{U} = 1$ m/s が維持される。定常層流なら simpleFoam で数百反復で収束する。
Fluent での設定
Fluent だとどうなりますか？
Fluent でも periodic boundary が使える。Setup → Boundary Conditions で inlet/outlet を periodic pair として定義し、Mass Flow Specification で目標流量を指定する。あるいは、パイプ入口に一様流を与えて十分長いパイプ（$L > 0.06 Re \cdot D$）を計算する方法もある。後者のほうが助走区間の発達も同時に検証できる。
「助走区間」を計算し忘れて大失敗
管内流れの数値解析でやりがちなミスの1つが「入口の助走区間を無視する」こと。Hagen-Poiseuille式が成立するのは流れが完全に発達した後の話で、実際には入口から直径の約50〜100倍の長さが必要です。工場の配管解析で「理論値と10%ずれる！」と騒いでいたら、計測点が入口から管径の20倍しか離れていなかった——というのは笑えない実話です。境界条件に「放物線プロファイル vs 一様流」を選べるのも、まさにこの助走区間問題が理由です。入口条件を雑に設定すると、解が「正しい計算式を使った間違い」になります。
層流管内流れ（Hagen-Poiseuille）の実務適用
実践的な解析手順
Hagen-Poiseuille 流れを使ってCFDコードを検証する手順を教えてください。
ステップバイステップで説明しよう。
1. パラメータ設定: $D = 1$ m, $L = 10D$, $\bar{U} = 1$ m/s, $\nu = 0.01$ m$^2$/s → $Re = 100$
2. 理論解の計算: $U_{max} = 2\bar{U} = 2$ m/s, $-dp/dx = 8\mu\bar{U}/R^2 = 0.32$ Pa/m
3. メッシュ生成: 径方向10, 20, 40セルの3水準。軸方向はセルアスペクト比 $\leq 5$ に
4. 定常計算実行: simpleFoam (OpenFOAM) または Pressure-Based Steady (Fluent)
5. 出口断面の速度プロファイル抽出: 理論放物線と重ねてプロット
6. $L_2$ 誤差の計算と収束オーダー確認: 2次精度なら $p \approx 2$
7. 摩擦係数の確認: 壁面せん断応力から $f$ を計算し $64/Re$ と比較
壁面せん断応力の確認
壁面せん断応力の理論値はいくつですか？
壁面での速度勾配から、
上の数値例（$\mu = \rho \nu = 0.01$ Pa$\cdot$s として $\rho = 1$ kg/m$^3$）では $\tau_w = 0.08$ Pa となる。CFDの結果がこの値と一致するか確認する。
メッシュ品質の注意点
パイプのメッシュで気をつけることはありますか？
円管メッシュには特有の問題がある。
- O型メッシュの中心特異性: 中心軸でセルが潰れる。対策として butterfly型（中心に正方形ブロック＋周囲にO型）を使用
- 壁面近傍の解像度: 放物線プロファイルは壁面で最も勾配が大きい。壁面に近いほどメッシュを細かくする（grading 比 $1.1\text{--}1.3$）
- 軸方向分割: 完全発達流では軸方向は均一分割でよいが、助走区間を含む場合は入口近傍を細かくする
butterfly 型メッシュってどういう形ですか？
管断面の中心に小さな正方形ブロックを配置し、その周囲をO型ブロックで囲む形だ。blockMeshDict では5つのブロック（中心1 + 周囲4）で構成する。これにより中心軸での特異性を回避できる。
代表的な検証結果
正しく計算できた場合のチェックポイントを教えてください。
Re = 100、パイプ長さ $10D$ の場合：
| 検証項目 | 理論値 | 許容誤差 |
|---|---|---|
| 出口中心速度 $U_{max}$ | $2\bar{U}$ | $< 0.1\%$（十分細かいメッシュ） |
| 壁面せん断応力 $\tau_w$ | $8\mu\bar{U}/D$ | $< 1\%$ |
| 圧力降下 $\Delta p$ | $128\mu L Q / (\pi D^4)$ | $< 0.5\%$ |
| 摩擦係数 $f$ | $64/Re = 0.64$ | $< 1\%$ |
| 速度プロファイルの $L_2$ 誤差 | — | 収束オーダー $\geq 1.9$（2次精度時） |
半導体工場の超精密配管——層流管理が命取り
半導体製造の洗浄工程では、超純水やフッ酸を管内で完全に層流状態で流すことが求められます。Re数が2000を超えると微細な粒子が舞い上がってウェハを汚染するため、配管設計はHagen-Poiseuille理論の守備範囲内に収めるのが鉄則です。あるファブでは流量増加のために配管径を少し細くしたところ、流速が上がってRe数が臨界値を超え、歩留まりが突然悪化した。原因特定に数週間かかったという話があります。「層流管内流れは地味な理論」と思いがちですが、数千億円の装置を動かす現場では死活問題です。
層流管内流れ（Hagen-Poiseuille）のソフトウェア比較
ツール別の実装
Hagen-Poiseuille 流れは簡単な問題ですけど、ツールによって違いはありますか？
基本的にどのツールでも解けるが、検証作業での使い勝手に差がある。
| ツール | 周期境界 | 理論解との比較機能 | 壁面量の出力 |
|---|---|---|---|
| Ansys Fluent | 対応（Periodic BC） | Custom Field Function で誤差計算可能 | Wall Shear Stress 直接出力 |
| STAR-CCM+ | 対応 | Field Function で理論解定義可能 | Wall Shear Stress モニター |
| OpenFOAM | cyclic/cyclicAMI | postProcess で sample + Python比較 | wallShearStress ユーティリティ |
| COMSOL | 対応 | 式で理論解を直接入力・プロット | 表面積分で直接計算 |
| FEniCS | 周期BC対応 | Python で直接比較 | 表面積分 |
COMSOL での実装
COMSOL はGUIで簡単に設定できそうですね。
COMSOL は教育目的に非常に向いている。Laminar Flow (spf) モジュールで、
1. 2D Axisymmetric モデルを選択
2. 矩形領域（$r: 0 \to R$, $x: 0 \to L$）を定義
3. 入口に Normal Inflow Velocity、出口に Pressure = 0
4. 壁面に No Slip、対称軸に Axis
5. メッシュを Free Triangular で生成し、壁面にBoundary Layer
6. 定常計算を実行
COMSOL の強みは、結果のプロット画面で直接 $u_{exact}(r) = 2 U_{avg} (1 - (r/R)^2)$ を入力して数値解と重ねてプロットできることだ。
配管系の設計計算への応用
実務では単純なパイプ流だけじゃなく、配管系全体を計算することもありますよね。
配管系の圧力損失計算では、Hagen-Poiseuille の法則が基本になる。直管部の摩擦損失は、
これにエルボ、バルブ、拡大管などの局所損失（$\Delta p_l = K \rho \bar{U}^2 / 2$）を加算する。CFDで配管系全体を解くよりも、1Dの圧力損失計算のほうが効率的な場合が多い。
CFDは局所的な流れの詳細を知りたいときに使うわけですね。
そうだ。例えば、バルブ近傍の剥離域やエルボの二次流れのような局所現象は、CFD以外では正確に予測できない。一方、100mの直管の圧力損失なら理論式で十分だ。
「HP式で計算できるのに、なぜソフトを買うの？」への答え
「Hagen-Poiseuille式は手計算できるのに、なんでFluent使うの？」という質問、新入社員から来ることがあります。正直な答えは「単純な直管はそうだけど、実際の配管はそんなじゃない」ということです。T字分岐、エルボ、バルブ、段差——これらが絡むと途端にHP式は使えなくなる。現場の配管は「HP式が成立しない部分の集合体」です。商用ツールは「この複雑な形状を丸ごと一発で計算したい」というニーズに応えるもので、HP式はその検証の出発点として使います。ツールの選定根拠を問われたとき、「理論値との一致を確認できた」というエビデンスは提案書の説得力を上げます。
層流管内流れ（Hagen-Poiseuille）の先端研究
層流-乱流遷移
先生、パイプ流の遷移って実はまだ完全には解明されていないんですよね？
その通り。パイプ流の遷移は流体力学の未解決問題の一つだ。Hagen-Poiseuille 流は全てのReで線形安定（線形安定性解析で不安定モードが見つからない）だが、実際にはRe $\approx 2000\text{--}2300$ で乱流に遷移する。これは「亜臨界遷移」と呼ばれる。
線形安定なのに乱流になるって矛盾してませんか？
有限振幅の擾乱が必要なんだ。線形安定性解析は微小擾乱しか扱わない。有限振幅の3D擾乱が加わると、非線形効果で「遷移的成長（transient growth）」を経て乱流に遷移する。
パフとスラグ
遷移のときに特徴的な構造が出るんですよね。
Re $\approx 2000$ 付近では「パフ（puff）」と呼ばれる局在化した乱流構造が現れる。
- パフ (Re $\approx 1700\text{--}2300$): 乱流の塊が上流に鋭い前面、下流になだらかな後面を持つ。流れ方向に移流されつつ、分裂や消滅を繰り返す
- スラグ (Re > $2300$ 程度): 乱流領域が上流・下流の両方向に拡大。最終的に管全体が乱流化
Avila et al. (2011, Science) はDNSと実験で、パフの分裂・消滅が統計的過程であることを示し、臨界Re $\approx 2040$ を決定した。これは統計力学の directed percolation と同じ普遍性クラスに属する。
DNS による遷移研究
DNS でパイプ遷移を計算するにはどうすればいいですか？
スペクトル法が標準的だ。Willis & Kerswell のコード「openpipeflow」がオープンソースで公開されている。円柱座標のFourier（軸方向・方位角方向）+ Chebyshev/有限差分（径方向）の組み合わせだ。
有限体積法のCFDコード（OpenFOAM等）でもDNS級の計算は可能だが、スペクトル法に比べて同じ精度を達成するのに $5\text{--}10$ 倍のメッシュ点数が必要になる。
マイクロ流路への応用
マイクロ流体デバイスでも Hagen-Poiseuille の法則は使えますか？
基本的には使える。マイクロ流路（$D \sim 10\text{--}100 \, \mu$m）では Re が低い（$Re \ll 1$ が典型的）ので、慣性項は完全に無視でき、Stokes 流れの世界だ。
ただし以下の効果が重要になる場合がある。
- スリップ流れ: Knudsen 数 $Kn = \lambda / D > 0.01$ で壁面スリップが顕著に。連続体仮定の限界
- 電気浸透流: 壁面の電荷による Debye 層の効果。Hagen-Poiseuille 流にプラグ流が重畳
- 非ニュートン効果: 血液などの非ニュートン流体では速度プロファイルが放物線から逸脱
マイクロ流路のCFDには COMSOL がよく使われますよね。
そうだ。COMSOL の Microfluidics モジュールでは、電気浸透流、電気泳動、拡散-反応を含むマルチフィジックス解析が可能だ。
マイクロ流体デバイスと層流の逆説的な利用
マクロスケールでは「乱流のほうが混合効率が高い」のが常識ですが、マイクロ流体チップ（LOC）の世界では層流が武器になります。幅100μmの流路ではRe数が1以下になることも珍しくなく、二種類の液体を流しても界面でしか混合が起きない。この「層流の混合しにくさ」を逆手にとって、DNA検査チップでは試薬をきれいに層状に維持して段階的に反応させます。もし乱流で混ざってしまったら試薬が無駄になる。HP式が支配するスケールで「意図的に混合させない設計」というのが面白い最前線です。
層流管内流れ（Hagen-Poiseuille）のトラブル対応
よくある問題と対策
こんな簡単な問題で失敗することあるんですか？
意外とあるんだ。Hagen-Poiseuille は「うまくいくはず」の問題だからこそ、合わないときに問題の切り分けがしやすい。
1. 速度プロファイルが放物線にならない
確認ポイント:
- パイプが十分長いか: $L > 0.06 Re \cdot D$ でなければ、出口断面で速度は完全発達していない。Re=100 なら $L > 6D$ が必要
- 入口条件: 一様流を入口に与えた場合、助走区間を経て放物線に発達する。入口断面で放物線を期待してはいけない
- メッシュの軸対称性: 非構造メッシュの場合、断面内の非対称性がプロファイルを歪ませることがある
2. 圧力降下が理論と合わない
圧力降下が理論値の $128\mu LQ / (\pi D^4)$ とずれるんですが。
確認ポイント:
| 原因 | 対策 |
|---|---|
| 助走区間の付加圧力降下 | 完全発達領域内の2点間で $\Delta p$ を計算 |
| 単位系の不整合 | $\mu$ と $\nu$ の単位を確認。$\mu = \rho \nu$ |
| 出口境界条件 | Pressure Outlet ($p = 0$) が出口にあるか確認 |
| 数値拡散 | 1次精度スキームでは圧力勾配の精度も低下 |
3. 収束が遅い / 収束しない
Hagen-Poiseuille は線形問題なので、適切に設定すれば数百反復以内に収束するはずだ。収束が遅い場合：
- 緩和係数が低すぎる: SIMPLE法の場合、圧力緩和 $\alpha_p = 0.3$、速度 $\alpha_U = 0.7$ が標準的。不必要に下げていないか確認
- メッシュの歪み: 非直交性が高いと収束が悪化。
checkMeshでnon-orthogonality を確認 - 境界条件の不整合: 入口で質量流入があるのに出口が壁面になっている等の致命的ミス
4. 壁面せん断応力がゼロ
wallShearStress がゼロと出力されるんですが。
壁面条件が slip になっていないか確認すること。noSlip（OpenFOAM では fixedValue uniform (0 0 0)）が正しい。
また、OpenFOAM の wallShearStress post-processing ユーティリティは、乱流モデルが有効でないと正しく計算されない場合がある。層流の場合は wallShearStress ではなく、wallGradU から手動で $\tau_w = \mu (\partial u / \partial n)_{wall}$ を計算するほうが確実だ。
5. 周期境界条件が機能しない
OpenFOAM で cyclic 境界を使う際の注意点：
- inlet/outlet パッチの面の法線ベクトルが反対向きで、面の形状が完全に一致している必要がある
createPatchユーティリティで周期パッチを作成した後、checkMeshでエラーがないことを確認meanVelocityForceの方向ベクトルが流れ方向と一致しているか確認
基本的な問題こそ丁寧に検証する姿勢が大事なんですね。
その通り。Hagen-Poiseuille 流が正しく解けないなら、それ以上複雑な問題を解く資格がないということだ。CFDの「健康診断」として定期的にベンチマークを回すことを勧める。
「Re数が正しいのに解が合わない」のよくある罠
層流管内流れの検証で「Re数は2000以下にしたのに速度プロファイルが放物線にならない」という相談がたまにあります。大抵の原因は「壁の粗さモデルが有効になっている」か「出口境界条件が逆流を許していない」のどちらかです。特に後者は盲点で、出口でゼロ勾配条件を使っても、数値的な揺らぎで流れが逆流すると収束が壊れます。もう一つの罠は「メッシュが粗すぎて放物線を十分解像できていない」ケース。断面方向に最低10要素は欲しいところを、5要素で「なんとなく形になってる」と思い込む。層流だから簡単、という先入観が落とし穴を招きます。
関連トピック
なった
詳しく
報告