こんにちは、JS2IIUです。Pythonで流体シミュレーションを攻略する連載の第3回です。
第1回では質量の保存（連続の式）を、第2回では運動量とエネルギーの保存を学びました。これらはいわば流体の「憲法」であり、あらゆる流体現象が従わなければならない絶対的なルールです。
しかし、現実の流体は非常に複雑です。空気や水の粘り気（粘性）、渦の発生、複雑な温度変化などが絡み合い、そのままでは数式を解くことも、直感的に理解することも困難です。そこで物理学の世界では、まず「理想的な状態」を想定して、本質的な構造を明らかにしようとします。
今回は、粘性がまったくない仮想的な流体である「完全流体（理想流体）」の世界へと足を踏み入れます。今回もよろしくお願いします。
1. 完全流体とオイラーの運動方程式
「完全流体」とは、以下の2つの理想的な性質を持つ流体を指します。
- 粘性がない（非粘性）: 流体内部に摩擦（粘性力）が全く存在せず、流体粒子同士が滑らかにすり抜ける。
- 熱伝導がない: 流体粒子間での熱のやり取り（熱伝導）が無視できる。
このような仮定を置くことで、現実の流体の複雑さ（摩擦や熱の拡散）を排除し、運動の本質的な構造を明らかにできます。
オイラーの運動方程式の導出と意味
第2回で学んだ運動量保存則（ニュートンの運動方程式）を、粘性ゼロ・熱伝導ゼロの完全流体に適用すると、オイラーの運動方程式が得られます。
$$
\frac{\partial \mathbf{u}}{\partial t} + (\mathbf{u} \cdot \nabla)\mathbf{u} = -\frac{1}{\rho}\nabla P + \mathbf{g}
$$
ここで：
- \(\mathbf{u}\)：速度ベクトル場（\(u, v, w\)）
- \(t\)：時間
- \(\rho\)：流体の密度
- \(P\)：圧力
- \(\mathbf{g}\)：重力加速度ベクトル（外力）
各項の物理的意味
- \(\frac{\partial \mathbf{u}}{\partial t}\)：非定常項（局所加速度）…その場での速度の時間変化
- \((\mathbf{u} \cdot \nabla)\mathbf{u}\)：移流項（対流項、慣性項）…流体粒子が空間的に移動することで生じる速度変化（流体自身が自分の速度を運ぶ）
- \(-\frac{1}{\rho}\nabla P\)：圧力勾配項…圧力の空間変化による加速度
- \(\mathbf{g}\)：外力項…重力などの体積力
移流項 \((\mathbf{u} \cdot \nabla)\mathbf{u}\) の詳細
この項は、流体力学の非線形性の本質を担っています。\((\mathbf{u} \cdot \nabla)\mathbf{u}\) は「流体粒子が空間を移動しながら、その移動先での速度場の変化を受ける」ことを意味します。これは、
- 速度場が空間的に一様ならゼロになる
- 速度場が空間的に変化していれば、粒子はその変化を“運ぶ”
という性質を持ちます。
この移流項の存在により、完全流体であっても流れは非常にダイナミックかつ複雑な挙動（渦や波、非線形現象）を示します。
オイラー方程式の成り立ちと仮定
オイラー方程式は、以下の仮定のもとで成り立ちます：
- 流体は連続体として扱う（分子レベルの離散性は無視）
- 粘性ゼロ（非粘性）
- 熱伝導ゼロ
- 圧縮性の有無は問わない（非圧縮性の場合は \(\nabla \cdot \mathbf{u} = 0\) を追加）
- 外力は体積力（重力など）のみ
オイラー方程式の意義
オイラー方程式は、粘性や熱伝導といった現実の複雑さを排除しつつも、流体の本質的な非線形ダイナミクスを捉えることができます。特に、移流項の非線形性が流体現象の多様性（渦、波、乱流の萌芽など）を生み出します。
しかし、オイラー方程式自体も一般には解析的に解くのが難しいため、さらに「渦がない」などの条件を加えることで、よりシンプルで美しい数学的構造（ポテンシャル流）が現れます。
2. ポテンシャル流れ：渦のない世界
流体の動きの中に「渦」がない状態（非回転：\(\nabla \times \mathbf{u} = 0\)）を想定してみましょう。
数学のベクトル解析には、「回転がゼロのベクトル場は、あるスカラ関数の勾配で表現できる」という便利な定理があります。このスカラ関数 \(\phi\)（ファイ）を速度ポテンシャルと呼びます。
$$
\mathbf{u} = \nabla \phi
$$
つまり、速度ベクトル \(\mathbf{u} = (u, v, w)\) をわざわざ3つの変数として追いかけなくても、たった1つのスカラ関数 \(\phi\) さえ分かれば、流体の動きがすべて記述できてしまうのです。これは計算負荷を劇的に減らす強力な武器になります。
ラプラス方程式の登場
ここで、第1回で学んだ「非圧縮性流体の連続の式（\(\nabla \cdot \mathbf{u} = 0\)）」を思い出してください。速度ベクトル \(\mathbf{u}\) に \(\nabla \phi\) を代入すると、以下の式が得られます。
$$
\nabla \cdot (\nabla \phi) = \nabla^2 \phi = 0
$$
これが有名なラプラス方程式です。
驚くべきことに、渦のない理想的な流れでは、速度ポテンシャルはラプラス方程式を満たすことがわかりました。ラプラス方程式は「滑らかな分布」を作る性質を持っており、重力ポテンシャルや静電ポテンシャルなど、物理学の至る所に現れる非常に美しい方程式です。
流体の運動をさらに単純化するため、「渦（vorticity）」が全く存在しない状態、すなわち非回転流（irrotational flow）を考えます。これは数学的には次の条件で表されます：
$$
\nabla \times \mathbf{u} = 0
$$
この条件は、流体粒子が微小なスケールで自転しない（回転成分を持たない）ことを意味します。現実の流体では渦が発生することが多いですが、理想化することで理論解析が大幅に容易になります。
ベクトル場のポテンシャル表現
ベクトル解析の定理（ポアンカレの補題）より、回転がゼロのベクトル場は必ず「あるスカラ関数の勾配」として表現できます。このスカラ関数 (\phi) を速度ポテンシャル（potential function）と呼びます。
$$
\mathbf{u} = \nabla \phi
$$
この式は、速度ベクトル \(\mathbf{u} = (u, v, w)\) を3つの成分で追いかける代わりに、たった1つのスカラ場 (\phi) だけで流体の運動を完全に記述できることを意味します。これは計算量や理論解析の観点で非常に大きな利点です。
また、速度ポテンシャルは物理的には「流体がどれだけ“流れたか”」の累積量のような意味合いも持ちます。
ラプラス方程式の導出と物理的意味
ここで、第1回で学んだ「非圧縮性流体の連続の式」
$$
\nabla \cdot \mathbf{u} = 0
$$
を思い出してください。\(\mathbf{u} = \nabla \phi\) を代入すると、
$$
\nabla \cdot (\nabla \phi) = \nabla^2 \phi = 0
$$
となります。
このラプラス方程式は、物理学・工学のさまざまな分野で現れる「滑らかさ」を保証する方程式です。ポテンシャル流では、速度ポテンシャルがこの方程式を満たすため、流れ場は非常に滑らかで連続的な分布となります。
ポテンシャル流の特徴と応用
- 線形性: ラプラス方程式は線形方程式なので、複数の基本解（点ソース、一様流、ダブルレットなど）を足し合わせて複雑な流れを構成できる（重ね合わせの原理）。
- 渦なし・摩擦なし: ポテンシャル流は渦も摩擦もない理想的な流れであり、現実の流体の一部現象（翼の揚力、流線形物体まわりの流れなど）の近似に使われる。
- 他分野との共通性: ラプラス方程式は重力場・静電場・熱伝導など多くの物理現象で現れるため、数学的な解析手法や直感が流体にも応用できる。
このように、ポテンシャル流の理論は流体力学の中でも特に美しく、解析的な解が得られる数少ないケースの一つです。
3. Pythonによる実装：速度ポテンシャルの視覚化
それでは、Pythonを使ってポテンシャル流れを計算してみましょう。
今回は、2つの「湧き出し（ソース）」と「吸い込み（シンク）」を配置したときに、どのような流れが生まれるかをシミュレーションします。
ラプラス方程式の解は、線形結合（足し算）ができるという素晴らしい性質を持っています。点ソースのポテンシャルを足し合わせるだけで、複雑な流れを簡単に作ることができるのです。
import torch
import matplotlib.pyplot as plt
def get_potential_source(x, y, x0, y0, strength):
"""
点ソース（湧き出し）による速度ポテンシャルを計算する。
phi = strength / (2 * pi) * ln(r)
"""
# 中心からの距離 r を計算 (数値安定のため小さな値 eps を加える)
r = torch.sqrt((x - x0)**2 + (y - y0)**2 + 1e-9)
return (strength / (2 * torch.pi)) * torch.log(r)
# --- パラメータ設定 ---
nx, ny = 100, 100
lx, ly = 2.0, 2.0
x = torch.linspace(-lx/2, lx/2, nx)
y = torch.linspace(-ly/2, ly/2, ny)
Y, X = torch.meshgrid(y, x, indexing='ij')
# --- 速度ポテンシャルの合成 ---
# (0.5, 0) に湧き出し(Source)、(-0.5, 0) に吸い込み(Sink)を配置
phi_source = get_potential_source(X, Y, x0=0.5, y0=0, strength=5.0)
phi_sink = get_potential_source(X, Y, x0=-0.5, y0=0, strength=-5.0)
# 重ね合わせの原理：単純に足し合わせる
phi = phi_source + phi_sink
# --- 速度場の計算 ---
# u = d_phi / dx, v = d_phi / dy
# PyTorchのスライシングを用いて中央差分で勾配を計算
dx = x[1] - x[0]
dy = y[1] - y[0]
u = torch.zeros_like(phi)
v = torch.zeros_like(phi)
# 内部点の勾配計算
u[:, 1:-1] = (phi[:, 2:] - phi[:, :-2]) / (2 * dx)
v[1:-1, :] = (phi[2:, :] - phi[:-2, :]) / (2 * dy)
# --- 可視化 ---
plt.figure(figsize=(10, 8))
# 速度ポテンシャルの等高線を表示
cp = plt.contour(X, Y, phi, levels=30, cmap='viridis', alpha=0.5)
plt.clabel(cp, inline=True, fontsize=8)
# 流線（Streamlines）を表示
# 流線は速度ベクトルに沿った線であり、ポテンシャルの等高線と直交する性質がある
plt.streamplot(x.numpy(), y.numpy(), u.numpy(), v.numpy(), color='black', linewidth=1, density=1.5)
# ソースとシンクの地点をプロット
plt.plot(0.5, 0, 'ro', label='Source (湧き出し)')
plt.plot(-0.5, 0, 'bo', label='Sink (吸い込み)')
plt.title("Potential Flow: Source and Sink Pair")
plt.xlabel("x")
plt.ylabel("y")
plt.legend()
plt.axis('equal')
plt.grid(alpha=0.3)
plt.show()
コードの解説
- 重ね合わせの原理:
phi = phi_source + phi_sinkというコードが、ポテンシャル流の最大の特徴を示しています。複雑な微分方程式を一から解かなくても、基本的な「部品（ソース、シンク、一様流など）」の解を足すだけで、新しい流れを構築できるのです。 - 速度の導出:
速度ポテンシャルphiを空間微分することで速度ベクトルu, vを得ています。これは物理的には「ポテンシャルの高い方から低い方へ流れる」という勾配の考え方に基づいています。 - 等高線と流線の関係:
グラフを見ると、速度ポテンシャルの等高線（カラーの線）と、黒い流線が常に直交していることに気づくはずです。これは地図の等高線（ポテンシャル）と、雨水が流れる方向（速度）の関係と同じです。
4. 理想流体の限界：ダランベールのパラドックス
ポテンシャル流の理論は非常に美しく、翼の揚力の計算や流線形物体まわりの流れ解析など、現代でも多くの応用があります。しかし、この「理想的な世界」には決定的な限界が存在します。それがダランベールのパラドックス（D’Alembert’s paradox）です。
ダランベールのパラドックスとは
18世紀にジャン・ル・ロン・ダランベールによって指摘されたこのパラドックスは、次のような現象です：
「完全流体（非粘性・非回転）の中に物体を置いた場合、物体には抗力（ドラッグ、空気抵抗）が全く働かない」
これはポテンシャル流の理論に従い、物体の前後で圧力分布を厳密に計算すると、前方で受ける圧力と後方で受ける圧力が完全に釣り合い、合力としての抗力がゼロになるためです。
数学的背景
ポテンシャル流では、粘性がないため流体粒子は物体表面に沿って滑らかに流れ、物体の後方で流れが分離したり渦が発生したりしません。そのため、物体の前後で圧力分布が対称となり、抗力が消失します。
現実との乖離とその理由
現実の流体では、どんなに粘性が小さくても、物体の後方で「剥離（separation）」や「渦（wake）」が発生し、これが大きな抗力の原因となります。これは、
- 物体表面近傍の境界層（boundary layer）で粘性の影響が無視できない
- 境界層が剥離し、後流（wake）に渦が発生する
といった現象によるものです。ポテンシャル流ではこれらの現象が一切表現できません。
パラドックスの意義と流体力学への影響
このパラドックスは、理論と現実のギャップを明確に示し、20世紀初頭の境界層理論（プラントル）やナビエ・ストークス方程式の発展につながりました。現実の流体抵抗や渦の発生を正しく記述するには、粘性の効果を数式に組み込む必要があるのです。
このように、理想流体の理論は美しい一方で、現実の流体現象を理解するためには「粘性」という“現実の複雑さ”を避けて通れないことを、ダランベールのパラドックスは教えてくれます。
5. まとめと次回の予告
第3回では、流体力学における「理想的な模型」について学びました。
- 粘性や熱伝導を無視した完全流体と、その運動方程式であるオイラー方程式。
- 渦がない場合に定義できる速度ポテンシャル。
- 速度ポテンシャルが満たすラプラス方程式と、解の重ね合わせの原理。
- Pythonによるソース・シンク流の視覚化。
理想的な世界を理解したことで、逆に「現実の流体の何が難しいのか」が明確になりました。
次回、第4回では、いよいよ流体力学の王様とも言える「ナビエ・ストークス方程式」が登場します。完全流体で切り捨てた「粘性」を数式に組み戻し、現実の流体が持つドロドロとした、しかし深みのある性質をどう記述するのかを解説します。
最後まで読んでいただきありがとうございます。