雷诺平均方法

湍流是瞬态的三维流动。工程师可能不关心所有细节，而只关心流动表现出的总体特性。比如，只关心某个物理量的平均值而不是瞬态值，因此用雷诺平均方程（Reynolds-Averaged Navier-Stokes, RANS）来模拟湍流。RANS方法将湍流中的所有物理量表示成时间和空间的函数，并将其分解为时间平均项和脉动项的叠加：

$$ \begin{equation} \phi(x_i, t) = \underbrace{\overline{\phi}(x_i)}_ % workaround for hugo {\text{时间平均项}} + \underbrace{\phi^\prime (x_i, t)}_{\text{脉动项}} \end{equation} $$

对于稳态问题，方程中的时间平均项可以表示为

$$ \begin{equation} \overline{\phi}(x_i) = \lim_{T \to \infty} \frac{1}{T} \int_0^T \phi(x_i, t) \mathrm{d}t \end{equation} $$

其中 $t$ 为时间，$T$ 为平均时间。$T$ 必须远大于湍流脉动的时间尺度，才能保证 $\phi$ 的值不受时间平均初始时刻的影响，因此这里取 $T \to \infty$。

而对于非稳态流动，时间平均项也是时间的函数

$$ \begin{equation} \phi (x_i, t) = \overline{\phi}(x_i, t) + \phi^\prime (x_i, t) \end{equation} $$

对N个试验进行系综平均（ensemble average），有

$$ \begin{equation} \overline{\phi}(x_i, t) = \lim_{N \to \infty} \frac{1}{N} \sum_{n=1}^N \phi(x_i, t) \end{equation} $$

系综平均有如下性质：

$$ \begin{equation} \overline{\overline{\phi}} = \overline{\phi}, \quad \overline{\phi^\prime} = 0 \end{equation} $$

$$ \begin{equation} \overline{\alpha \phi + \beta \psi} = \alpha \overline{\phi} + \beta \overline{\psi} \quad (\alpha 和 \beta 为常数) \end{equation} $$

$$ \begin{equation} \begin{aligned} \overline{\phi \psi} &= \overline{(\overline{\phi}+\phi^\prime)(\overline{\psi}+\psi^\prime)} \newline &= \overline{(\overline{\phi},\overline{\psi} + \overline{\psi} \phi^\prime + \overline{\phi}\psi^\prime + \phi^\prime \psi^\prime)} \newline &= \overline{\overline{\phi},\overline{\psi}} + \overline{\overline{\psi} \phi^\prime} + \overline{\overline{\phi}\psi^\prime} + \overline{\phi^\prime \psi^\prime} \newline &= \overline{\phi},\overline{\psi} + \overline{\phi^\prime \psi^\prime} \end{aligned} \end{equation} $$

$$ \begin{equation} \frac{\overline{\partial \phi}}{\partial t} = \frac{\partial \overline{\phi}}{\partial t}, \quad \frac{\overline{\partial \phi}}{\partial x_i} = \frac{\partial \overline{\phi}}{\partial x_i} \end{equation} $$

动量方程

动量方程的形式如下

$$ \begin{equation} \underbrace{\frac {\partial (\rho {u}_ % workaround for hugo i)}{\partial t}}_ % workaround for hugo {\text{时间项}} + \underbrace{\frac{\partial (\rho u_j u_i)}{\partial x_j}}_ % workaround for hugo {\text{对流项}} - \underbrace{\frac{\partial}{\partial x_j} \left[ \mu \left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) \right]}_{\text{粘性项}} = - \frac{\partial p}{\partial x_i} \end{equation} $$

对动量方程雷诺平均后会出现雷诺应力（Reynolds stresses）：

$$ \begin{equation} \frac {\partial (\rho \overline{u}_i)}{\partial t} + \frac{\partial}{\partial x_j}( \rho \overline{u}_i \overline{u}_j + \rho \overline{u_i^\prime u_j^\prime}) - \frac{\partial}{\partial x_j}\left[\mu\left( \frac{\partial \overline{u}_i}{\partial x_j} + \frac{\partial \overline{u}_j}{\partial x_i} \right) \right] = - \frac{\partial \overline{p}}{\partial x_i} \end{equation} $$

通常将雷诺应力放到粘性项中，即

$$ \begin{equation} \frac {\partial (\rho \overline{u}_ % workaround for hugo i)}{\partial t} + \frac{\partial}{\partial x_j}( \rho \overline{u}_ % workaround for hugo i \overline{u}_ % workaround for hugo j) - \frac{\partial}{\partial x_j}\left[\mu\left( \frac{\partial \overline{u}_ % workaround for hugo i}{\partial x_j} + \frac{\partial \overline{u}_ % workaround for hugo j}{\partial x_i} \right) \underbrace{-\rho \overline{u_i^\prime u_j^\prime}}_{\text{雷诺应力}} \right] = - \frac{\partial \overline{p}}{\partial x_i} \end{equation} $$

对不可压缩流体，常写作

$$ \begin{equation} \frac {\partial \overline{u}_ % workaround for hugo i}{\partial t} + \frac{\partial}{\partial x_j}\left( \overline{u}_ % workaround for hugo i \overline{u}_ % workaround for hugo j \right) - \frac{\partial}{\partial x_j}\left[\nu\left( \frac{\partial \overline{u}_ % workaround for hugo i}{\partial x_j} + \frac{\partial \overline{u}_ % workaround for hugo j}{\partial x_i} \right) \underbrace{- \overline{u_i^\prime u_j^\prime}}_{\text{雷诺应力}} \right] = - \frac{1}{\rho}\frac{\partial \overline{p}}{\partial x_i} \label{eq:ranse} \end{equation} $$

上式被称为 RANS 方程，其中的 $-\overline{u_i^\prime u_j^\prime}$ （或 $-\rho \overline{u_i^\prime u_j^\prime}$ ）就是雷诺应力，通常记为 $R_{ij}$。写作分量形式如下

$$ \begin{equation} R_{ij} \equiv - \overline{u_i^\prime u_j^\prime} = -\left ( \begin{matrix} \overline{u^\prime u^\prime} & \overline{u^\prime v^\prime} & \overline{u^\prime w^\prime} \newline \overline{v^\prime u^\prime} & \overline{v^\prime v^\prime} & \overline{v^\prime w^\prime} \newline \overline{w^\prime u^\prime} & \overline{w^\prime v^\prime} & \overline{w^\prime w^\prime} \newline \end{matrix} \right ) \end{equation} $$

$R_{ij}$ 是对角对称的，包含六个未知量：三个正应力分量（对角）和三个切应力分量（非对角）。

封闭问题

要求解 RANS 方程，直观的做法是求解各应力分量的对流输运方程。对应力分量的对流输运方程进行雷诺平均，会得到更高阶的未知量 $\overline{u_i^\prime u_j^\prime u_k^\prime}$ 。实际上这样下去会得到无限多的未知量，整个方程组永远无法封闭。这就需要采用一些近似的模型来封闭方程组，比如用某些平均量来近似表示雷诺应力。这些不同的近似方法则被统称为湍流模型（turbulence models）。

涡粘模型

从数学形式上看，雷诺应力很像二阶应力张量。因此 Boussinesq 仿照分子粘性应力与速度变形率的关系，提出了关于涡粘性（eddy viscosity）的 Boussinesq 假设，认为雷诺应力与平均速度变形率成线性关系：

$$ \begin{equation} -\overline{u_i^\prime u_j^\prime} = \nu_t \left( \frac{\partial \overline{u}_ % workaround for hugo i}{\partial x_j} + \frac{\partial \overline{u}_ % workaround for hugo j}{\partial x_i} - \frac{2}{3} \frac{\partial \overline{u}_ % workaround for hugo k}{\partial x_k} \delta_{ij} \right) - \frac{2}{3} k \delta_{ij} \label{eq:boussinesq} \end{equation} $$

式 $\eqref{eq:boussinesq}$ 中的 $k$ 为湍动能（turbulence kinetic energy），定义为 $k=\frac{1}{2}(\overline{u^\prime u^\prime}+\overline{v^\prime v^\prime}+\overline{w^\prime w^\prime})$，$\delta_{ij}$ 为 Kronecker delta。$\nu_t$ 被称为湍流粘度或涡粘系数，湍流粘度和一般流体粘度不同，它不是流体的固有物理属性，而和当地的湍流强度有很大关系。对于流场中不同位置的点，$\nu_t$ 的值都不尽相同。

注意对不可压缩流体，有$\frac{\partial \overline{u}_k}{\partial x_k}=0$ 。将式 $\eqref{eq:boussinesq}$ 代入 RANS 方程 $\eqref{eq:ranse}$，简化后可以得到

$$ \begin{equation} \frac {\partial \overline{u}_ % workaround for hugo i}{\partial t} + \frac{\partial}{\partial x_j}\left( \overline{u}_ % workaround for hugo i \overline{u}_ % workaround for hugo j \right) - \frac{\partial}{\partial x_j}\left[\nu_\mathit{eff}\left( \frac{\partial \overline{u}_i}{\partial x_j} + \frac{\partial \overline{u}_j}{\partial x_i} \right) \right] = - \frac{1}{\rho}\frac{\partial p^*}{\partial x_i}\label{eq:ranse_final} \end{equation} $$

其中

  • $\nu_\mathit{eff} = \nu+\nu_t$
  • $p^* = \overline{p} + \frac{2}{3} \rho k \delta_{ij}$

注意 RANS 方程 $\eqref{eq:ranse_final}$ 中的压力 $p^*$ 包含湍动能 $k$ 的影响 。

基于以上假设的模型被称为涡粘模型（Eddy Viscosity Model, EVM）。涡粘模型一般通过求解额外的偏微分方程得到 $\nu_t$ ，并根据需要求解的方程数量将涡粘模型分为：零方程模型、一方程模型、两方程模型等。工业界常用的涡粘模型有一方程 Spalart-Allmaras 模型，两方程 $k-\varepsilon$ 及 $k-\omega$ 模型等。

雷诺应力模型

涡粘模型基于湍流各项同性假设，对于各向异性的复杂流动（如大曲率流动、强旋流、冲击流等）效果不佳。这时可采用雷诺应力模型（Reynolds Stress Model, RSM）在更高阶的应力上封闭 RANS 方程。具体做法通常是求解六个未知应力分量的方程及其对应的耗散率方程。雷诺应力模型可以较好的模拟具有强烈非线性特征的复杂流动，然而由于它需要求解许多额外的方程，因此计算量相比涡粘模型大得多，在工程中应用不多。常见的雷诺应力模型有 LRR （Launder-Reece-Rodi）和 SSG（Speziale-Sarkar-Gatski） 模型。