キャビティ流れ（蓋駆動）
キャビティ流れ（蓋駆動）の理論基礎
概要
先生、lid-driven cavity ってCFDで一番最初にやるベンチマーク問題ですよね？
その通り。正方形キャビティの上面（蓋）を一定速度で水平にスライドさせる問題だ。幾何形状が単純で境界条件も明確、しかもRe数を変えると豊かな流れ構造が現れる。Ghia et al. (1982) のベンチマークデータが40年以上にわたって使われ続けている。
問題設定
問題の定式化を教えてください。
正方形キャビティ（一辺 $L$）の上壁面が速度 $U$ で $x$ 方向に移動する。他の3壁面は静止。支配方程式は非圧縮Navier-Stokes方程式だ。
無次元パラメータは $Re = UL/\nu$ の一つだけだ。境界条件は、
- 上壁面: $u = U$, $v = 0$
- 下壁面・左壁面・右壁面: $u = 0$, $v = 0$
流れ構造のRe依存性
Re数によってどう変わるんですか？
以下のように整理できる。
| Re | 主渦の位置（中心座標） | コーナー渦 | 流れの性質 |
|---|---|---|---|
| 100 | (0.6189, 0.7344) | 下隅に微小渦 | 定常、主渦がやや上方・右寄り |
| 400 | (0.5547, 0.6055) | 下両隅に渦 | 定常、主渦が中心に近づく |
| 1000 | (0.5313, 0.5625) | 下両隅＋左上に渦 | 定常、主渦がほぼ中心 |
| 5000 | (0.5117, 0.5352) | 全コーナーに渦 | 定常（2D）、3Dでは不安定 |
| 10000 | (0.5117, 0.5313) | 渦の階層構造 | 2Dでは定常/弱非定常、3Dでは乱流化 |
コーナー渦って、角に小さな渦ができるんですよね。
そうだ。Moffatt (1964) が理論的に示したように、鋭角コーナーでは無限の渦列が存在する。各渦の強度は幾何級数的に減少する。CFDでは最初の2〜3個の渦が解像できれば十分だが、メッシュをコーナーに向けて細かくする必要がある。
渦度-流れ関数定式化
2Dの場合は渦度-流れ関数で解くことが多いんですか？
2D非圧縮流では渦度-流れ関数定式化が連続の式を自動的に満たすため効率的だ。
ここで渦度 $\omega = \partial v/\partial x - \partial u/\partial y$、流れ関数 $\psi$ は $u = \partial \psi / \partial y$, $v = -\partial \psi / \partial x$ だ。
ただし3Dや可変密度流れへの拡張が困難なので、汎用CFDコードでは速度-圧力定式化（primitive variable formulation）が標準だ。
蓋の角の特異性
蓋の端で壁面速度が不連続になりますよね。これは問題にならないんですか？
良い指摘だ。蓋の角（$(0,L)$ と $(L,L)$）では速度が不連続になる。$u = U$（上壁面）と $u = 0$（側壁面）が同時に成立しない。これは数学的な特異点であり、メッシュを細かくしても解が収束しない。
実務的な対処法は、
- 角から数セル離れた領域で結果を評価する
- 正則化条件を使う（蓋の速度を角に向かってスムーズにゼロに遷移させる）
- 特異性を承知の上で、十分なメッシュ密度で計算し、角以外の領域の精度を確保する
Ghia（1982年）のデータが40年以上引用される理由
キャビティ流れの「正解データ」として世界中で使われているのが、GhiaらがJournal of Computational Physicsに1982年に発表した論文です。Re数1000、3200、5000、10000の各ケースでの速度プロファイルが表データで掲載されており、今も年間数百〜千本の論文で引用されています。なぜ40年以上通用し続けているのか——それは「問題設定がシンプルで再現性が高く、格子収束されたデータが当時の最高水準で示されているから」です。新しいCFDソルバーや数値スキームを開発したとき、「Ghia基準」に通るかどうかが最初のチェックポイントになっています。
キャビティ流れ（蓋駆動）の数値計算手法
数値手法の選択
キャビティ流れを解くにはどんな手法がいいですか？
キャビティ流れは閉じた領域で入口も出口もないという特殊性がある。圧力の絶対値が一意に決まらないため、参照圧力を1点で固定する必要がある。
| 手法 | 適用範囲 | 備考 |
|---|---|---|
| 有限差分法（等間隔格子） | Re < $10^4$ | Ghia et al. の元データもこの手法 |
| 有限体積法（FVM） | 汎用 | 商用CFDの標準。SIMPLE/PIMPLEで連成 |
| 有限要素法（FEM） | 汎用 | Taylor-Hood要素（P2/P1）でLBB条件充足 |
| スペクトル法 | 高精度計算 | Chebyshev多項式基底。指数的収束 |
| 格子ボルツマン法（LBM） | Re < $10^4$ 程度 | 壁面バウンスバックで実装が容易 |
FVM での実装（SIMPLE法）
SIMPLE法の手順を簡単に教えてください。
キャビティ流れでのSIMPLE反復は以下の通りだ。
1. 圧力場 $p^*$ を仮定
2. 運動量方程式を $p^$ で解いて仮速度 $\mathbf{u}^$ を得る
3. 圧力補正方程式 $\nabla \cdot (\frac{1}{a_P} \nabla p') = \nabla \cdot \mathbf{u}^*$ を解く
4. 速度と圧力を補正: $p = p^ + \alpha_p p'$, $\mathbf{u} = \mathbf{u}^ - \frac{1}{a_P} \nabla p'$
5. 収束判定。未収束なら1に戻る
キャビティ流れでは入口出口がないため、圧力の基準を1点（例えば左下角のセル）で固定する。OpenFOAM では p の fvSolution で reference cell と reference value を指定する。
格子ボルツマン法
格子ボルツマン法（LBM）でも解けるんですか？
LBM は cavity flow の入門として非常に人気がある。D2Q9 モデル（2次元9速度）で、BGK衝突演算子を使えば簡単に実装できる。
緩和時間 $\tau$ と動粘度の関係は $\nu = c_s^2 (\tau - 0.5) \Delta t$ で、$c_s = 1/\sqrt{3}$（格子音速）だ。上壁面の移動は Zou-He 境界条件で実装する。
LBM だと圧力-速度連成を解く必要がないんですね。
そうだ。LBM は陽的手法なので反復が不要で、並列化も容易だ。ただし、高Re数では数値安定性に注意が必要で、MRT（Multi-Relaxation Time）モデルや正則化手法が使われる。
Ghia ベンチマークとの比較方法
Ghia et al. のデータとどう比較すればいいですか？
比較対象は2つだ。
1. キャビティ中央の垂直線（$x = 0.5$）上の $u$ 速度: 上壁面で $u = 1$、底壁面で $u = 0$、中間で再循環による負の $u$ が現れる
2. キャビティ中央の水平線（$y = 0.5$）上の $v$ 速度: 左壁面近くで正、右壁面近くで負
Ghia et al. は $129 \times 129$ の均一格子で Re = 100, 400, 1000, 3200, 5000, 7500, 10000 のデータを tabulate している。これらの離散点の値と自分のCFD結果を比較して、相対誤差が $1\%$ 以内であれば良好だ。
コーナー渦が「圧力-速度連成の鍵穴」になる話
蓋駆動キャビティの数値計算では、4つのコーナー近傍に微小な二次渦・三次渦が入れ子状に生まれます。これらは理論的には無限に小さな渦が続く「Moffatt渦」として知られていますが、数値計算では格子解像度に限界があるため必ず途中で消えます。面白いのは「この入れ子渦を正確に再現できるかどうか」が圧力-速度連成スキームの品質試験になること。SIMPLE法・SIMPLEC法・PISO法などを比較するベンチマークとして、コーナー渦の数と形状が使われることがあります。「コーナーをちゃんと計算できるか」は意外と難しい問題です。
キャビティ流れ（蓋駆動）の実務適用
実践手順
キャビティ流れのベンチマーク計算を実際にやる手順を教えてください。
OpenFOAM を例に説明しよう。
1. メッシュ生成: blockMeshDict で $N \times N$ の均一格子を定義。$N = 64, 128, 256$ の3水準
2. 境界条件: 上壁面 movingWall に $U = (1, 0, 0)$、他3壁面は fixedValue (0 0 0)
3. 物性値: $\nu = 1/Re$（$U = 1$, $L = 1$ の無次元化）
4. ソルバー: icoFoam（層流・非定常）または simpleFoam（定常）
5. 後処理: 中央線上の速度を sampleDict で抽出し、Ghia データと比較
メッシュ収束性
何分割あれば十分ですか？
Re数によって必要な解像度が変わる。
| Re | 最低必要分割数 | Ghia の分割数 | 推奨分割数 |
|---|---|---|---|
| 100 | $32 \times 32$ | $129 \times 129$ | $64 \times 64$ |
| 1000 | $64 \times 64$ | $129 \times 129$ | $128 \times 128$ |
| 5000 | $128 \times 128$ | $257 \times 257$ | $256 \times 256$ |
| 10000 | $256 \times 256$ | $257 \times 257$ | $512 \times 512$ |
均一格子じゃなくて、壁面近くを細かくしたほうがいいんじゃないですか？
その通り。壁面近傍とコーナー部にメッシュを集中させると、同じセル数でも精度が向上する。blockMeshDict の simpleGrading で壁面方向にストレッチをかける（grading ratio 5〜10）のが効果的だ。
Ghia データとの比較プロット
Re = 1000 のときの代表的な値を教えてください。
Ghia et al. (1982) の Re = 1000 データから主要な値を示す。
$x = 0.5$ 上の $u$ 速度:
| $y$ | $u$ |
|---|---|
| 1.0000 | 1.00000 |
| 0.9766 | 0.65928 |
| 0.9688 | 0.57492 |
| 0.5000 | -0.06080 |
| 0.0547 | -0.24533 |
| 0.0000 | 0.00000 |
特に $y \approx 0.05$ 付近の負の $u$（$\approx -0.25$）は、底壁面近くの再循環を表している。この値が合えば、下部コーナー渦の構造が正しく解像できている証拠だ。
3D キャビティ
3Dのキャビティ流れも研究されているんですか？
3Dではスパン方向（$z$方向）の端壁面効果が加わる。Koseff & Street (1984) の実験が有名だ。2Dでは定常な Re = 3200 が、3Dでは非定常になることが知られている。
3D cavity の臨界 Re はスパン方向のアスペクト比に依存するが、SAR (Span Aspect Ratio) = 3 の場合、Re $\approx 800$ で3D不安定性が発現する（Albensoeder & Kuhlmann 2006）。
塗料の撹拌タンクとキャビティ流れの共通点
工場の塗料撹拌タンクは「上蓋を動かす代わりにインペラで液体を動かす」という意味でキャビティ流れに近い物理を持っています。タンク内に「デッドゾーン（流れが淀んで撹拌されない領域）」ができると塗料が不均一になったり沈殿が起きたりします。このデッドゾーンはキャビティ流れでいう「コーナー近傍の低速域」にあたります。蓋駆動キャビティのCFDで「副渦の位置と大きさ」を正確に計算できれば、タンク設計でのデッドゾーン予測精度も上がります。学術的なベンチマーク問題が、現場の製品品質改善に直結するのは珍しいことではありません。
キャビティ流れ（蓋駆動）のソフトウェア比較
ツール別の実装
キャビティ流れはどのCFDツールでも解けますよね。違いはありますか？
ベンチマーク問題としてはどのツールでも解けるが、教育・学習目的での使いやすさに差がある。
| ツール | 設定の容易さ | 精度確認 | 教育利用 |
|---|---|---|---|
| Ansys Fluent | チュートリアルあり | Custom Field Function | 大学ライセンスあり |
| COMSOL | GUI で数分で設定完了 | 理論解・参照解の直接比較可能 | 教育に最適 |
| OpenFOAM | cavity チュートリアルが初学者向け | Python/GNUplot で比較 | 無償で学習に最適 |
| FEniCS/Firedrake | Pythonスクリプトで実装 | 高次FEM、自動微分 | 研究・教育向け |
| Palabos | LBMベース。C++ | cavity サンプルあり | LBM学習に最適 |
OpenFOAM の cavity チュートリアル
OpenFOAM の有名な cavity チュートリアルについて教えてください。
OpenFOAM の入門チュートリアルとして $FOAM_TUTORIALS/incompressible/icoFoam/cavity/cavity が提供されている。$20 \times 20$ の粗い均一格子で Re = 10 のキャビティ流れを解く。
これをベンチマーク用に拡張するには、
- メッシュを $128 \times 128$ 以上に変更（
blockMeshDictの hex の分割数を変更） - $\nu$ を変更して Re = 100, 400, 1000 を計算
sampleDictで中央線データを抽出- Ghia データ（テキストファイルで入手可能）と Python matplotlib で比較プロット
Ansys Fluent での設定
Fluent だとどう設定しますか？
Fluent の場合は以下の手順だ。
1. Meshing: 正方形領域に Face Meshing で $128 \times 128$ の均一格子を生成
2. Setup: Viscous Model → Laminar、Fluid Material → 密度と粘度を設定
3. Boundary Conditions: 上壁面を Moving Wall（$U = 1$ m/s, $x$ 方向）、他は Stationary Wall
4. Solution Methods: SIMPLE、Second Order Upwind
5. Reference Values: 圧力の Operating Point を設定（浮動圧力を避ける）
6. Run: 3000〜5000 iterations で残差が $10^{-6}$ 以下に収束
7. Post: 中央線上のXY Plot で速度プロファイルを確認
COMSOL での設定
COMSOL なら初学者でも簡単そうですね。
COMSOL では数クリックで設定できる。
1. 2D モデル → Laminar Flow (spf) を選択
2. 正方形ジオメトリ（$1 \times 1$ m）を作成
3. 上辺を Wall → Sliding Wall（$U_w = 1$ m/s）に設定
4. 他の3辺は No Slip（デフォルト）
5. Fluid Properties で $\mu$ と $\rho$ を設定（Re に対応）
6. Mapped メッシュで均一分割
7. Stationary Study で求解
8. Line Graph で中央線上の速度を Ghia データと比較
ツールの違いを理解するのに、同じ問題を複数のツールで解いてみるのも勉強になりますよね。
その通り。同じ問題を OpenFOAM, Fluent, COMSOL で解いて結果を比較するのは、各ツールの特性を理解する良い方法だ。離散化スキームの影響もよく分かる。
「Ghia論文との比較」がソルバー入門の登竜門になっている理由
CFDを学び始めた人が最初に取り組む「自分のコードが動いているか確認する」課題として、蓋駆動キャビティとGhia論文との比較は世界標準になっています。理由は単純で「入口・出口境界条件が不要」「3D計算が不要」「実験装置がなくても答えが確認できる」の3つが揃っているから。OpenFOAMのチュートリアルにも含まれていますし、商用ソルバーの入門トレーニングでも必ず登場します。「キャビティでGhia論文と合いました」は「CFD計算の基礎が動いている」の証明書として業界で通用します。シンプルな問題だからこそ、教育・検証・ソルバー評価の3役を一手に担えます。
キャビティ流れ（蓋駆動）の先端研究
安定性解析と分岐
キャビティ流れの安定性解析ってどうやるんですか？
定常解を基底状態として微小擾乱の成長率を調べる。時間発展法（DNS で擾乱のパワー法的成長を追跡）またはグローバル安定性解析（一般化固有値問題を解く）が使われる。
2Dキャビティの Hopf 分岐は Re $\approx 8000\text{--}8500$ で起こる（Fortin et al. 1997、Bruneau & Saad 2006）。この値は格子解像度に敏感で、Re = 8000 とも 8500 とも報告されている。3Dの場合は、centrifugal instability（Taylor-Goertler 型）が Re $\approx 800\text{--}1000$ で現れる。
正則化 lid-driven cavity
蓋の角の特異性を避ける方法はあるんですか？
正則化（regularized）lid-driven cavity では、蓋の速度を角に向かってスムーズにゼロに遷移させる。例えば、
この正則化により角の特異性が消え、メッシュ収束性が改善する。高精度ベンチマーク（Botella & Peyret 1998 のスペクトル法解など）ではこの正則化が使われている。
高精度参照解
Ghia のデータより精度が高い参照解はありますか？
以下の研究がより高精度な参照解を提供している。
| 研究者 | 手法 | Re | 格子/次数 | 精度 |
|---|---|---|---|---|
| Ghia et al. (1982) | FDM | 100-10000 | $257 \times 257$ | 工学的精度 |
| Botella & Peyret (1998) | スペクトル法 | 1000 | $N = 96$ | 12桁以上 |
| Erturk et al. (2005) | FDM（高次） | 100-21000 | $601 \times 601$ | 高精度 |
| Bruneau & Saad (2006) | マルチグリッドFDM | 1000-10000 | $2048 \times 2048$ | 非定常解含む |
Botella & Peyret は12桁の精度ですか。すごいですね。
スペクトル法の指数的収束の威力だ。ただし正則化された lid-driven cavity での結果なので、特異点がある場合のGhia データとは直接比較できないことに注意。
機械学習への応用
キャビティ流れは機械学習の研究でも使われているんですか？
サロゲートモデルの構築やPINN（Physics-Informed Neural Networks）の検証に頻繁に使われている。
- PINN: Raissi et al. (2019) は cavity flow を含むNS方程式のPINN解法を実証
- オートエンコーダ: 異なるRe数の流れ場をlatent spaceで補間
- 強化学習: キャビティ内の能動制御（壁面吹出し/吸込み）の最適化
PINN ってデータなしでNS方程式の解を学習できるんですか？
損失関数にNS方程式の残差を含めるので、原理的にはデータなしで解ける。ただし、実際には境界条件のデータが必要だし、高Re数では収束が困難だ。cavity flow はRe $\leq 1000$ 程度でPINNの動作を確認するのに適した問題だ。
Re=8000を超えると突然3Dになる——キャビティ流れの「転換点」
2Dのキャビティ流れは、Re数を上げていくとある時点で解が3次元化します。実験とLESの研究によると、Re≈8000〜10000程度で流れがスパン方向（奥行き方向）に不均一になり始め、真の3D乱流に移行します。2DのRANS計算ではこの変化は検知できません。建築や電子機器冷却の実設計では筐体内部の「空洞」でまさにこの現象が起きることがあり、「2Dで計算したら均一な混合が出た、でも実機では偏りがあった」という事例の裏に3D遷移が潜んでいることがあります。Ghia論文のデータが2Dでよく使われる一方、3D拡張は今も研究の最前線です。
キャビティ流れ（蓋駆動）のトラブル対応
よくある問題と対策
キャビティ流れの計算で起きるトラブルを教えてください。
1. Ghia データと合わない
確認ポイント:
- メッシュ密度: Re = 1000 なら最低 $64 \times 64$ 以上。$32 \times 32$ だと $5\text{--}10\%$ の誤差が出る
- 空間離散化スキーム: 1次風上（First Order Upwind）だと数値拡散で渦が弱まる。最低2次精度を使う
- 参照点のずれ: Ghia のデータは $129 \times 129$ の格子点上の値。自分のメッシュとの内挿誤差に注意
- Re数の定義: $Re = UL/\nu$ であることを確認。$L$ をキャビティの一辺の長さとしているか
2. 圧力が振動する
圧力のコンター図がチェッカーボード模様になるんですが。
原因: collocated配置（速度と圧力を同じ格子点に配置）で Rhie-Chow 補間が不適切な場合に起こる。
対策:
- Rhie-Chow momentum interpolation が有効になっているか確認（OpenFOAM では標準で有効）
- staggered grid（速度と圧力を異なる格子点に配置）を使う手法もある（古典的FDMコード）
- 圧力の under-relaxation を $0.2\text{--}0.3$ に設定
3. 高Re数で収束しない
Re > 5000 では非線形性が強く、SIMPLE法の収束が遅くなる。
対策:
- 擬似時間進行法: 非定常ソルバー（PIMPLE）で十分長い時間走らせて定常解に収束させる
- continuation法: 低Re（例: Re=100）の収束解を初期条件にして、Re を段階的に上げる
- マルチグリッド: AMG前処理を圧力方程式に使用。OpenFOAMでは
GAMGソルバー
Re を段階的に上げるのは良い方法ですね。
いわゆる parameter continuation だ。Re=100 → 400 → 1000 → 3200 → 5000 → 10000 と段階的に上げる。各段階で前段の解を初期条件にする。
4. コーナー渦が解像できない
Moffatt 渦は幾何級数的に小さくなるので、均一格子では解像が困難だ。
対策: コーナーに向かって格子を幾何級数的に細かくする。OpenFOAM の blockMesh なら simpleGrading で壁面への grading ratio を $10\text{--}20$ に設定。
5. 3D計算で2Dと異なる結果が出る
3Dで計算したら2Dの結果と違うんですが。
Re > 800 程度では3D不安定性（Taylor-Goertler渦）が発現するので、3Dの結果が2Dと異なるのは物理的に正しい。
確認事項:
- スパン方向に周期境界を使っているか（端壁面がある場合は端壁効果が加わる）
- スパン方向の解像度が十分か（不安定波長の $1/10$ 以下のメッシュ幅）
- 2D計算との差が不安定性によるものか、数値誤差によるものかを区別するため、低Re（Re=100）で2Dと3Dが一致することを先に確認
段階的に検証を進める姿勢が大事なんですね。
その通りだ。cavity flow は「簡単な問題」と思われがちだが、高Re数での分岐現象、コーナー特異性、3D不安定性など、数値流体力学の本質的な問題を含んでいる。だからこそ40年以上にわたってベンチマーク問題として使われ続けているんだ。
「中央渦の中心位置がずれる」のはメッシュのせいか数値拡散か
キャビティ流れで「Ghia論文の渦中心位置と自分の計算が合わない」場合、真っ先に疑うのはメッシュ解像度ですが、実はそれ以外にも原因があります。移流項の離散化スキームが低精度（1次精度の上流差分）の場合、数値拡散で渦が「実際より下方向にぼやける」効果が出ます。少なくとも2次精度以上のスキームを使わないと中心位置はずれがちです。また、高Re数（Re=10000）ではキャビティが定常解に落ち着かず、弱い振動解になることがある。収束判定を「残差が10⁻⁶」だけで行うと「振動しながら残差が下がっている」状態を正常終了と誤認することがあります。速度の時系列も確認することが大切です。
関連トピック
なった
詳しく
報告