電影工業中的流體模擬（四）- 納維斯托克斯方程（下）

01-26

上一節從動量守恆角度推導了NS方程的動量方程形式，由物質守恆，我們還可以得到連續方程形式：

考慮空間中的一塊區域 $Omega$ 中質量為 $M$ 密度為 $ho$ 的流體，我們有

$egin{equation}M = iiint_{Omega} ho,end{equation}$

並且動量 $mathbf{P} = Mmathbf{u}$ 滿足

$egin{equation}mathbf{P} = iiint_{Omega} ho mathbf{u}.end{equation}$

對這一塊區域 $Omega$ ，由於區域內部的流體是不變的，那麼流經該區域的流體的變化量可以用其邊緣（ $partial Omega$ ）流量的變化量來表示，也即

$egin{equation}frac{partial M}{partial t} = - iint_{partial Omega} ho mathbf{u} cdot mathbf{n},end{equation}$

這裡我們假定 $mathbf{n}$ 指向外部，所以有負號。

由高斯定理，上式又可以寫成

$egin{equation}frac{partial M}{partial t} = - iiint_{Omega} abla cdot ( ho mathbf{u} ).end{equation}$

注意到M的表達式，利用雷諾定理，我們可以得到

$egin{equation}frac{partial M}{partial t} = frac{partial}{partial t}iiint_{Omega} ho = iiint_{Omega}frac{partial ho}{partial t} = -iiint_{Omega} abla cdot ( ho mathbf{u}).end{equation}$

同樣地，上式對於任何區域 $Omega$ 均成立，所以

$egin{equation}frac{partial ho}{partial t} + abla cdot ( ho mathbf{u}) = 0.end{equation}$

這就是NS方程由物質守恆得到的連續方程。

上一節我們從動量的角度推導了描述流體運動的NS方程，其中我們引入了一個記號 D ，對於一個流體粒子的任何量（例如速度，位置） $mathbf{q}$ ，記作

$egin{equation}frac{Dmathbf{q}}{Dt} = frac{partial mathbf{q}}{partial t} + mathbf{u} cdot abla mathbf{q} .end{equation}$

這裡 $mathbf{u}$ 為粒子速度。這個導數我們把他稱之為物質導數（material derivative）。

引入物質導數是直觀的。考慮下面一個例子

考慮以恆定流速（例如上游有一以恆定速度泵水的水泵）流經上圖管口區域的流體，顯然如果我們隨著x軸移動，對於流體的任何一個粒子，該粒子的速度（相對管壁）都是越來越快的。但是由於前文所述，這是一個「流速恆定」的流體，那麼怎麼區分「恆定「和「越來越快」」呢？這就是物質導數的由來：我們「跟隨」流體粒子運動，流體粒子速度是越來越快的。但是如果我們考察空間中的任何一個定點，那麼流過該定點的粒子的速度都是不變的。前者就是物質導數（「越來越快」），後者就是一般的全導數（「恆定」）。

對於不可壓縮流體，顯然任何一個包含固定數量流體粒子的流體塊所佔空間體積都是不變的，也就是該流體任何處的物質導數滿足 $D ho = 0$ 。代入到連續方程，得到

$egin{equation} abla cdot mathbf{u} = 0,end{equation}$

這便是不可壓縮條件。

接下來我們考察流體的內部壓強 $p$ 需要滿足的關係。

對動量方程

$egin{equation}frac{D mathbf{u}}{D t} = mathbf{g} - frac{1}{ ho} abla p + u abla cdot abla mathbf{u}.end{equation}$

兩邊求散度，得到

$egin{equation} abla cdot frac{partial mathbf{u}}{partial t} + abla cdot (mathbf{u} cdot abla mathbf{u}) + abla cdot frac{1}{ ho} abla p= abla cdot (mathbf{g} + u abla cdot abla mathbf{u}).end{equation}$

由不可壓縮條件得到

$egin{equation} abla cdot frac{partial mathbf{u}}{partial t} = frac{partial}{partial t} abla cdot mathbf{u} = 0.end{equation}$

這樣，動量方程轉化成

$egin{equation} abla cdot frac{1}{ ho} abla p = abla cdot (- mathbf{u} cdot abla mathbf{u} + mathbf{g} + u abla cdot abla mathbf{u}),end{equation}$

流體模擬中一個很重要的過程就是求解滿足上式的 $p$ 。由於 $abla^2$ 離散 $p$ 以後可以得到對稱正定矩陣，考慮到邊界條件，上式的解是存在且唯一的。

後面的章節我們將陸續介紹流體模擬的演算法分類，並討論屬於拉氏(Lagrangian)方法的SPH（平滑粒子流體動力學）。

Reference:

Bridson, Robert. Fluid simulation for computer graphics. CRC Press, 2015.

Reynolds transport theorem. Wikipedia

\_(:3」∠)\_ \_(?ω?」∠)\_ \_(:з)∠)\_ ∠( ? 」∠)＿ \_(:зゝ∠)\_

請毫不猶豫地關注我們：

我們的網站：GraphiCon

知乎專欄：GraphiCon圖形控 - 知乎專欄

公眾號：GraphiCon

如果你有什麼想法，建議，或者想加入我們，你可以：

給我們發郵件：[email protected]

加入我們的QQ群：SIQGRAPH（342086343）

加入我們的slack群：GraphiCon

本作品採用知識共享署名-非商業性使用-禁止演繹 4.0 國際許可協議進行許可。