测试----

[李东岳]

不只是N-S方程

谈到CFD，大部分人最开始想到的就是Navier-Stokes方程（N-S方程）。N-S方程从守恒定律推导而来：
$$\begin{array}{}\text{(1)}& \frac{\mathrm{\partial }\rho }{\mathrm{\partial }t}+\mathrm{\nabla }\cdot (\rho \mathbf{U})=0,\end{array}$$
$$\begin{array}{}\text{(2)}& \frac{\mathrm{\partial }\rho \mathbf{U}}{\mathrm{\partial }t}+\mathrm{\nabla }\cdot (\rho \mathbf{U}\mathbf{U})=-\mathrm{\nabla }p+\mathrm{\nabla }\cdot \tau ,\end{array}$$
其中$\rho$为密度，$\mathbf{U}$为速度，$p$为压力，$\tau$为剪切应力。N-S方程具有以下特点：

方程$\text{(2)}$中左边第二项是关于$\mathbf{U}$的偏导数，这种未知量和未知量乘积的问题构成非线性问题，CFD对非线性问题需要特殊处理。另一方面，非线性的双曲问题的解可能会存在间断（如激波）。激波通常存在于高超声速的欧拉问题求解中。同时，非线性项也是湍流在数学方程中的体现；

方程$\text{(2)}$的数学特征为抛物线。不同数学特征的问题需要调用不同的时间/空间离散格式，隐性时间格式更有利于求解抛物线问题。若方程$\text{(2)}$中省略若干项则会改变方程的数学特征，例如若将方程右侧置为$0$，则变为双曲特征的欧拉方程。欧拉方程得益于其双曲特性，可采用迎风类显性算法推进，各种基于有限体积法的高分辨率格式因此而生（交错网格中心格式、中心-迎风格式等）。同时，动量方程的对流项和扩散项的相对强弱，也会影响边界条件的设置（如DNS的出口无反射边界条件）；

在马赫数较大时（如大于0.3），方程$\text{(1)}$可用来求解密度，方程$\text{(2)}$可用来求解速度，同时附加能量方程求解温度以及状态方程求解压力，即密度基求解器。在马赫数较小时，并没有单独的压力方程，并且方程$\text{(1)}$缺少主要求解变量。这导致压力的求解需要特殊的策略。这也是CFD中压力基SIMPLE/PISO算法要处理的问题。在构建离散矩阵的情况下，对速度和压力整合处理还是单独处理，是分离求解与耦合算法的重要问题；

N-S方程之前传：N-S方程为宏观方程，调用了宏观假定，其可从玻尔兹曼方程附加Chapman-Enskog展开推导而来。N-S方程也可以看做从介尺度模型演化的宏观二阶矩模型。在无压力无粘性的条件下具备弱双曲特征。由于失去了高阶矩的统计学特征，因此N-S方程在某些情况下是不适用的，比如稀薄流、微流动、以及气固多相流中的颗粒轨迹交叉现象都不能用N-S方程来描述。其数学本源在于N-S方程假定流体为一种近似均衡状态，在某一个网格点存在一个单值的分布。然而，在分子/颗粒数量非常小的情况下，这些分子/颗粒基本无碰撞，会导致多值问题的产生。因此，在这种情况下，需要对统计学分布的高阶矩进行追踪，N-S方程作为二阶矩模型是往往不够的。

正如上面讨论的，N-S方程中存在大量的数学问题，初学者也看的不明不白。这无关紧要。本章仅仅抛出一块砖头，让大家了解隐藏在CFD中的方程的美妙含义。这些方程中充满火星符号（如$\mathrm{\nabla }$，$\mathrm{\partial }$），也让人摸不清头脑。\href{https://www.cfd-china.com/topic/2397}{东岳流体CFD课程}因此分为两步，第一步是基本的N-S方程入门，第二步是N-S方程求解。第一步需要同学们通过本笔记进行预习。第二步将在课程上讲授。

学习CFD理论首先遇到的就是各种各样的偏微分算符。能看懂、拆分CFD方程是研究算法的最基本步骤。给出一个CFD方程（如动量方程），可以不知道是怎么推出来的，但要能看懂。本章介绍CFD方程的各种写法，不涉及到任何的CFD算法。

矢量标识法

本笔记中主要采用矢量标识法讨论CFD方程，在矢量标识法中，标量全部采用斜体，如压力$p$。矢量采用正体加粗，如速度矢量$\mathbf{U}$，其具有三个分量$u_1$，$u_2$，$u_3$或$u$，$v$，$w$。二阶张量也采用正体\footnote{如果难以理解二阶张量的含义，可以这样尝试：矢量是一阶张量，具有三个分量，二阶张量则具有9个分量。}，比如应力张量$\mathit{\tau }$，其具备9个分量，其可表示为：
$$\begin{array}{}\text{(1)}& \mathit{\tau }=\left[\begin{array}{ccccccc}{\tau }_{11}& {\tau }_{12}& {\tau }_{13}{\tau }_{21}& {\tau }_{22}& {\tau }_{23}{\tau }_{31}& {\tau }_{32}& {\tau }_{33}\end{array}\right]\end{array}$$
下面粘度$\nu$为$1$的不可压缩流体动量方程为例对其进行拆分，这个方程若采用矢量标识法可以写为：
$$\begin{array}{}\text{(2)}& \frac{\mathrm{\partial }\mathbf{U}}{\mathrm{\partial }t}+\mathrm{\nabla }\cdot (\mathbf{U}\mathbf{U})=-\mathrm{\nabla }\frac{p}{\rho }+\mathrm{\nabla }\cdot (\mathrm{\nabla }\mathbf{U})\end{array}$$
其中的$\mathbf{U}$为速度矢量，$p$为压力，$\rho$为密度。一般来讲，CFD文献中通常采用方程$\text{(2)}$的形式，而并不进行展开进而更加紧凑。下面介绍如何将其展开为3个方程。
\begin{itemize}
\item 方程$\text{(2)}$第一项表示$\mathbf{U}$对时间的偏导数，因为$\mathbf{U}$为矢量，故其导数的分量形式为：
$$\begin{array}{}\text{(3)}& \frac{\mathrm{\partial }\mathbf{U}}{\mathrm{\partial }t}=\left[\begin{array}{c}\frac{\mathrm{\partial }{u}_1}{\mathrm{\partial }t}\frac{\mathrm{\partial }{u}_2}{\mathrm{\partial }t}\frac{\mathrm{\partial }{u}_3}{\mathrm{\partial }t}\end{array}\right]\end{array}$$
其中$u_1$表示$x$方向速度，$u_2$表示$y$方向速度，$u_3$表示$z$方向速度。这样拆分之后的方程，即为各个方向的速度针对时间的偏导数。如果理解方程$\text{(3)}$有困难，那么有必要预习一下《高等数学》（同济大学版）第二章。
\item 方程$\text{(2)}$第二项$\mathrm{\nabla }\cdot (\mathbf{U}\mathbf{U})$中的$\mathbf{U}\mathbf{U}$是一种简写，完整形式为$\mathbf{U}\otimes \mathbf{U}$，$\otimes$是一个张量运算符。依据$\otimes$的定义，$\mathbf{U}\mathbf{U}$可以写为：
$$\begin{array}{}\text{(4)}& \mathbf{U}\otimes \mathbf{U}=\mathbf{U}\mathbf{U}=\mathbf{U}\cdot {\mathbf{U}}^T=\left[\begin{array}{c}{u}_1{u}_2{u}_3\end{array}\right][u_1,u_2,u_3]=\left[\begin{array}{ccccccc}{u}_1{u}_1& u_1{u}_2& u_1{u}_3{u}_2{u}_1& u_2{u}_2& u_2{u}_3{u}_3{u}_1& u_3{u}_2& u_3{u}_3\end{array}\right]\end{array}$$
接下来看 $\mathrm{\nabla }\cdot$，其为散度算符，有时用$\mathrm{d}\mathrm{i}\mathrm{v}$来表示。对一个矢量（$1$阶张量）做散度的结果为一个标量（$0$阶张量），对一个$2$阶张量做散度的结果为矢量（$1$阶张量）。因此，对任意$n$阶张量做散度操作之后，结果为$n-1$阶张量。举例，对一个矢量$\mathbf{U}$做散度有：$\mathrm{\nabla }\cdot \mathbf{U}=\frac{\mathrm{\partial }{u}_1}{\mathrm{\partial }x}+\frac{\mathrm{\partial }{u}_2}{\mathrm{\partial }y}+\frac{\mathrm{\partial }{u}_3}{\mathrm{\partial }z}$。因此，方程$\text{(2)}$中的第二项$\mathrm{\nabla }\cdot (\mathbf{U}\mathbf{U})$即为对一个$2$阶张量做散度：
$$\begin{array}{}\text{(5)}& \mathrm{\nabla }\cdot \left(\mathbf{U}\mathbf{U}\right)=\mathrm{\nabla }\cdot \left[\begin{array}{ccccccc}{u}_1{u}_1& u_1{u}_2& u_1{u}_3{u}_2{u}_1& u_2{u}_2& u_2{u}_3{u}_3{u}_1& u_3{u}_2& u_3{u}_3\end{array}\right]=\left[\begin{array}{c}\frac{\mathrm{\partial }{u}_1{u}_1}{\mathrm{\partial }x}+\frac{\mathrm{\partial }{u}_2{u}_1}{\mathrm{\partial }y}+\frac{\mathrm{\partial }{u}_3{u}_1}{\mathrm{\partial }z}\frac{\mathrm{\partial }{u}_1{u}_2}{\mathrm{\partial }x}+\frac{\mathrm{\partial }{u}_2{u}_2}{\mathrm{\partial }y}+\frac{\mathrm{\partial }{u}_3{u}_2}{\mathrm{\partial }z}\frac{\mathrm{\partial }{u}_1{u}_3}{\mathrm{\partial }x}+\frac{\mathrm{\partial }{u}_2{u}_3}{\mathrm{\partial }y}+\frac{\mathrm{\partial }{u}_3{u}_3}{\mathrm{\partial }z}\end{array}\right]\end{array}$$

[warnerchang]

@李东岳 东岳老师，方程双曲性质具体是什么？不是特别懂，有时是双曲，有时又是抛物线，椭圆的，他们之间是什么关系呢？

其次，为啥式（2）中的运动粘度为1呢？式（2）是openfoam中求解方程的形式吗？

[李东岳]

@warnerchang 方程2是简写，省略粘度符号，要不看起来太麻烦

[warnerchang]

@李东岳 噢噢，好的