能量守恒方程推导

能量守恒方程的推导可能有很多方式，笔者的推导方式应该更容易理解，供大家参考。

1. 简介

以控制体为研究对象，能量守恒定律可以表述为：

控制体内的能量变化率 = 流入控制体的能量 - 流出控制体的能量 + 表面力的做功功率 + 体积源项的功率 + 体积力的做功功率。

其中，流入、流出控制体的能量包括流动携带的能量的和热传导携带的能量，能量是一个标量，流入控制体为正，流出控制体为负。力的功率等于力的大小乘以速度大小，另外我们规定阻力做功为负，动力做功为正，阻力的方向和速度相反，动力的方向和速度相同。

2. 推导

流体的能量包括内能和动能，控制体内的能量计算式为（其中，括弧内的表达式为单位质量的能量 J /kg）：

\[E=\rho\left(\mathrm{e}+\frac{1}{2} V^{2}\right) d x d y d z=\rho\left(\mathrm{e}+\frac{1}{2}\left(v_{x}^{2}+v_{y}^{2}+v_{z}^{2}\right)\right) d x d y d z\]

需要注意的是，上式的密度、内能、速度分量的值均是位于控制体的质心。

控制体内的能量变化率为：

\[\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) d x d y d z\right]}{\partial t}=\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right)\right]}{\partial t} d x d y d z\]

接下来，计算通过介质流动而流入、流出控制体的能量。

根据以下示意图，X 方向上通过介质流动而流入、流出控制体的净能量为：

\[\rho\left(e+\frac{1}{2} V^{2}\right) v_{x} d y d z-\left\{\rho\left(e+\frac{1}{2} V^{2}\right) v_{x}+\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{x}\right]}{\partial x} d x\right\} d y d z=-\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{x}\right]}{\partial x} d x d y d z\left({}^{*}\right)\]

我们注意到，上式（*) 的速度分量 vx 是位于控制体的角落（原点）的，而单位质量能量的参考点是控制体的质心。若我们将密度、能量、速度分量的参考点也调整到质心，则 X 方向上单位面积通过介质流动而流入（下式蓝色，前半部分）、流出（下式红色，后半部分）控制体的净能量也可以表示为：

\[\begin{array}{l} \textcolor{blue}{\left(\rho-\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e-\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{y}-\frac{\partial v_{y}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{z}-\frac{\partial v_{z}}{\partial x} \frac{d x}{2}\right)^{2}\right]\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\right.} \\ \textcolor{red}{-\left(\rho+\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e+\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{y}+\frac{\partial v_{y}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{z}+\frac{\partial v_{z}}{\partial x} \frac{d x}{2}\right)^{2}\right]\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\right.} \end{array}\]

对该式进行简化，略去二次方以上的高价微量（红色部分）：

\[\begin{array}{l} \left(\rho-\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e-\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{y}-\frac{\partial v_{y}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{z}-\frac{\partial v_{z}}{\partial x} \frac{d x}{2}\right)^{2}\right]\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\right.\\ -\left(\rho+\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e+\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{y}+\frac{\partial v_{y}}{\partial x} \frac{d x}{2}\right)^{2}+\left(v_{z}+\frac{\partial v_{z}}{\partial x} \frac{d x}{2}\right)^{2}\right]\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\right.\\ =\left(\rho-\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e-\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[\begin{array}{c} v_{x}^{2}-v_{x} \frac{\partial v_{x}}{\partial x} d x+\textcolor{red}{\left(\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)^{2}}+v_{y}^{2}-v_{y} \frac{\partial v_{y}}{\partial x} d x+ \\ \textcolor{red}{\left(\frac{\partial v_{y}}{\partial x} \frac{d x}{2}\right)^{2}}+v_{z}^{2}-v_{z} \frac{\partial v_{z}}{\partial x} d x+\textcolor{red}{\left(\frac{\partial v_{z}}{\partial x} \frac{d x}{2}\right)^{2}} \end{array}\right]\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\right.\\ -\left(\rho+\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e+\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[\begin{array}{c} v_{x}^{2}+v_{x} \frac{\partial v_{x}}{\partial x} d x+\textcolor{red}{\left(\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)^{2}}+v_{y}^{2}+v_{y} \frac{\partial v_{y}}{\partial x} d x+ \\ \textcolor{red}{\left(\frac{\partial v_{y}}{\partial x} \frac{d x}{2}\right)^{2}}+v_{z}^{2}+v_{z} \frac{\partial v_{z}}{\partial x} d x+\textcolor{red}{\left(\frac{\partial v_{z}}{\partial x} \frac{d x}{2}\right)^{2}} \end{array}\right]\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\right.\\ =\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\left(\rho-\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e-\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[V^{2}-\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)\right]\right)\\ -\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\left(\rho+\frac{\partial \rho}{\partial x} \frac{d x}{2}\right)\left(e+\frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2}\left[V^{2}+\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)\right]\right) \end{array}\] \[\begin{array}{l} =\left(v_{x}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\left(\begin{array}{c} \rho e-\rho \frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2} \rho V^{2}-\rho \frac{1}{2}\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)-\frac{\partial \rho}{\partial x} \frac{d x}{2} e \\ +\textcolor{red}{\frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{\partial e}{\partial x} \frac{d x}{2}}-\frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2} V^{2}+\textcolor{red}{\frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2}\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)} \end{array}\right)\\ -\left(v_{x}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2}\right)\left(\begin{array}{c} \rho e+\rho \frac{\partial e}{\partial x} \frac{d x}{2}+\frac{1}{2} \rho V^{2}+\rho \frac{1}{2}\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)+\frac{\partial \rho}{\partial x} \frac{d x}{2} e \\ +\textcolor{red}{\frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{\partial e}{\partial x} \frac{d x}{2}}+\frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2} V^{2}+\textcolor{red}{\frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2}\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)} \end{array}\right)\\ =\left[\begin{array}{c} v_{x} \rho e-v_{x} \rho \frac{\partial e}{\partial x} \frac{d x}{2}+v_{x} \frac{1}{2} \rho V^{2}-v_{x} \rho \frac{1}{2}\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)-v_{x} \frac{\partial \rho}{\partial x} \frac{d x}{2} e-v_{x} \frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2} V^{2} \\ -\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \rho e+\textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \rho \frac{\partial e}{\partial x} \frac{d x}{2}}-\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \frac{1}{2} \rho V^{2}+\textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \rho\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)}+ \\ \textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \frac{\partial \rho}{\partial x} \frac{d x}{2} e}+\textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2} V^{2}} \end{array}\right]\\ \left[\begin{array}{c} v_{x} \rho e+v_{x} \rho \frac{\partial e}{\partial x} \frac{d x}{2}+v_{x} \frac{1}{2} \rho V^{2}+v_{x} \rho \frac{1}{2}\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)+v_{x} \frac{\partial \rho}{\partial x} \frac{d x}{2} e+v_{x} \frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2} V^{2} \\ +\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \rho e+\textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \rho \frac{\partial e}{\partial x} \frac{d x}{2}}+\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \frac{1}{2} \rho V^{2}+\textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \rho\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)}+ \\ \textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \frac{\partial \rho}{\partial x} \frac{d x}{2} e}+\textcolor{red}{\frac{\partial v_{x}}{\partial x} \frac{d x}{2} \frac{\partial \rho}{\partial x} \frac{d x}{2} \frac{1}{2} V^{2}} \end{array}\right]\\ =-\left[\begin{array}{c} v_{x} \rho \frac{\partial e}{\partial x} d x+v_{x} \rho\left(v_{x} \frac{\partial v_{x}}{\partial x} d x+v_{y} \frac{\partial v_{y}}{\partial x} d x+v_{z} \frac{\partial v_{z}}{\partial x} d x\right)+v_{x} e \frac{\partial \rho}{\partial x} d x+v_{x} \frac{1}{2} V^{2} \frac{\partial \rho}{\partial x} d x \\ +\rho e \frac{\partial v_{x}}{\partial x} d x+\frac{1}{2} \rho V^{2} \frac{\partial v_{x}}{\partial x} d x \end{array}\right] \end{array}\]

将上式乘以 dydz 得到 X 方向上通过介质流动而流入、流出控制体的净能量：

\[\begin{array}{l} \left[\begin{array}{c} \textcolor{blue}{v_{x} \rho \frac{\partial e}{\partial x}}+\textcolor{red}{v_{x} \rho\left(v_{x} \frac{\partial v_{x}}{\partial x}+v_{y} \frac{\partial v_{y}}{\partial x}+v_{z} \frac{\partial v_{z}}{\partial x}\right)}+\textcolor{blue}{v_{x} e \frac{\partial \rho}{\partial x}}+\textcolor{red}{v_{x} \frac{1}{2} V^{2} \frac{\partial \rho}{\partial x}}+\textcolor{blue}{\rho e \frac{\partial v_{x}}{\partial x}} \\ +\textcolor{red}{\frac{1}{2} \rho V^{2} \frac{\partial v_{x}}{\partial x}} \end{array}\right] d x d y d z\\ =-\left[\textcolor{blue}{\frac{\partial\left(\rho e v_{x}\right)}{\partial x}}+\textcolor{red}{\frac{\partial\left(\rho v_{x} \frac{1}{2} V^{2}\right)}{\partial x}}\right] d x d y d z\\ =-\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{x}\right]}{\partial x} d x d y d z \end{array}\]

可以看出，与式（*）的结果是一致的。
根据以下示意图，Y 方向上通过介质流动而流入、流出控制体的净能量为：

\[\rho\left(e+\frac{1}{2} V^{2}\right) v_{y} d y d z-\left\{\rho\left(e+\frac{1}{2} V^{2}\right) v_{y}+\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{y}\right]}{\partial y} d y\right\} d x d z=-\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{y}\right]}{\partial y} d x d y d z\]

据以下示意图，Z 方向上通过介质流动而流入、流出控制体的净能量为：

\[\rho\left(e+\frac{1}{2} V^{2}\right) v_{z} d x d y-\left\{\rho\left(e+\frac{1}{2} V^{2}\right) v_{z}+\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{z}\right]}{\partial z} d z\right\} d x d y=-\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{z}\right]}{\partial z} d x d y d z\]

于是，通过介质流动而流入、流出控制体的总净能量为：

\[-\left\{\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{x}\right]}{\partial x}+\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{y}\right]}{\partial y}+\frac{\partial\left[\rho\left(e+\frac{1}{2} V^{2}\right) v_{z}\right]}{\partial z}\right\} d x d y d z\]

接着，计算通过热传导而流入、流出控制体的能量。根据傅里叶热传导公式，热流密度为：

\[-\lambda \frac{\partial T}{\partial n} \]

根据以下示意图，X 方向上通过热传导而流入、流出控制体的净能量为：

\[-\lambda \frac{\partial T}{\partial x} d y d z-\left[-\lambda \frac{\partial T}{\partial x}+\frac{\partial\left(-\lambda \frac{\partial T}{\partial x}\right)}{\partial x} d x\right] d y d z=\frac{\partial\left[\partial\left(\lambda \frac{\partial T}{\partial x}\right)\right]}{\partial x} d x d y d z\]

根据以下示意图，Y 方向上通过热传导而流入、流出控制体的净能量为：

\[-\lambda \frac{\partial T}{\partial y} d x d z-\left[-\lambda \frac{\partial T}{\partial y}+\frac{\partial\left(-\lambda \frac{\partial T}{\partial y}\right)}{\partial y} d y\right] d x d z=\frac{\partial\left[\partial\left(\lambda \frac{\partial T}{\partial y}\right)\right]}{\partial y} d x d y d z\]

根据以下示意图，Z 方向上通过热传导而流入、流出控制体的净能量为：

\[-\lambda \frac{\partial T}{\partial z} d x d y-\left[-\lambda \frac{\partial T}{\partial z}+\frac{\partial\left(-\lambda \frac{\partial T}{\partial z}\right)}{\partial z} d z\right] d x d y=\frac{\partial\left[\partial\left(\lambda \frac{\partial T}{\partial z}\right)\right]}{\partial z} d x d y d z\]

于是，通过热传导而流入、流出控制体的总净能量为：

\[\left\{\frac{\partial}{\partial x}\left[\partial\left(\lambda \frac{\partial T}{\partial x}\right)\right]+\frac{\partial}{\partial y}\left[\partial\left(\lambda \frac{\partial T}{\partial y}\right)\right]+\frac{\partial}{\partial z}\left[\partial\left(\lambda \frac{\partial T}{\partial z}\right)\right]\right\} d x d y d z\]

接下来，计算表面力的做功功率。

回顾一下如下的控制体受力分析示意图：

根据以下示意图，法向为 X 的表面上的表面力做功功率为（注意前文所述的做功正负号）：

正应力：\[\begin{array}{l} \sigma_{x x}:-\sigma_{x x} d y d z \cdot v_{x}+\left(\sigma_{x x}+\frac{\partial \sigma_{x x}}{\partial x} d x\right) d y d z \cdot\left(v_{x}+\frac{\partial v_{x}}{\partial x} d x\right) \\ =-\sigma_{x x} d y d z \cdot v_{x}+\sigma_{x x} d y d z v_{x}+\sigma_{x x} d y d z \frac{\partial v_{x}}{\partial x} d x+\frac{\partial \sigma_{x x}}{\partial x} d x d y d z v_{x}+\textcolor{red}{\frac{\partial \sigma_{x x}}{\partial x} d x \frac{\partial v_{x}}{\partial x} d x} \\ =\sigma_{x x} \frac{\partial v_{x}}{\partial x} d x d y d z+v_{x} \frac{\partial \sigma_{x x}}{\partial x} d x d y d z \end{array}\]

切应力：\[\begin{array}{l} \tau_{x y}:-\tau_{x y} d y d z \cdot v_{y}+\left(\tau_{x y}+\frac{\partial \tau_{x y}}{\partial x} d x\right) d y d z \cdot\left(v_{y}+\frac{\partial v_{y}}{\partial x} d x\right) \\ =-\tau_{x y} d y d z \cdot v_{y}+\tau_{x y} d y d z v_{y}+\tau_{x y} d y d z \frac{\partial v_{y}}{\partial x} d x+\frac{\partial \tau_{x y}}{\partial x} d x d y d z v_{y}+\textcolor{red}{\frac{\partial \tau_{x y}}{\partial x} d x \frac{\partial v_{y}}{\partial x} d x }\\ =\tau_{x y} \frac{\partial v_{y}}{\partial x} d x d y d z+v_{y} \frac{\partial \tau_{x y}}{\partial x} d x d y d z \\ \tau_{x z}:-\tau_{x z} d y d z \cdot v_{z}+\left(\tau_{x z}+\frac{\partial \tau_{x z}}{\partial x} d x\right) d y d z \cdot\left(v_{z}+\frac{\partial v_{z}}{\partial x} d x\right) \\ =-\tau_{x z} d y d z \cdot v_{z}+\tau_{x z} d y d z v_{z}+\tau_{x z} d y d z \frac{\partial v_{z}}{\partial x} d x+\frac{\partial \tau_{x z}}{\partial x} d x d y d z v_{z}+\textcolor{red}{\frac{\partial \tau_{x z}}{\partial x} d x \frac{\partial v_{z}}{\partial x} d x} \\ =\tau_{x z} \frac{\partial v_{z}}{\partial x} d x d y d z+v_{z} \frac{\partial \tau_{x z}}{\partial x} d x d y d z \end{array}\]

总应力：\[\begin{array}{l} \sigma_{x x} \frac{\partial v_{x}}{\partial x} d x d y d z+v_{x} \frac{\partial \sigma_{x x}}{\partial x} d x d y d z+\tau_{x y} \frac{\partial v_{y}}{\partial x} d x d y d z+v_{y} \frac{\partial \tau_{x y}}{\partial x} d x d y d z \\ +\tau_{x z} \frac{\partial v_{z}}{\partial x} d x d y d z+v_{z} \frac{\partial \tau_{x z}}{\partial x} d x d y d z \end{array}\]

作者：鱼花生