本周，我们讨论一下这一个语句：

fvm::laplacian(nu, U);

这一项是动量方程 $$ \frac{\partial U}{\partial t} + \nabla\cdot(UU) = -\nabla p + \nabla\cdot\tau $$ 的右边第二项。

考虑不可压缩牛顿流体，则有 $$ \tau = \nu ( \nabla U +(\nabla U)^T ) $$ 由于不可压缩牛顿流体有 $$ \nabla U = (\nabla U)^T,\, \nabla\cdot U =0 $$ 因此，$\tau$ 的散度为 $$ \nabla\cdot\tau = \nu \nabla^2 U $$ 这一步具体推导可见 wiki N-S 方程推导。

现在，我们主要来看一下这一项是如何离散的。把 $U$ 看成未知量 $\varphi$，设 $\nu=1$ $$ \iiint_V\nabla\cdot(\nabla \varphi) d V = \oint_{\partial V}\nabla \varphi d S $$ 可以看到，这一项是对 $\varphi$ 的梯度（是一向量）进行求散度，高斯积分后成为对 $\varphi$ 的梯度求面积分，离散后为 $$ \iiint_V\nabla\cdot(\nabla \varphi) d V = \sum \nabla \varphi d\vec{S} $$

设该项等于零，并以 $x$ 方向为例进行说明，网格间距均为 $\Delta x$，有 $$ \sum \nabla \varphi d\vec{S} = \left(\frac{d \varphi}{d x}\right)_A \vec{S}_A + \left(\frac{d \varphi}{d x}\right)_B \vec{S}_B=0 $$

考虑单个 cell，以 cell 4 进行解释，考虑 $S_A = S_B = S$

$$ \sum \nabla \varphi d\vec{S} = S\left(\frac{\varphi_5-\varphi_4}{\Delta x} - \frac{\varphi_4 - \varphi_3}{\Delta x}\right) \ = \frac{S}{\Delta x} ( -2\varphi_4 +\varphi_5 + \varphi_3)=0 $$

设左右两面边界上的值分别为 $\varphi^l_0$、$\varphi^r_0$，按此方法，对其余网格进行离散，形成矩阵乘法，有 $$ \left[ \begin{array}% -3 & 1 & 0 & 0\ 1 & -2 & 1 & 0 \ 0 & 1 & -2 & 1\ 0 & 0 & 1 & -3\ \end{array} \right] \times \left[ \begin{array}% \varphi_1\\ \varphi_2\\ \varphi_3\\ \varphi_4\\ \end{array}

\right]

\left[ \begin{array}% -2\varphi^l_0\\ 0\\ 0\\ -2\varphi^r_0\\ \end{array} \right] $$ 上式与动量方程的其他项相加，就形成动量方程的离散形式。

CFDwired Forum 社区讨论版：code of the week #3: laplacian(nu,U)