? 本文内容主要是不可压的雷诺平均方程（Reynolds-averaged Navier-Stokes equations，RANS）的推导。在推导雷诺方程之前先来总结一些关于时间平均的运算规律。
1. 关于时间平均的运算规律证明 ? 对于瞬时量 ? = Φ + ? ′ \phi=\Phi+\phi^\prime ?=Φ+?′和ψ = Ψ + ψ ′ \psi=\Psi+\psi^\prime ψ=Ψ+ψ′， Φ \Phi Φ和 Ψ \Psi Ψ已经是与时间无关的量，所以根据时间平均的定义有
? ′ ￣ = ψ ′ ￣ = 0 (1.1) \overline{\phi^\prime}=\overline{\psi^\prime}=0 \tag{1.1} ?′?=ψ′?=0(1.1)
Φ ￣ = 1 Δ t ∫ 0 Δ t Φ d t = Φ 1 Δ t ∫ 0 Δ t d t = Φ (1.2) \overline{\Phi}=\frac{1}{\Delta t} \int_0^{\Delta t} \Phi dt =\Phi \frac{1}{\Delta t} \int_0^{\Delta t} dt=\Phi \tag{1.2} Φ=Δt1?∫0Δt?Φdt=ΦΔt1?∫0Δt?dt=Φ(1.2)
? ? ? s ￣ = 1 Δ t ∫ 0 Δ t ? ? ? s d t = ? ? s ( 1 Δ t ∫ 0 Δ t ? d t ) = ? ? ˉ ? s = ? Φ ? s (1.3) \begin{aligned} \overline{\frac{\partial \phi}{\partial s}} &= \frac{1}{\Delta t} \int_0^{\Delta t} \frac{\partial \phi}{\partial s} dt = \frac{\partial}{\partial s} \left(\frac{1}{\Delta t} \int_0^{\Delta t} \phi dt \right) = \frac{\partial \bar{\phi}}{\partial s} = \frac{\partial \Phi}{\partial s} \end{aligned} \tag{1.3} ?s?????=Δt1?∫0Δt??s???dt=?s??(Δt1?∫0Δt??dt)=?s??ˉ??=?s?Φ??(1.3)
∫ ? d s ￣ = 1 Δ t ∫ 0 Δ t ( ∫ ? d s ) d t = ∫ ( 1 Δ t ∫ 0 Δ t ? d t ) d s = ∫ ? ˉ d s = ∫ Φ d s (1.4) \begin{aligned} \overline{\int \phi ds} &= \frac{1}{\Delta t}\int_0^{\Delta t}\left( \int \phi ds \right)dt =\int \left( \frac{1}{\Delta t}\int_0^{\Delta t} \phi dt\right)ds =\int \bar{\phi}ds = \int\Phi ds \end{aligned} \tag{1.4} ∫?ds??=Δt1?∫0Δt?(∫?ds)dt=∫(Δt1?∫0Δt??dt)ds=∫?ˉ?ds=∫Φds?(1.4)
? + ψ ￣ = 1 Δ t ∫ 0 Δ t ( ? + ψ ) d t = 1 Δ t ∫ 0 Δ t ? d t + 1 Δ t ∫ 0 Δ t ψ d t = Φ + Ψ (1.5) \begin{aligned} \overline{\phi + \psi} &= \frac{1}{\Delta t} \int_0^{\Delta t} (\phi + \psi)dt = \frac{1}{\Delta t} \int_0^{\Delta t} \phi dt + \frac{1}{\Delta t} \int_0^{\Delta t} \psi dt = \Phi + \Psi \end{aligned} \tag{1.5} ?+ψ??=Δt1?∫0Δt?(?+ψ)dt=Δt1?∫0Δt??dt+Δt1?∫0Δt?ψdt=Φ+Ψ?(1.5)
? Ψ ￣ = 1 Δ t ∫ 0 Δ t ( ? Ψ ) d t = Ψ 1 Δ t ∫ 0 Δ t ? d t = Φ Ψ (1.6) \overline{\phi\Psi} = \frac{1}{\Delta t } \int_0^{\Delta t} (\phi \Psi) dt = \Psi \frac{1}{\Delta t} \int_0^{\Delta t} \phi dt = \Phi\Psi \tag{1.6} ?Ψ?=Δt1?∫0Δt?(?Ψ)dt=ΨΔt1?∫0Δt??dt=ΦΨ(1.6)
? ′ Ψ ￣ = 1 Δ t ∫ 0 Δ t ( ? ′ Ψ ) d t = Ψ 1 Δ t ∫ 0 Δ t ? ′ d t = Ψ ? ′ ˉ = 0 (1.7) \overline{\phi^\prime\Psi} = \frac{1}{\Delta t} \int_0^{\Delta t} (\phi^\prime \Psi) dt = \Psi \frac{1}{\Delta t} \int_0^{\Delta t} \phi^\prime dt = \Psi \bar{\phi^\prime} = 0 \tag{1.7} ?′Ψ?=Δt1?∫0Δt?(?′Ψ)dt=ΨΔt1?∫0Δt??′dt=Ψ?′ˉ?=0(1.7)
? ψ ￣ = ( Φ + ? ′ ) ( Ψ + ψ ′ ) ￣ = Φ Ψ + Φ ψ ′ + ? ′ Ψ + ? ′ ψ ′ ￣ = Φ Ψ ￣ + Φ ψ ′ ￣ + ? ′ Ψ ￣ + ? ′ ψ ′ ￣ = Φ Ψ ￣ + 0 + 0 + ? ′ ψ ′ ￣ = Φ Ψ + ? ′ ψ ′ ￣ (1.8) \begin{aligned} \overline{\phi \psi} = \overline{(\Phi + \phi^\prime)(\Psi + \psi^\prime)} &= \overline{\Phi\Psi + \Phi\psi^\prime + \phi^\prime \Psi + \phi^\prime \psi^\prime} \\ &= \overline{\Phi\Psi} + \overline{\Phi\psi^\prime} + \overline{ \phi^\prime \Psi} + \overline{\phi^\prime \psi^\prime} \\ &= \overline{\Phi\Psi} + 0 + 0 + \overline{\phi^\prime \psi^\prime} \\ & = \Phi\Psi + \overline{\phi^\prime \psi^\prime} \end{aligned} \tag{1.8} ?ψ?=(Φ+?′)(Ψ+ψ′)??=ΦΨ+Φψ′+?′Ψ+?′ψ′?=ΦΨ+Φψ′?+?′Ψ?+?′ψ′?=ΦΨ+0+0+?′ψ′?=ΦΨ+?′ψ′??(1.8)
? 以上均为标量的运算，矢量的散度和梯度计算可以通过上述公式扩展而来。
? 本文中均以粗体字母代表矢量，如速度矢量 u \bold u u和其他任意矢量 a \bold a a等，以 ? ? u \nabla \cdot \bold u ??u表示矢量 u \bold u u的散度，以 ? ? \nabla \phi ??表示标量 ? \phi ?的梯度。根据雷诺分解，瞬时速度矢量 u = U + u ′ \bold u = \bold U + \bold u^\prime u=U+u′和瞬时标量 ? = Φ + ? ′ \phi = \Phi + \phi^\prime ?=Φ+?′有如下运算规律：
u = ( u , v , w ) T \bold u = (u,v,w)^T u=(u,v,w)T
u ˉ = U = ( u ˉ , v ˉ , w ˉ , ) = ( U , V , W ) T \bar \bold u = \bold U =(\bar u,\bar v,\bar w,)= (U,V,W)^T uˉ=U=(uˉ,vˉ,wˉ,)=(U,V,W)T
u ′ = ( u ′ , v ′ , w ′ ) T \bold u^\prime = (u^\prime,v^\prime,w^\prime)^T u′=(u′,v′,w′)T
? ? u ￣ = ? u ? x + ? v ? y + ? w ? z ￣ = ? u ˉ ? x + ? v ˉ ? y + ? w ˉ ? z = ? U ? x + ? V ? y + ? W ? z = ? ? U = ? ? u ˉ (1.9) \begin{aligned} \overline{\nabla \cdot \bold u} &= \overline{\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z} } = \frac{\partial \bar u}{\partial x} + \frac{\partial \bar v}{\partial y} + \frac{\partial \bar w}{\partial z}\\ &= \frac{\partial U}{\partial x} + \frac{\partial V}{\partial y} + \frac{\partial W}{\partial z} =\nabla \cdot \bold U = \nabla \cdot \bar \bold u \end{aligned} \tag{1.9} ??u?=?x?u?+?y?v?+?z?w??=?x?uˉ?+?y?vˉ?+?z?wˉ?=?x?U?+?y?V?+?z?W?=??U=??uˉ?(1.9)
? u ￣ = ( ? u ￣ , ? v ￣ , ? w ￣ ) T = [ Φ U + ? ′ u ′ ￣ Φ V + ? ′ v ′ ￣ Φ W + ? ′ w ′ ￣ ] = Φ U + ? ′ u ′ ￣ (1.10) \overline{\phi \bold u} = (\overline{\phi u},\overline{\phi v}, \overline{\phi w})^T= \left[ \begin{matrix} \Phi U + \overline{\phi^\prime u^\prime } \\ \\ \Phi V +\overline{\phi^\prime v^\prime } \\ \\ \Phi W +\overline{\phi^\prime w^\prime } \end{matrix} \right] = \Phi \bold U +\overline{\phi^\prime \bold u^\prime} \tag{1.10} ?u?=(?u?,?v?,?w?)T=???????ΦU+?′u′?ΦV+?′v′?ΦW+?′w′?????????=ΦU+?′u′?(1.10)
? ? ( ? u ) ￣ = ? ? ? u ￣ = ? ? ( Φ U ) + ? ? ( ? ′ u ′ ￣ ) (1.11) \overline{\nabla \cdot (\phi \bold u)} = \nabla \cdot \overline{\phi \bold u} =\nabla \cdot (\Phi \bold U) + \nabla \cdot (\overline{\phi^\prime \bold u^\prime}) \tag{1.11} ??(?u)?=???u?=??(ΦU)+??(?′u′?)(1.11)
? ? ￣ = ( ? ? ? x , ? ? ? y , ? ? ? z ) T ￣ = ( ? ? ˉ ? x , ? ? ˉ ? x , ? ? ˉ ? x ) T = ( ? Φ ? x , ? Φ ? x , ? Φ ? x ) T = ? Φ (1.12) \overline{\nabla \phi} =\overline{\left(\frac{\partial \phi}{\partial x},\frac{\partial \phi}{\partial y},\frac{\partial \phi}{\partial z}\right)^T} = \left( \frac{ \partial \bar \phi}{\partial x} ,\frac{ \partial \bar \phi}{\partial x},\frac{ \partial \bar \phi}{\partial x} \right)^T =\left( \frac{\partial \Phi}{\partial x} ,\frac{\partial \Phi}{\partial x},\frac{\partial \Phi}{\partial x}\right)^T = \nabla \Phi \tag{1.12} ???=(?x???,?y???,?z???)T?=(?x??ˉ??,?x??ˉ??,?x??ˉ??)T=(?x?Φ?,?x?Φ?,?x?Φ?)T=?Φ(1.12)
? ? ? ? ￣ = ? ? ? ? ￣ = ? ? ? Φ (1.13) \overline{\nabla \cdot \nabla \phi} = \nabla \cdot \overline{\nabla \phi} = \nabla \cdot \nabla \Phi \tag{1.13} ?????=?????=???Φ(1.13)
?? ? ? \nabla \cdot \nabla ???也可写作 ? 2 \nabla^2 ?2,即拉普拉斯算子。
文章图片
2. 雷诺平均方程的推导过程 ? 有了上述准备，推导雷诺平均方程就容易多了，这里考虑三维笛卡尔坐标系下瞬态不可压的连续性方程和N-S方程，忽略体积力，方程如下：
? ? u = 0 (2.1) \nabla \cdot \bold u = 0 \tag{2.1} ??u=0(2.1)
【湍流模型（2）——雷诺平均方程】 ? u ? t + ? ? ( u u ) = ? 1 ρ ? p ? x + ν ? 2 u (2.2a) \frac{\partial u}{\partial t} + \nabla \cdot (u \bold u) = -\frac{1}{\rho} \frac{\partial p}{\partial x} + \nu \nabla^2 u \tag{2.2a} ?t?u?+??(uu)=?ρ1??x?p?+ν?2u(2.2a)
? v ? t + ? ? ( v u ) = ? 1 ρ ? p ? y + ν ? 2 v (2.2b) \frac{\partial v}{\partial t} + \nabla \cdot (v \bold u) = -\frac{1}{\rho} \frac{\partial p}{\partial y} + \nu \nabla^2 v \tag{2.2b} ?t?v?+??(vu)=?ρ1??y?p?+ν?2v(2.2b)
? w ? t + ? ? ( w u ) = ? 1 ρ ? p ? z + ν ? 2 w (2.2c) \frac{\partial w}{\partial t} + \nabla \cdot (w \bold u) = -\frac{1}{\rho} \frac{\partial p}{\partial z} + \nu \nabla^2 w \tag{2.2c} ?t?w?+??(wu)=?ρ1??z?p?+ν?2w(2.2c)
? 雷诺平均方程的思想就是先将方程中的流动变量分解成平均量和脉动量之和，然后对方程两边同时取时间平均。根据雷诺分解，方程中的各瞬时量分解如下：
u = U + u ′ u = U + u ′ v = V + v ′ w = W + w ′ p = P + p ′ (2.3) \bold u = \bold U + \bold u^\prime \quad u=U+u^\prime \quad v=V+v^\prime \quad w=W+w^\prime \quad p=P + p^\prime \tag{2.3} u=U+u′u=U+u′v=V+v′w=W+w′p=P+p′(2.3)
将 ( 2.3 ) (2.3) (2.3)各式代入到式 ( 2.1 ) (2.1) (2.1)~ ( 2.2 c ) (2.2c) (2.2c)中，然后对方程取时间平均，再根据第1节中的运算规律对方程进行整理。
文章图片
先来推导方程 ( 2.1 ) (2.1) (2.1)，根据式 ( 1.9 ) (1.9) (1.9)
? ? u ￣ = ? ? u ￣ = ? ? U \overline{\nabla \cdot \bold u} = \nabla \cdot \overline \bold u = \nabla \cdot \bold U ??u=??u=??U
? 于是得到平均流的连续性方程：
? ? U = 0 (2.4) \nabla \cdot \bold U = 0 \tag{2.4} ??U=0(2.4)
? 其中 U \bold U U是平均流速度矢量，写成分量形式为
? U ? x + ? V ? y + ? W ? z = 0 \frac{\partial U}{\partial x} + \frac{\partial V}{\partial y} +\frac{\partial W}{\partial z} = 0 ?x?U?+?y?V?+?z?W?=0
推导方程 ( 2.2 a ) (2.2a) (2.2a)：
? 方程左边第一项，根据式 ( 1.3 ) (1.3) (1.3)
? u ? t ￣ = ? U ? t \overline{\frac{\partial u}{\partial t}} = \frac{\partial U}{\partial t} ?t?u??=?t?U?
? 左边第二项，根据式 ( 1.11 ) (1.11) (1.11)
? ? ( u u ) ￣ = ? ? ( U U ) + ? ? ( u ′ u ′ ￣ ) \overline{\nabla \cdot (u\bold u)} = \nabla \cdot (U\bold U) + \nabla \cdot (\overline{u^\prime \bold u^\prime}) ??(uu)?=??(UU)+??(u′u′)
? 方程右边第一项，
? 1 ρ ? p ? x ￣ = ? 1 ρ ? p ˉ ? x = ? 1 p ? P ? x \overline{-\frac{1}{\rho} \frac{\partial p}{\partial x}} = -\frac{1}{\rho} \frac{\partial \bar p}{\partial x} = -\frac{1}{p} \frac{\partial P}{\partial x} ?ρ1??x?p??=?ρ1??x?pˉ??=?p1??x?P?
? 右边第二项，根据式 ( 1.13 ) (1.13) (1.13)
ν ? 2 u ￣ = ν ? 2 U \overline{\nu \nabla^2 u} = \nu \nabla^2 U ν?2u=ν?2U
? 综上，对方程 ( 2.2 a ) (2.2a) (2.2a)取时间平均得到 x x x分量的时均动量方程：
? U ? t + ? ? ( U U ) + ? ? ( u ′ u ′ ￣ ) = ? 1 ρ ? P ? x + ν ? 2 U ( 1 ) ( 2 ) ( 3 ) ( 4 ) ( 5 ) (2.5a) \begin{aligned} \frac{\partial U}{\partial t} + \nabla \cdot (U \bold U) + \nabla \cdot (\overline{u^\prime \bold u^\prime}) = -\frac{1}{\rho} \frac{\partial P}{\partial x} + \nu \nabla^2U \tag{2.5a} \\ (1) \quad \qquad (2) \qquad \qquad (3) \qquad \qquad (4) \qquad \quad(5) \end{aligned} ?t?U?+??(UU)+??(u′u′)=?ρ1??x?P?+ν?2U(1)(2)(3)(4)(5)?(2.5a)
同样的推导过程可以得到 y y y和 z z z分量的时均动量方程：
? V ? t + ? ? ( V U ) + ? ? ( v ′ u ′ ￣ ) = ? 1 ρ ? P ? y + ν ? 2 V (2.5b) \frac{\partial V}{\partial t} + \nabla \cdot (V \bold U) + \nabla \cdot (\overline{v^\prime \bold u^\prime}) = -\frac{1}{\rho} \frac{\partial P}{\partial y} + \nu \nabla^2V \tag{2.5b} ?t?V?+??(VU)+??(v′u′)=?ρ1??y?P?+ν?2V(2.5b)
? W ? t + ? ? ( W U ) + ? ? ( w ′ u ′ ￣ ) = ? 1 ρ ? P ? z + ν ? 2 W (2.5c) \frac{\partial W}{\partial t} + \nabla \cdot (W \bold U) + \nabla \cdot (\overline{w^\prime \bold u^\prime}) = -\frac{1}{\rho} \frac{\partial P}{\partial z} + \nu \nabla^2W \tag{2.5c} ?t?W?+??(WU)+??(w′u′)=?ρ1??z?P?+ν?2W(2.5c)
以方程 ( 2.5 a ) (2.5a) (2.5a)为例，可以看到第（2）、（3）项都是来自原方程的对流项，第（3）项是湍流脉动量的乘积，是时间平均后多出来的项，其代表的物理含义是湍流涡的对流所引起动量输运，这就类似于流体在粘性剪应力的作用下又附加了湍流剪应力。所以习惯上，这一项会移到方程的右边和粘性项放在一起，以反应其物理属性。为了展示湍流应力的结构，把第（3）项中的脉动速度矢量写出分量形式，并移到等号右边，则有
? U ? t + ? ? ( U U ) = ? 1 ρ ? P ? x + ν ? 2 U + 1 ρ [ ? ( ? ρ u ′ 2 ￣ ) ? x + ? ( ? ρ u ′ v ′ ￣ ) ? y + ? ( ? ρ u ′ w ′ ￣ ) ? z ] (2.6a) \begin{aligned} \frac{\partial U}{\partial t} + \nabla \cdot (U \bold U) =& -\frac{1}{\rho} \frac{\partial P}{\partial x}+\nu \nabla^2U \\ &+\frac{1}{\rho}\left[ \frac{\partial (-\rho \overline{u^{\prime 2}})}{\partial x} + \frac{\partial (-\rho \overline{u^{\prime }v^\prime})}{\partial y} +\frac{\partial (-\rho \overline{u^{\prime }w^\prime})}{\partial z} \right] \tag{2.6a} \end{aligned} ?t?U?+??(UU)=??ρ1??x?P?+ν?2U+ρ1?[?x?(?ρu′2)?+?y?(?ρu′v′)?+?z?(?ρu′w′)?]?(2.6a)
? V ? t + ? ? ( V U ) = ? 1 ρ ? P ? y + ν ? 2 V + 1 ρ [ ? ( ? ρ u ′ v ′ ￣ ) ? x + ? ( ? ρ v ′ 2 ￣ ) ? y + ? ( ? ρ v ′ w ′ ￣ ) ? z ] (2.6b) \begin{aligned} \frac{\partial V}{\partial t} + \nabla \cdot (V \bold U) =& -\frac{1}{\rho} \frac{\partial P}{\partial y}+\nu \nabla^2V \\ & +\frac{1}{\rho}\left[ \frac{\partial (-\rho \overline{u^{\prime} v^\prime})}{\partial x} + \frac{\partial (-\rho \overline{v^{\prime2} })}{\partial y} +\frac{\partial (-\rho \overline{v^{\prime }w^\prime})}{\partial z} \right] \tag{2.6b} \end{aligned} ?t?V?+??(VU)=??ρ1??y?P?+ν?2V+ρ1?[?x?(?ρu′v′)?+?y?(?ρv′2)?+?z?(?ρv′w′)?]?(2.6b)
? W ? t + ? ? ( W U ) = ? 1 ρ ? P ? z + ν ? 2 W + 1 ρ [ ? ( ? ρ u ′ w ′ ￣ ) ? x + ? ( ? ρ v ′ w ′ ￣ ) ? y + ? ( ? ρ w ′ 2 ￣ ) ? z ] (2.6c) \begin{aligned} \frac{\partial W}{\partial t} + \nabla \cdot (W \bold U) =& -\frac{1}{\rho} \frac{\partial P}{\partial z}+\nu \nabla^2W \\ & +\frac{1}{\rho}\left[ \frac{\partial (-\rho \overline{u^{\prime} w^\prime})}{\partial x} + \frac{\partial (-\rho \overline{v^{\prime }w^\prime})}{\partial y} +\frac{\partial (-\rho \overline{w^{\prime2 }})}{\partial z} \right] \tag{2.6c} \end{aligned} ?t?W?+??(WU)=??ρ1??z?P?+ν?2W+ρ1?[?x?(?ρu′w′)?+?y?(?ρv′w′)?+?z?(?ρw′2)?]?(2.6c)
湍流应力一共有6个独立分量，包括3个法向应力分量：
τ x x = ? ρ u ′ 2 ￣ τ y y = ? ρ v ′ 2 ￣ τ z z = ? ρ w ′ 2 ￣ \tau_{xx}=-\rho \overline{u^{\prime 2}} \quad \tau_{yy}=-\rho \overline{v^{\prime 2}} \quad \tau_{zz}=-\rho \overline{w^{\prime 2}} τxx?=?ρu′2τyy?=?ρv′2τzz?=?ρw′2
和3个切应力分量：
τ x y = τ y x = ? ρ u ′ v ′ ￣ τ x z = τ z x = ? ρ u ′ w ′ ￣ τ z y = τ y z = ? ρ w ′ v ′ ￣ \tau_{xy}=\tau_{yx}=-\rho \overline{u^\prime v^\prime} \quad \tau_{xz}=\tau_{zx}=-\rho \overline{u^\prime w^\prime} \quad \tau_{zy}=\tau_{yz}=-\rho \overline{w^\prime v^\prime} \quad τxy?=τyx?=?ρu′v′τxz?=τzx?=?ρu′w′τzy?=τyz?=?ρw′v′
这些湍流应力称为雷诺应力，方程式 ( 2.4 ) (2.4) (2.4)和 ( 2.6 a ) (2.6a) (2.6a)~ ( 2.6 c ) (2.6c) (2.6c)称为雷诺平均N-S方程(RANS)。
参考资料：
- Versteeg H K , Malalasekera W . An introduction to computational fluid dynamics : the finite volume method = 计算流体动力学导论[M]. 世界图书出版公司, 2010.
推荐阅读
- 有限体积法（6）——离散格式的特性
- 有限体积法（5）——对流-扩散方程的离散
- 有限体积法（4）——一维扩散方程数值求解（第二类边界条件）
- CFD|有限体积法（3）——一维扩散方程数值求解（第一类边界条件）
- 【小白的CFD之旅】08 CFD速成之道
- Star|Star CCM+（RPM要被淘汰了）
- 有限体积法（9）——高阶差分格式（QUICK格式）
- 有限体积法（10）——格式精度与待定系数法
- 有限体积法（8）——混合差分格式
- 有限体积法（7）——迎风格式