标签：tau 方程 partial 推导 equation 体积 rho delta frac

3 有限体积法：推导方程

基本原理和目标

（注意：这一节看不懂没关系，在后面的推导中会慢慢用到）

• 质量、动量和能量的守恒

  1. 流体的质量守恒
  2. 动量改变的速度 = 一个流体粒子上受到的力的总和（牛顿第二定律）
  3. 能量改变的速度 = 一个流体粒子吸收的热量，和作用在其上的功的总和（热力学第一定律）

• 推导出控制流动的偏微分方程

• 通过牛顿模型描述粘性应力，推导出纳维-斯托克斯方程

• 控制流体行为的方程和传输方程是相似的

• 把传输方程的积分形式应用于有限时间间隔和有限控制体积

• 把物理现象分为椭圆型、抛物线型、双曲线型

• 将流体看作一个连续体（continuum）

• 考虑一个小流体元素，其边长为\(\delta x\)，\(\delta y\)，\(\delta z\)​

三维物理量守恒

质量守恒

一个流体元素中质量的增加速度 = 流入流体元素的净流动速度

质量增加速度：

\[\begin{align}\frac{\partial}{\partial t}(\rho\delta x\delta y\delta z)=\frac{\partial\rho}{\partial t}(\delta x\delta y\delta z)\end{align} \]

注意：这里\(\rho\delta x\delta y\delta z\)可理解为\(\rho \delta V\)，也就是代表了该流体微元的总质量。

然后求质量关于时间的一阶导，就是质量的净增加速度了。

穿过流体元素边界，流入流体元素的净质量流动速率为：

​​

其中，约定流入流体微元的质量流动标正号，流出的标负号。

注意这里是如何运用泰勒级数展开的（以x方向，即图中横着的方向为例）：

1. 在点x处的质量流动速率密度（即单位面积的质量流率）为\(\rho u\)；

2. 使用一阶泰勒级数展开：

   \[\begin{align}f(x+\delta x)=f(x)+\frac{df(x)}{dx} \delta x\end{align} \]

3. 在\(x-\frac{1}{2} \delta x\)（即流入处）：

   \[\begin{align} f(x-\frac{1}{2}\delta x)&=f(x)-\frac{\partial f(x)}{\partial x}\cdot \frac{1}{2}\delta{x}\\ &=\rho u-\frac{\partial (\rho u)}{\partial x}\cdot \frac{1}{2}\delta{x} \end{align} \]

   在\(x+\frac{1}{2} \delta x\)（即流出处）同理，此处不再展示。

4. x方向流入的真正质量流动速率为（即乘上个截面面积）：

   \[\begin{align}\left(\rho u-\frac{\partial(\rho u)}{\partial x}\frac12\delta x\right)\delta y\delta z\end{align} \]

把六个方向（x，y，z方向的流入和流出）加一下，可得：

\[\begin{align} \begin{gathered} \left(\rho u-\frac{\partial(\rho u)}{\partial x}\frac12\delta x\right)\delta y\delta z-\left(\rho u+\frac{\partial(\rho u)}{\partial x}\frac12\delta x\right)\delta y\delta z \\ +\left(\rho v-\frac{\partial(\rho v)}{\partial y}\frac{1}{2}\delta y\right)\delta x\delta z-\left(\rho v+\frac{\partial(\rho v)}{\partial y}\frac{1}{2}\delta y\right)\delta x\delta z \\ +\left(\rho w-\frac{\partial(\rho w)}{\partial z}\frac{1}{2}\delta z\right)\delta y\delta x-\left(\rho w+\frac{\partial(\rho w)}{\partial z}\frac{1}{2}\delta z\right)\delta y\delta x \end{gathered}\end{align} \]

化简可得最终穿过流体元素边界，流入流体元素的净质量流动速率：

\[\begin{align}\left(\frac{\partial(\rho u)}{\partial x}+\frac{\partial(\rho v)}{\partial y}+\frac{\partial(\rho w)}{\partial z}\right)\delta x\delta y\delta z\end{align} \]

那么，我们可以得到质量守恒方程：

\[\begin{align} \frac{\partial\rho}{\partial t}(\delta x\delta y\delta z)+\left(\frac{\partial(\rho u)}{\partial x}+\frac{\partial(\rho v)}{\partial y}+\frac{\partial(\rho w)}{\partial z}\right)\delta x\delta y\delta z=0\end{align} \]

把\(\delta x \delta y \delta z\)消掉，可得：

\[\begin{align}&\frac{\partial\rho}{\partial t}+\frac{\partial(\rho u)}{\partial x}+\frac{\partial(\rho v)}{\partial y}+\frac{\partial(\rho w)}{\partial z}=0 \\&\frac{\partial\rho}{\partial t}+\nabla\cdot(\rho\boldsymbol{u})=0\end{align} \]

注意，这里的Nabla算子：

\[\begin{align} \nabla=\frac{\partial}{\partial x}\vec{i}+\frac{\partial}{\partial y}\vec{j}+\frac{\partial}{\partial z}\vec{k} \end{align} \]

上述方程可表示在可压缩流体中的某一点，非稳态、三维流动下的质量守恒或连续性方程。

而对于不可压缩流动，密度是守恒的，此时：

\[\begin{align} &\nabla\cdot\boldsymbol{u}=0\\ \quad&\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}+\frac{\partial w}{\partial z}=0 \end{align} \]

推广：物理量守恒

现设有一物理量在每单位质量的值为\(\phi\)。

沿着一个流体粒子的轨迹，对函数\(\phi\)的总（total）或实质（substantive）导数可表示为：

\[\begin{align} \frac{D\phi}{Dt}=\frac{\partial\phi}{\partial t}+\frac{\partial\phi}{\partial x}\frac{dx}{dt}+\frac{\partial\phi}{\partial y}\frac{dy}{dt}+\frac{\partial\phi}{\partial z}\frac{dz}{dt}\end{align} \]

注意，流体粒子会随着流动而移动，故：

\[\begin{align} \frac{dx}{dt}=u,~\frac{dy}{dt}=v,~\frac{dz}{dt}=w \end{align} \]

由此可得，\(\phi\)的实质导数可写作：

\[\begin{align} \frac{D\phi}{Dt}&=\frac{\partial\phi}{\partial t}+u\frac{\partial\phi}{\partial x}+v\frac{\partial\phi}{\partial y}+w\frac{\partial\phi}{\partial z}\\ &=\frac{\partial\phi}{\partial t}+\boldsymbol{u}\cdot\nabla\phi \end{align} \]

该方程表示，该物理量随着时间的变化+流体粒子在空间中的运动与物理量在空间中变化的叠加效应=单位质量物质的\(\phi\)属性沿着一个流体微元运动轨迹的总变化率。

那么，对于一个流体粒子而言，其物理量\(\phi\)每单位体积的变化速度为：

\[\begin{align} \rho\frac{D\phi}{Dt}=\rho\left(\frac{\partial\phi}{\partial t}+\boldsymbol{u}\cdot\nabla\phi\right) \end{align} \]

下一步，我们向任意守恒的物理量（arbitrary conserved property）推广（generalize）出其守恒方程：

（注意：别直接用这个方程，这个是来自PPT的，请参见下面的方程27）

\[\begin{align} \frac{\partial(\rho\phi)}{\partial t}+\nabla\cdot(\rho\phi\boldsymbol{u})=\rho\left[\frac{\partial\phi}{\partial t}+\boldsymbol{u}\cdot\nabla\phi\right]+\phi\left[\frac{\partial\rho}{\partial t}+\nabla\cdot(\rho\boldsymbol{u})\right]=\rho\frac{D\phi}{Dt} \end{align} \]

其具体推导过程为：

1. 原实质导数参见方程17：\(\frac{D\phi}{Dt}=\frac{\partial\phi}{\partial t}+\boldsymbol{u}\cdot\nabla\phi\)​

2. 首先，推广到密度和物理量一起，即考虑\(\rho \phi\)随时间的变化率，此时使用乘积法则：

   \[\begin{align} d(uv)=udv+vdu \end{align} \]

   利用该法则，有：

   \[\begin{equation} \frac{\partial(\rho\phi)}{\partial t}=\phi\frac{\partial\rho}{\partial t}+\rho\frac{\partial\phi}{\partial t} \end{equation} \]

3. 下一步，考虑\(\rho \phi\)在空间上的变化，加入流体运动速度\(\boldsymbol{u}\)（这里加粗代表x，y，z三个方向的速度矢量），即对\(\rho \phi \boldsymbol{u}\)取散度，三个物理量都可能随着位置变化：

   \[\begin{equation} \nabla\cdot(\rho\phi\vec{u})=\phi[\nabla\cdot(\rho\vec{u})]+(\rho\vec{u})\cdot\nabla\phi \end{equation} \]

4. 将过程2和过程3的结果，即时间导数和空间上的散度项结合起来：

   \[\begin{align} \frac{\partial(\rho\phi)}{\partial t}+\nabla\cdot(\rho\phi\vec{u})=\left(\phi\frac{\partial\rho}{\partial t}+\rho\frac{\partial\phi}{\partial t}\right)+(\phi(\nabla\cdot(\rho\vec{u}))+\rho\vec{u}\cdot\nabla\phi) \end{align} \]

5. 根据连续性方程，质量守恒要求：

   \[\begin{equation} \frac{\partial\rho}{\partial t}+\nabla\cdot(\rho\vec{u})=0 \end{equation} \]

   注意到方程23中：

   \[\begin{equation} \phi \frac{\partial\rho}{\partial t}+\phi \cdot \nabla\cdot(\rho\vec{u})=0 \end{equation} \]

   故结合后的方程可写作：

   \[\begin{equation} \frac{\partial(\rho\phi)}{\partial t}+\nabla\cdot(\rho\phi\vec{u})=\rho\frac{\partial\phi}{\partial t}+\rho\vec{u}\cdot\nabla\phi \end{equation} \]

6. 注意到方程26中右侧项恰等于\(\rho\)乘以\(\frac{D\phi}{Dt}\)。故可将推广后的守恒方程写作：

   \[\begin{equation} \frac{\partial(\rho\phi)}{\partial t}+\nabla\cdot(\rho\phi\vec{u})=\rho\left(\frac{\partial\phi}{\partial t}+\vec{u}\cdot\nabla\phi\right)=\rho\frac{D\phi}{Dt} \end{equation} \]

   其物理意义为：

   \[物理量\phi随时间的变化率+净流出速率=总变化率 \]

动量守恒

牛顿第二定律表明，一个流体粒子动量的改变=作用在流体粒子上力的总和。

对于一个流体粒子，在x, y, z方向，每单位体积动量增加的速度为：

\[\begin{equation} \rho\frac{Du}{Dt},~\rho\frac{Dv}{Dt},~\rho\frac{Dw}{Dt} \end{equation} \]

注意这里：动量\(P=mv\)（这里v不代表特定于y方向的速度，而是泛指速度），每单位体积动量即为\(\rho v\)，则每单位体积动量增加的速度即为密度\(\times\)速度关于时间的导数。

作用于流体元上的力：

• 表面力：作用于流体表面，由于流体元与其周围流体或边界的接触产生。

  1. 压力力（pressure forces）

     • 由于流体中的压力梯度作用于流体元表面产生的力；
     • 垂直于流体元表面，指向外界（不过这里好像也不一定，要看具体流动是出微元还是进）；
     • 导致流体加减速，是推动流体流动的重要因素。

  2. 黏性力（viscous forces）

     • 由于流体的黏性（流体分子之间的内摩擦），在流体元表面产生的剪切应力；
     • 在流体层之间起阻力作用，特别是在层流和接触边界处；
     • 造成速度梯度（摩擦-相对速度），使得流体从高速区域向低速区域传递动量，最终影响整个流体的流动特性；
     • 补充：摩擦力是分子间的吸引力（如范德华力）、原子间的电磁排斥力、原子间电子交换或共享形成弱键合（如氢键）等等电磁力在宏观下的体现。

• 体力（body forces，不是HP）：外界力场直接作用于流体的每个体积单元内部。

  1. 重力，影响流体压力分布、密度梯度等。
  2. 离心力（centrifugal force）：在旋转系统（如地球自转引起的大气、海洋环流）中，由于旋转参考系下的惯性效应产生的力，指向远离旋转轴的方向。
  3. 科里奥利力（Coriolis force）：旋转参考系下的虚拟力，在地球表面上主要由地球自转引起。在北半球，流体运动向右偏；在南半球，流体运动向左偏。
  4. 电磁力（electromagnetic force）：带电流体元在电场和磁场中所受的洛伦兹力，在工程力学中常用于控制电解液、等离子体的运动。

一般地，在动量方程中，我们将表面力（压力和黏性力）视为主要项并将它们分开考虑，而将体力视为源项。

通过压力和9个粘性应力分量，我们定义了流体元素的应力状态。

​​

该图是为x方向动量方程准备的，故只画出了6个表面沿x方向走的力。

注意，图中所示的物理量变化速率仍然运用了一阶泰勒级数展开，具体推导过程见质量守恒部分，方程2-5：

\[\begin{align}f(x+\delta x)=f(x)+\frac{df(x)}{dx} \delta x\end{align} \]

考虑x方向，由于压力P和应力分量\(\tau_{xx}\)，\(\tau_{yx}\)和\(\tau_{zx}\)造成的总力。

如上图所示，考虑一对面（E，W）（东/右侧壁面，西/左侧壁面，被x轴穿过的一对，上北下南左西右东），把图示力加一下，得：

\[\begin{equation} \begin{aligned} &\left[\left(p-\frac{\partial p}{\partial x}\frac{1}{2} \delta x\right)-\left(\tau_{xx}-\frac{\partial\tau_{xx}}{\partial x}\frac{1}{2} \delta x\right)\right]\delta y\delta z \\ &+\left[-\left(p+\frac{\partial p}{\partial x}\frac{1}{2}\delta x\right)+\left(\tau_{xx}+\frac{\partial\tau_{xx}}{\partial x}\frac{1}{2}\delta x\right)\right]\delta y\delta z\\&=\left(-\frac{\partial p}{\partial x}+\frac{\partial\tau_{xx}}{\partial x}\right)\delta x\delta y\delta z \end{aligned} \end{equation} \]

注意，这里考虑了4个力：两个表面上的压力和黏性力，而下面两个方程是不考虑压力的，因为压力垂直于流体元表面 。

类似地，对于另一对面（N，S）（上面和底下面，被y轴穿过），有：

\[\begin{equation} -\left(\tau_{yx}-\frac{\partial\tau_{yx}}{\partial y}\frac{1}{2}\delta y\right)\delta x\delta z+\left(\tau_{yx}+\frac{\partial\tau_{yx}}{\partial y}\frac{1}{2}\delta y\right)\delta x\delta z=\frac{\partial\tau_{yx}}{\partial y}\delta x\delta y\delta z \end{equation} \]

而对于（T，B）（前面和后面，被z轴穿过），有：

\[\begin{equation} -\left(\tau_{zx}-\frac{\partial\tau_{zx}}{\partial z}\frac{1}{2}\delta z \right)\delta x\delta y+\left(\tau_{zx}+\frac{\partial\tau_{zx}}{\partial z}\frac{1}{2}\delta z \right)\delta x\delta y=\frac{\partial\tau_{zx}}{\partial z}\delta x\delta y\delta z \end{equation} \]

把上述方程30、31、32相加，可得由上述表面力导致的，该流体微元每单位体积受到的总力，也就是该流体微元在x方向的每单位体积动量变化速度：

\[\begin{equation} \rho \frac{Du}{Dt}=\frac{\partial(-p+\tau_{xx})}{\partial x}+\frac{\partial\tau_{yx}}{\partial y}+\frac{\partial\tau_{zx}}{\partial z} \end{equation} \]

类似地，y方向的动量方程为：

\[\begin{equation} \rho \frac{Dv}{Dt}=\frac{\partial\tau_{xy}}{\partial x}+\frac{\partial\left(-p+\tau_{yy}\right)}{\partial y}+\frac{\partial\tau_{zy}}{\partial z} \end{equation} \]

z方向的动量方程：

\[\begin{equation} \rho\frac{Dw}{Dt}=\frac{\partial\tau_{xz}}{\partial x}+\frac{\partial\tau_{yz}}{\partial y}+\frac{\partial(-p+\tau_{zz})}{\partial z} \end{equation} \]

能量守恒

基于热力学第一定律：

\[\begin{equation} \begin{aligned} 流体微元能量增加速度=\\净热量增加速度+作用于流体微元上的净功增加速度 \end{aligned} \end{equation} \]

定义一个流体微元/流体粒子每单位体积的能量变化速率为：

\[\rho \frac{DE}{Dt} \]

做功

接下来，以x方向为例：

\[\begin{equation} \begin{aligned} &作用于流体粒子上的表面力所做功的速率=\\ &力\times 力作用方向上的速度 \end{aligned} \end{equation} \]

这样，我们把表示x方向合力的方程30乘以x方向速度分量\(u\)，得：

\[\begin{equation} \begin{aligned} &\left[\left(pu-\frac{\partial(pu)}{\partial x}\frac{1}{2} \delta x\right)-\left(\tau_{xx}u-\frac{\partial(\tau_{xx}u)}{\partial x}\frac{1}{2} \delta x\right)-\left(pu+\frac{\partial(pu)}{\partial x}\frac{1}{2} \delta x\right)\right] \\ &\left.+\left(\tau_{xx}u+\frac{\partial(\tau_{xx}u)}{\partial x}\frac{1}{2}\delta x\right)\right]\delta y\delta z \\ &+\left[-\left(\tau_{yx}u-\frac{\partial\left(\tau_{yx}u\right)}{\partial y}\frac{1}{2}\delta y\right)+\left(\tau_{yx}u+\frac{\partial\left(\tau_{yx}u\right)}{\partial y}\frac{1}{2}\delta y\right)\right]\delta x\delta z \\ &+\left[-\left(\tau_{zx}u-\frac{\partial(\tau_{zx}u)}{\partial z}\frac{1}{2} \delta z\right)+\left(\tau_{zx}u+\frac{\partial(\tau_{zx}u)}{\partial z}\frac{1}{2} \delta z\right)\right]\delta x\delta y \end{aligned} \end{equation} \]

加减乘除一下，可得x方向表面力的净功率为：

\[\begin{equation} \left[\frac{\partial\left(u(-p+\tau_{xx})\right)}{\partial x}+\frac{\partial\left(u\tau_{xy}\right)}{\partial y}+\frac{\partial(u\tau_{xz})}{\partial z}\right]\delta x\delta x\delta y \end{equation} \]

类似地，有y和z方向表面力对x方向流动的净作用功率：

\[\begin{equation} \left[\frac{\partial\left(v\tau_{xy}\right)}{\partial x}+\frac{\partial\left(v\left(-p+\tau_{yy}\right)\right)}{\partial y}+\frac{\partial\left(v\tau_{yz}\right)}{\partial z}\right]\delta x\delta y\delta z \end{equation} \]

\[\begin{equation} \left[\frac{\partial(v\tau_{xz})}{\partial x}+\frac{\partial(v\tau_{xy})}{\partial y}+\frac{\partial\left(v(-p+\tau_{zz})\right)}{\partial z}\right]\delta x\delta y\delta z \end{equation} \]

将上述方程39、40、41相加，可得表面力对流体粒子作用的总功率为：

\[\begin{equation}\begin{aligned} &-\nabla\cdot(p\boldsymbol{u}) \\ &+\frac{\partial(u\tau_{xx})}{\partial x}+\frac{\partial\big(u\tau_{xy}\big)}{\partial y}+\frac{\partial(u\tau_{xz})}{\partial z}\\&+\frac{\partial\big(v\tau_{xy}\big)}{\partial x}+\frac{\partial\big(v\tau_{yy}\big)}{\partial y}+\frac{\partial\big(v\tau_{yz}\big)}{\partial z}\\&+\frac{\partial(w\tau_{xz})}{\partial x}+\frac{\partial\big(w\tau_{yz}\big)}{\partial y} +\frac{\partial(w\tau_{zz})}{\partial z} \end{aligned}\end{equation} \]

热

由于热流穿过流体微元边界，流体微元每单位体积的净热增加速度为：

\[\begin{equation} -\frac{\partial q_x}{\partial x}-\frac{\partial q_y}{\partial y}-\frac{\partial q_z}{\partial z}=-\nabla\cdot \vec{q} \end{equation} \]

此时，我们引入傅里叶定律：

\[\begin{equation} \vec{q}=-k\nabla T \end{equation} \]

这里：

• \(\vec{q}\)是热流密度向量，表示热量通过单位面积的速率；
• k是材料的热导率；
• \(\nabla T\)是温度梯度。

这样，我们能重写方程44的推导结果为：

\[\begin{equation} -\nabla \cdot \vec{q}=\nabla \cdot (k\nabla T) \end{equation} \]

动能

为便于进一步推导内能i或温度T的方程，我们建立动能变化的方程。

该方程实际上就是把动量方程再乘上速度\(\vec{u}\)。

\[\begin{equation} \begin{aligned} &\rho\frac{DK}{Dt} \\ &=- \boldsymbol{u}\cdot\nabla p+u\left(\frac{\partial\tau_{xx}}{\partial x}+\frac{\partial\tau_{xy}}{\partial y}+\frac{\partial\tau_{xz}}{\partial z}\right)+v\left(\frac{\partial\tau_{yx}}{\partial x}+\frac{\partial\tau_{yy}}{\partial y}+\frac{\partial\tau_{yz}}{\partial z}\right) \\ &+w\left(\frac{\partial\tau_{zx}}{\partial x}+\frac{\partial\tau_{zy}}{\partial y}+\frac{\partial\tau_{zz}}{\partial z}\right)+\boldsymbol{u}\cdot S_E \end{aligned} \end{equation} \]

注意：为得到此方程，需要把方程33-35加起来（这里只考虑了压力和黏性力），四则运算一下，再补一个源项\(S_E\)，最后再乘以\(\boldsymbol{u}\)。

总能量方程与内能

总能量方程应为：

\[\begin{equation} \rho\frac{De}{Dt}=-\nabla\cdot(p\boldsymbol{u})+\nabla\cdot(\boldsymbol{u}\cdot\boldsymbol{\tau})+\nabla\cdot(k\nabla T)+S_i++\boldsymbol{u}\cdot\boldsymbol{S}_E \end{equation} \]

其中：

• \(-\nabla\cdot(p\boldsymbol{u})\)是压力对流体做功的速率，\(\nabla\cdot(\boldsymbol{u}\cdot\boldsymbol{\tau})\)是黏性力对流体做功的贡献，二者共同组成描述表面力对流体微元作用的总功率，即方程42；
• \(\nabla \cdot (k\nabla T)\)是热传导的贡献，见方程43-45；
• 新加入的元素：\(S_i\)是单位体积的内部热源项，\(\boldsymbol{u}\cdot\boldsymbol{S}_E\)表示体力对流体微元的作用。
  （PS：这个体力的项我没搞明白，它要不要乘以\(\rho\)？还有它到底该不该存在，为什么就突然引入动能方程了？虽然最后是要被抵消的）

而总动能方程简写为：

\[\begin{equation} \rho\frac{DK}{Dt}=-\boldsymbol{u}\cdot\nabla p+\boldsymbol{u}\cdot\nabla\cdot\boldsymbol{\tau}+\boldsymbol{u}\cdot\boldsymbol{S}_E \end{equation} \]

总能量-总动能=总内能i，对每一小项分别计算：

• 压力项：
  总能量方程中的为\(-\nabla\cdot (\rho \boldsymbol{u})\)，动能方程中的为\(-\boldsymbol{u}\cdot \nabla p\)，根据链式法则可得：

  \[\begin{equation} -\nabla\cdot(p\boldsymbol{u})-(-\boldsymbol{u}\cdot\nabla p)=-p\nabla\cdot\boldsymbol{u} \end{equation} \]

• 黏性力项：
  总能量方程中的为\(\nabla \cdot (\boldsymbol{u}\cdot \boldsymbol{\tau})\)，而动能方程中的为\(\boldsymbol{u} \cdot \nabla \cdot \boldsymbol{\tau}\)，二者相减的结果为：

  \[\begin{equation} \begin{aligned} (\nabla\boldsymbol{u}):\boldsymbol{\tau}&=\tau_{xx}\frac{\partial u}{\partial x}+\tau_{xy}\frac{\partial u}{\partial y}+\tau_{xz}\frac{\partial u}{\partial z}\\&+\tau_{yx}\frac{\partial v}{\partial x}+\tau_{yy}\frac{\partial v}{\partial y}+\tau_{yz}\frac{\partial v}{\partial z}\\&+\tau_{zx}\frac{\partial w}{\partial x}+\tau_{zy}\frac{\partial w}{\partial y}+\tau_{zz}\frac{\partial w}{\partial z} \end{aligned} \end{equation} \]

• 热传导项和内部热源项保留。

将上述方程49和50合并进来，得到内能方程：

\[\begin{equation} \rho\frac{Di}{Dt}=-p\nabla\cdot\boldsymbol{u}+(\nabla\boldsymbol{u}):\boldsymbol{\tau}+\nabla\cdot(k\nabla T)+S_i \end{equation} \]

展开之后的形式为：

\[\begin{equation} \begin{aligned} &\rho\frac{Di}{Dt} \\ &=- p\nabla\cdot\boldsymbol{u}+\nabla\cdot(k\nabla T)\\&+\tau_{xx}\frac{\partial u}{\partial x}+\tau_{yx}\frac{\partial u}{\partial y}+\tau_{zx}\frac{\partial u}{\partial z}\\&+\tau_{xy}\frac{\partial v}{\partial x}+\tau_{yy}\frac{\partial v}{\partial y}+\tau_{zy}\frac{\partial v}{\partial z} \\ &+\tau_{xz}\frac{\partial w}{\partial x}+\tau_{yz}\frac{\partial w}{\partial y}+\tau_{zz}\frac{\partial w}{\partial z}+S_{i} \end{aligned} \end{equation} \]

进一步地，对于不可压缩流体，有：

• 密度恒定；
• 连续性方程，\(\nabla \cdot \boldsymbol{u}=0\);
• 有内能与温度的关系：\(i=cT\)，其中c为比热容。

那么此时由于连续性，不存在压力梯度，压力项消失，内能方程简化为：

\[\begin{equation} \rho\frac{Di}{Dt}=(\nabla\boldsymbol{u}):\boldsymbol{\tau}+\nabla\cdot(k\nabla T)+S_i \end{equation} \]

\[\begin{equation} \begin{aligned} &\rho c\frac{DT}{Dt} \\ &=\nabla\cdot(k\nabla T)\\&+\tau_{xx}\frac{\partial u}{\partial x}+\tau_{yx}\frac{\partial u}{\partial y}+\tau_{zx}\frac{\partial u}{\partial z}\\&+\tau_{xy}\frac{\partial v}{\partial x}+\tau_{yy}\frac{\partial v}{\partial y}+\tau_{zy}\frac{\partial v}{\partial z} \\ &+\tau_{xz}\frac{\partial w}{\partial x}+\tau_{yz}\frac{\partial w}{\partial y}+\tau_{zz}\frac{\partial w}{\partial z}+S_{i} \end{aligned} \end{equation} \]

‍

标签：tau,方程,partial,推导,equation,体积,rho,delta,frac
From： https://www.cnblogs.com/yukibvvd/p/18528647/3-limited-volume-method-derivative-equation-18