湍流模型和离散格式在汽车外流场计算中的对比

摘要

以35° Ahmed汽车标准模型为对象，基于OpenFOAM软件研究Realizable k⁃ε （RKE）， Spalart⁃Allmaras （SA）和Shear Stress Transport k⁃ω （SST）3种湍流模型以及Linear Upwind （LU）， Linear UpwindV （LUV）， Localblend， Total Variation Diminishing （TVD）4种速度对流离散格式在Ahmed模型外流场计算中的表现。同时，对比了OpenFOAM、 FLUENT和STARCCM+这3个软件平台的差异。结果表明：OpenFOAM中RKE和SST湍流模型预测的阻力系数较SA模型更接近实验结果；LUV和Localblend等更稳定的格式可以一定程度上减弱计算不稳定性带来的尾迹振荡问题，但只有TVD格式能够完全的抑制这种震荡行为。3种软件平台计算的阻力系数均与实验吻合较好，其中OpenFOAM采用TVD格式的计算误差最小，仅为0.7%。此外， OpenFOAM和STARCCM+可以得到较为合理和相近的三维尾迹流场结构，而FLUENT则会高估流向C柱涡的强度。

关键词

汽车空气动力学; OpenFOAM软件; 湍流模型; 对流格式; 定常计算

作为汽车空气动力学研究中的标准模型，Ahmed模型^［1］被广泛用于计算流体力学（CFD）算法的验证开发和风洞试验标定。随着斜面后倾角的变化，Ahmed模型尾部流场会呈现出C柱涡（流向涡）、斜面分离泡和马蹄涡等流动结构的相互耦合^［2］。对于Ahmed模型的气动阻力，压差阻力占比多达85%，并且大部分是来源于模型背部和尾部斜面。在斜面后倾角0°~30°变化时，由于C柱涡、斜面分离泡和马蹄涡等流动结构随倾角的不断变化，使得气动阻力系数先减少后增加；后倾角大于30°后，气流在斜面前缘发生分离，C柱涡破裂，斜面和背部压力明显回升，使得气动阻力系数显著下降。

对于25° Ahmed模型，刘训等人^［3］使用了Standard k⁃ε （SKE）、 Renormalization group k⁃ε （RNG）， Realizable k⁃ε （RKE）这3种湍流模型进行了数值仿真，结果表明RKE模型计算的气动阻力系数和实验最接近，但是并没针对其他流场变量进行对比。贾志浩^［4］基于半模25° Ahmed模型， 对比了SKE、 RNG、 Standard k⁃ω （SKW）和SST这4种湍流模型计算的阻力值， 发现SST模型的仿真结果较好，平均误差在7%以内。Guilmineau 等^［2］使用Explicit Algebraic Stress Model （EARSM）、 SST模型也计算了这一模型， 发现这两种模型得到的平均速度场和中截面流线图等均无法和实验结果取得很好的一致性。其他一些学者也采用了不同的湍流模型去计算该模型，但均没有得到很好的结果^［2， 5-7］ 。对于35° Ahmed模型，田思等人^［8］对比了不同网格策略和湍流模型对外流场计算的影响，结果表明RKE模型在计算阻力上最为准确，但是其缺乏对流场信息的对比。在Guilmineau等^［2］的研究结果中，两种湍流模型的结果均能和实验结果取得较好的一致性，其中EARSM计算得到的阻力系数偏小，而SST计算的结果偏大。Ashton等^［5］对比了SST、RKE、Elliptic Blending k⁃ε v²⁃f （B-EVM）、Elliptic Blending Reynolds Stress Model （EB-RSM） 的湍流模型在计算时的表现，发现各种模型的结果也都比较合理，各种模型的结果只有一些细微的差别。另外，国内的学者们多是仅从阻力系数这一量来衡量计算的优劣，但评价流场准确度的量应该是多角度的，这样才能得到更合理的结论。

另一方面，当前的数值研究多是基于Fluent、STAR CCM+等商业软件，而基于开源软件OpenFOAM对Ahmed模型进行数值研究的较少。而OpenFOAM作为一种使用非常广泛的开源软件平台，比商业软件更适合进行数值方法的探索。崔阳阳^［9］基于OpenFOAM建立了基于欧拉双流体方法的多相流模型，实现了侧风条件下的单相空气流动和降雪条件下的多相流动的数值模拟。吴德久^［10］以OpenFOAM自带的Smagorinsky-Lilly Model （SLM） 亚格子模型为基础，建立了Dynamic Smagorinsky-Lilly Model （DLS） 和Wall-Adapting Local Eddy-Viscosity Model （WALEM） 亚格子模型，实现了3种亚格子模型（SLM、 DLSM、 WALEM）的大涡模拟方法并且添加壁面函数kqRWallFunction。目前，针对开源软件中的湍流模型，速度对流离散格式等在地面交通工具的外流场计算中的适用性研究还比较少，仅有李田等^［11］研究了不同数值格式和湍流模型^［12］对计算高速列车外流场的影响，发现除了一阶迎风格式，其他二阶格式得到的结果都十分相似，而SST模型更适合计算高铁侧风。此外，不同于OpenFOAM开源软件， 商业软件为了让软件具有更强的稳定性和更好的收敛性，其通常会选择更稳定的计算设置，这会在一定程度上损失精度。目前，针对汽车外流场计算，不同软件平台计算结果对比的研究鲜有学者涉及。

基于上述分析，本文以35° Ahmed模型为研究对象，重点对比分析OpenFOAM软件平台提供的常见湍流模型以及一些典型对流格式的差别，同时对比不同软件平台（OpenFOAM、STARCCM+、FLUENT）的计算结果差异。在湍流模型方面，本文选择了使用较为广泛的RKE、SA及SST这3个湍流模型。在速度对流离散格式方面，本文选择Linear Upwind （LU）、 Linear UpwindV （LUV）、 Localblend、 Total Variation Diminishing （TVD）这4种数值离散格式。

1 数值方法

1.1　湍流模型

目前，工程上仍较多采用雷诺平均法（Reynold Average Navier⁃Stokes，RANS）对瞬时的湍流控制方程进行处理，但由于方程中多出了与脉动量有关的项（雷诺应力），研究者需通过对雷诺应力作出某种假定，即引入新的湍流模型，来封闭方程组。目前常用的湍流模型有两大类：雷诺应力模型和涡粘模型。涡粘模型是基于Boussinesq方程提出的涡粘假设，该假设建立了雷诺应力相对于平均速度梯度的关系。依据确定涡粘系数的微分方程数目多少，涡粘模型包括：零方程模型、一方程模型、两方程模型，如SA^［13］、RKE ^［14］、SST^［15］等湍流模型。以下简单介绍本文比较的SA、RKE和SST这3个湍流模型，关于模型的详细介绍可参考相关文献^［13-15］。

1.1.1　SA模型

SA模型^［13］是基于经验和量纲分析来对湍流粘度直接进行建模的，其假设生成项与当地涡量成正比，最终得到的涡粘系数$\tilde{v}$的控制方程如下：

$$\begin{array}{c}\frac{\partial \tilde{v}}{\partial t}+{\dot{u}}_j\frac{\partial \tilde{v}}{\partial x_j}=c_{b\mathrm{1}}\tilde{S}\tilde{v}-c_{w\mathrm{1}}{f}_w{\left(\frac{\tilde{v}}{d}\right)}^{\mathrm{2}}+\\ \frac{\mathrm{1}}{\sigma }\frac{\partial }{\partial x_k}(\nu +\tilde{\nu })\frac{\partial \tilde{\nu }}{\partial x_j}+\frac{c_{b\mathrm{2}}}{\sigma }\frac{\partial \tilde{\nu }}{\partial x_k}\frac{\partial \tilde{\nu }}{\partial x_k}\end{array}$$

（1）

为了保证近壁处湍流粘度系数快速衰减到0，引入衰减函数：

$$f_{v\mathrm{1}}=\frac{{\chi }^{\mathrm{3}}}{{\chi }^{\mathrm{3}}+c_{v\mathrm{1}}^{\mathrm{3}}}$$

（2）

其中：$\chi =\tilde{v}/v$。相比于两方程模型，SA模型计算量小且稳定性高。

1.1.2　RKE模型

由于标准的k⁃ε模型对时均应变率特别大的情形，有可能导致负的正应力。为了使流动符合湍流的物理定律，将湍流粘度计算式中的系数$c_{\mu }$与应变率联系起来^［14］。得到如下的控制方程：

$$\frac{\partial k}{\partial t}+{\dot{u}}_j\frac{\partial k}{\partial x_j}=-\overline{u{\text{'}}_{i}u{\text{'}}_j}\frac{\partial {\dot{u}}_i}{\partial x_j}-\epsilon +\frac{\partial }{\partial x_j}\left[\left(\nu +{\nu }_t/{\sigma }_k\right)\frac{\partial k}{\partial x_j}\right]$$

（3）

$$\frac{D\epsilon }{Dt}=-c_{\epsilon \mathrm{1}}\frac{\epsilon }{k}\overline{{u^{\text{'}}}_i{u^{\text{'}}}_l}\frac{\partial {\dot{u}}_i}{\partial x_l}-c_{\epsilon \mathrm{2}}\frac{{\epsilon }^{\mathrm{2}}}{k}+\frac{\partial }{\partial x_l}\left[(\nu +{\nu }_t/{\sigma }_{\epsilon })\frac{\partial \epsilon }{\partial x_l}\right]$$

（4）

式中： ${\nu }_t=c_{\mu }\frac{k^{\mathrm{2}}}{\epsilon }$， $c_{\mu }=\frac{\mathrm{1}}{A_{\mathrm{0}}+A_S{U}^{*}k/\epsilon }$。

RKE模型已经被有效的应用于不同类型的流动模拟，在汽车空气动力学领域也得到了广泛的应用。

1.1.3　SST模型

SST模型^［15］在近壁处使用k⁃ω模型，在外层使用k⁃ε模型。ε方程和ω方程在形式上十分相似，仅有系数差异，因此可以将其统一写成如下统一的k⁃ω形式：

$$\begin{array}{c}\frac{\partial k}{\partial t}+{\dot{u}}_j\frac{\partial k}{\partial x_j}=-\overline{u{\text{'}}_{i}u{\text{'}}_j}\frac{\partial {\dot{u}}_i}{\partial x_j}-{\beta }^{*}\omega k+\\ \frac{\partial }{\partial x_j}\left[\left(\nu +{\nu }_t/{\sigma }_k\right)\frac{\partial k}{\partial x_j}\right]\end{array}$$

（5）

$$\begin{array}{c}\frac{\partial \omega }{\partial t}+{\dot{u}}_j\frac{\partial \omega }{\partial x_j}=P_{\omega }-\beta {\omega }^{\mathrm{2}}+\mathrm{2}(\mathrm{1}-F_{\mathrm{1}})\frac{{\sigma }_{\omega \mathrm{2}}}{\omega }\frac{\partial k}{\partial x_j}\frac{\partial \omega }{\partial x_j}+\\ \frac{\partial }{\partial x_j}\left[\left(\nu +{\nu }_t/{\sigma }_w\right)\frac{\partial k}{\partial x_j}\right]\end{array}$$

（6）

其中： ${\nu }_t=\frac{{\alpha }_{\mathrm{1}}k}{\mathrm{m}\mathrm{a}\mathrm{x}(a_{\mathrm{1}}\omega ,\Omega F_{\mathrm{2}})}$

$$F_{\mathrm{2}}=\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{h}{\left[\mathrm{m}\mathrm{a}\mathrm{x}\left(\mathrm{2}\frac{\sqrt[]{k}}{\mathrm{0.99}\omega y},\frac{\mathrm{500}\mu }{y^{\mathrm{2}}\omega }\right)\right]}^{\mathrm{2}}$$

$$F_{\mathrm{1}}=\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{h}{\left\{\mathrm{m}\mathrm{i}\mathrm{n}\left[\mathrm{m}\mathrm{a}\mathrm{x}\left(\mathrm{2}\frac{\sqrt[]{k}}{\mathrm{0.99}\omega y},\frac{\mathrm{500}\mu }{y^{\mathrm{2}}\omega }\right),\frac{\mathrm{4}{\sigma }_{\omega \mathrm{2}}k}{C{D}_{k\omega }{y}^{\mathrm{2}}}\right]\right\}}^{\mathrm{2}}$$

SST模型考虑了湍流剪切应力的输运，不但能够对各种来流进行准确的预测，还能再各种压力梯度下精确的模拟分离现象，有效的结合了k⁃ω模型在近壁处和k⁃ε模型在外层流场的优势，是当前应用较为广泛的湍流模型之一。

1.2　离散格式

对于不可压缩流体，可以使用有限体积法^［16］在三维控制体（control volume， CV）上建立如下方程：

$${\int }_{\mathrm{C}\mathrm{V}}\underset{\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{p}\mathrm{o}\mathrm{r}\mathrm{a}\mathrm{l}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{m}}{\underbrace{\frac{\partial (\rho \varphi )}{\partial t}\mathrm{d}V}}+{\int }_{\mathrm{C}\mathrm{V}}\underset{\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{v}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{m}}{\underbrace{\mathrm{d}\mathrm{i}\mathrm{v}\left(\rho \varphi \mathrm{u}\right)\mathrm{d}V}}={\int }_{\mathrm{C}\mathrm{V}}\underset{\mathrm{d}\mathrm{i}\mathrm{f}\mathrm{f}\mathrm{u}\mathrm{s}\mathrm{i}\mathrm{o}\mathrm{n}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{m}}{\underbrace{\mathrm{d}\mathrm{i}\mathrm{v}\left(\mathrm{\Gamma }\mathrm{g}\mathrm{r}\mathrm{a}\mathrm{d}\varphi \right)\mathrm{d}V}}+{\int }_{\mathrm{C}\mathrm{V}}\underset{\mathrm{s}\mathrm{o}\mathrm{u}\mathrm{r}\mathrm{c}\mathrm{e}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{m}}{\underbrace{S_{\varphi }\mathrm{d}V}}$$

（7）

这个过程中很重要的一步是将控制体积界面上的物理量及其导数通过节点物理量插值求出，引入插值方式的目的就是为了建立离散方程，这种插值方式常称为离散格式。而在众多物理量的离散中，以速度对流项离散格式的选择最为困难。下面针对本文采用的LU、 LUV、 Localblend、 TVD这4种速度对流格式分别进行介绍。

1.2.1　LU

二阶迎风格式可以描述为：

$${\varphi }_{f,\mathrm{L}\mathrm{U}}=\left\{\begin{array}{c}{\varphi }_P+\left(\nabla {\varphi }_P\right)\cdot {\mathit{d}}_{fP},F\ge \mathrm{0}\\ {\varphi }_N+\left(\nabla {\varphi }_N\right)\cdot {\mathit{d}}_{fN},F\le \mathrm{0}\end{array}\right.$$

（8）

其中：下标P和N分别表示当前网格的相邻网格中心；矢量d_fP和d_fN分别表示网格中心距P点和N点的距离，$\nabla {\varphi }_P$和$\nabla {\varphi }_N$分别代表P点和N点的对应物理量梯度。

1.2.2　LUV

LUV格式由LU格式演化而来，与LU不同的是，它计算一个适用于矢量所有分量的单一限制器，不是为矢量的每个分量计算单独的限制器，这个限制的方向基于梯度变化最大的方向，因此限制的作用会更强，从而保证计算的稳定性。

1.2.3　Localblend

混合格式由一个权重因子$f$决定，其定义如下：

$$\mathrm{混合}\mathrm{格式}=f\times \mathrm{格式}\mathrm{1}+(\mathrm{1}-f)\times \mathrm{格式}\mathrm{2}$$

（9）

其中，权重因子$f$由6个准则决定，分别是mesh non-orthogonality， magnitude of cell centres gradient， convergence rate of residuals， faceWeight， skeness， Co number。具体的表达式如下：

$$f=\mathrm{m}\mathrm{a}\mathrm{x}(f\mathrm{N}\mathrm{o}\mathrm{n},f\mathrm{M}\mathrm{a}\mathrm{g},f\mathrm{C}\mathrm{o}\mathrm{n},f\mathrm{F}\mathrm{a}\mathrm{c},f\mathrm{S}\mathrm{k}\mathrm{e},f\mathrm{C}\mathrm{o})$$

（10）

一般来说，格式1设置为高稳定性格式，格式2设置为高精度格式，从而实现在网格质量较差较差或者收敛性较差的地方使用高稳定性的格式，保证计算结果的稳定性，在其余区域使用高精度格式保证计算结果的准确性。

1.2.4　TVD

在计算的过程中，格式的有界性可以保证计算的准确性和稳定性，但只有一阶迎风格式是严格有界的，因此TVD格式既能够保证有界性，又具有高精度的特点^［16］。其格式的定义如下：

$${\mathrm{\phi }}_e={\mathrm{\phi }}_p+\frac{\mathrm{1}}{\mathrm{2}}\psi (r)({\mathrm{\phi }}_P-{\mathrm{\phi }}_w)$$

（11）

其中：$r=\frac{{\mathrm{\phi }}_E-{\mathrm{\phi }}_P}{{\mathrm{\phi }}_P-{\mathrm{\phi }}_w}$，对于本文选择的Limited linear 1 格式，$\psi \left(r\right)=\mathrm{m}\mathrm{a}\mathrm{x}(\mathrm{m}\mathrm{i}\mathrm{n}\left(\mathrm{2}r,\mathrm{1}\right),\mathrm{0})$。$\psi$下标的3个相邻的网格中心，r可以判断梯度是否存在拐点，当r小于0时，表明梯度存在拐点，此时格式退化为一阶迎风保证有界性；当r大于0时，格式为二阶精度。

2 仿真设置

2.1　模型计算域及网格划分

本文计算的对象是全尺寸35° Ahmed模型，模型特征长度取模型高度H=0.288 m，计算域及边界条件设置如图1所示。计算域尺寸为36H（长）×7H（宽）×5H（高）。车头离入口为8H，车尾离出口为24H，车模两侧离计算域侧面为3H。计算域入口设置为速度入口，来流风速为40 m/s，对应的雷诺数（基于特征长度H）为7.68×10⁵；出口设置为压力出口；计算域侧面及顶面设置为对称边界条件；地面及模型表面为无滑移边界条件。

图1  计算域尺寸

Fig.1  The computational domain

网格划分采用的是OpenFOAM的SnappyHexMesh，网格布局如图2所示。y+控制在30左右，边界层为5层，最终生成3套疏密不同的网格，网格数分别为320万（粗网格），660万（中等网格）和1 310万（细网格）。由于OpenFOAM对网格质量要求较高，故将所有网格正交角控制在75°以下，扭曲率控制在3以下，以保证OpenFOAM计算的稳定性和准确性。

图2  网格分布

Fig.2  The meshing of the entire domain

2.2　求解器设置

流场采用potentialFoam求解器初始化流场，采用simpleFoam 进行迭代求解计算。求解器具体参数设置见表1，其中壁面模型方面的湍动能选择kqR壁面函数，湍流耗散率选择epsilon壁面函数，比耗散率选择omega壁面函数，湍流粘度选择nutUsplading壁面函数。速度对流项的离散格式采用一阶迎风计算2 000步后再分别采用4种不同格式（LU、 LUV、 Localblend、TVD）计算至收敛，以此来保证计算稳定性和收敛性。其他标量对流项则采用一阶迎风计算2 000步后再采用TVD（Limited linear 1）格式， 以保证标量的有界性。湍流模型分别采用SA、RKE和SST这3种进行对比分析，在网格无关性部分只针对LU的速度对流离散格式和RKE湍流模型开展研究。

表1  计算设置

Tab.1  The simulation setting

类别	设置
壁面函数	湍动能： kqR壁面函数
	湍流耗散率：epsilon壁面函数
	比耗散率：omega壁面函数
	湍流粘度：nutUsplading壁面函数
湍流模型	RKE、SA、SST
对流项离散格式	U：一阶迎风—二阶迎风
	k/ε/ω：一阶迎风—Limitedlinear 1
算法	potentialFoam—simpleFoam

3 结果与讨论

3.1　网格无关性验证

表2展示了3套网格数量下的阻力系数C_d值及风洞实验C_d值^［17］（与本文仿真采用的Ahmed模型几何尺寸和雷诺数均相同）。由表可知，相较于实验结果，粗、中、细网格的偏差分别为7.2%、3.6%、2.7%，其中中等网格和密网格的结果较接近，并与实验结果的误差均在5%以内。同时，从中截面的流向速度分布对比来看（见图3），3套网格的速度剖面结果和实验^［18］吻合均较好，彼此间差异较小。综上可知，目前采用的中等网格已经能够足够解析时均流场，为兼顾计算效率，本文选用中等网格进行后续不同湍流模型和离散格式的对比分析。

表2  不同网格阻力系数

Tab.2  The drag coefficient of different mesh

网格	粗网格	中等网格	细网格	实验 ［17］
Cd	0.297	0.289	0.285	0.279

图3  不同网格速度场对比

Fig.3  The velocity U of different mesh

3.2　不同湍流模型的对比

不同湍流模型的C_d值对比见表3。由表可以发现，SA的C_d值最大，远超实验值；其次是SST的C_d值，误差也达到了10%，但均大于实验值；RKE的C_d值与实验最为吻合，误差为3.6%。从中截面速度剖面分布（见图4）来看，在斜背处3个湍流模型的预测相当，均与实验结果^［18］较为接近。而在近尾迹区（x/h=0~0.5），SA预测的更好，在x/h=0.5~1.5区内3个湍流模型预测的都不太准确。

表3  不同湍流模型阻力系数

Tab.3  the drag coefficient of different software

模型	RKE	SA	SST	实验［17］
Cd	0.289	0.354	0.307	0.279

图4  不同湍流模型结果的流向速度

Fig.4  The velocity U of different turbulence model

观察图5所示的中截面湍动能分布图（由于SA模型不包含湍动能项，故无法展示SA模型计算结果的湍动能分布），RKE和SST模型预测的结果仍然较实验偏低很多，但从分布来看，RKE模型更接近实验结果的湍动能分布。

图5  不同湍流模型结果的中截面湍动能

Fig.5  The TKE of different turbulence model in symmetry plane

（a） RKE （b） SST （c） 实验^[18]

由图6所示的中截面速度流线图来看，可以发现SA模型预测的回流区（x≈0.93 H）较其他两个模型更短（x≈1.25 H），这里回流区长度的定义为尾迹流向速度为0的x位置离模型背部的距离，而更短的回流区意味着压力恢复能力弱，因此其阻力值偏大。同时可以观察到，SA模型捕捉的上下分离涡和实验的相当，而RKE的中截面结果上没有明显的下分离涡。SST模型的中截面结果上虽然有分离的涡的存在，但是其位置靠后，与实验偏差相对较大。

图6  不同湍流模型结果的中截面速度分布

Fig.6  The streamlines of different turbulence model in the symmetry plane

(a) RKE (b) SA (c) SST (d) 实验^[18]

此外，图7展示了3个湍流模型计算的Q⁃Criterion=5 000 s^-2的等值面图，由图可以清晰地反映出尾涡的三维时均结构。从图中可以观察到RKE模型的结果中有类似非定常的涡结构，可能是计算结果的震荡造成的。而SA和SST的结果均没有发现类似非定常的涡结构，这表明这两个模型的计算过程中没有发生震荡的情况，且SST模型的计算结果和文献^［2，6］中的SST模型结果基本一致。另外可以观察到，在直背下面的分离涡各有不同，在RKE结果中，这个涡结构在中部断开，因此在RKE的中截面流线图中没有观察到下部的分离涡。而SA模型结果的涡结构则是完全连续且健壮的，SST的这一涡结构更小且更靠后。

图7  不同湍流模型结果的Q （Q=5000 s^-2）等值面

Fig.7  The iso surface (Q=5000 s^-2) of different turbulence model

(a) RKE (b) SA (c) SST

图8展示了这3种湍流模型的背部压图和表面流线图，由图可以发现在SA模型和SST模型均捕捉到了斜背尾端处的二次分离涡，但比前人^［2］用非定常方法计算的二次分离涡尺寸更小，而在RKE的结果中，斜背处并没有产生二次分离涡。

图8  不同湍流模型结果的表面流线图

Fig.8  The surface streamline visualization of different turbulence model

(a) RKE (b) SA (c) SS

从各个湍流模型的诞生背景来看，SA模型是基于经验直接对湍流粘度进行建模，因而结果会有很大的不确定性，虽然其在航空领域得到了广泛的应用，但在汽车空气动力学领域，应用并不广泛。主要是由于SA模型是基于航空领域的工况进行建模和设定的，所以在汽车外流场这种典型三维钝体低速空气动力学的计算中并不完全适用，尤其是其C_d值的预测与实验偏差较大。RKE模型是SKE模型的改进版本，其将湍流粘度计算式中的C_u与应变率联系起来，使得RKE模型可以应对更复杂的分离旋转流动。RKE模型在汽车领域有很广泛的应用，而从上文的结果来看，除了其流场结果出现了震荡，其它的结果是比较合理的，而这一震荡主要是由于速度对流格式造成的（下文将进行详细探究）。SST模型结合了k⁃ω模型在近壁处的优势和k⁃ε模型在远离壁面处的优势，因而也取得了很好的应用，其计算的结果也能和实验较好的吻合。值得一提的是，SA模型和SST均捕捉到了斜背处的二次分离涡，虽然尺寸和前人的结果有所差别，但至少捕捉到了这一流动结构，而RKE模型则没有捕捉到这一结构，这一结果可能得益于SA模型和SST模型在近壁处的优势^［19］。

3.3　速度对流离散格式

速度对流格式的不合理选择有可能造成流场计算的不稳定，因此这部分基于发生震荡的RKE模型，利用OpenFOAM进行LU、LUV、Localblend、TVD （Limited linear 1）速度对流格式和原来的LU格式进行比较，不同速度对流格式的C_d值对比见表4。由表可以观察到，相较于LU格式，其他3种格式的C_d值均有所降低，更接近实验值，误差均小于1.1%。

表4  不同速度对流离散格式的C_d值

Tab.4  The drag coefficient of different velocity convection schemes

格式	LU	LUV	Localblend	TVD	实验［14］
Cd	0.289	0.277	0.280	0.281	0.279

图9展示了LU、Localblend、TVD这3种速度对流格式的流向速度剖面分布和实验结果的对比，由图可以观察到3种格式的速度分布几乎一致，这说明不同格式对中截面的尾迹流场影响相对较小。图10展示了4种格式计算的中截面速度流线图，总体来看流线分布十分相似，但是在下方分离涡的捕捉上略有不同。除了LU格式，其他格式都捕捉到了下方的小分离涡，但尺寸比较小，其中TVD格式捕捉的最为清晰。从Q⁃Criterion=5 000 s^-2等值面（见图11）和表面流线图（见图12）来看，LUV和Localblend这两种格式均能一定程度降低尾迹的震荡，尤其是斜背上方，异常的压力分布都能够被抑制，但是在侧边仍然有这种类似非定常的涡结构，而TVD格式则几乎完全抑制了这一震荡，得到了较为合理的结果。这表明在格式上的妥协在一定程度上能够增强计算稳定性。

图9  不同速度对流格式结果的流向速度

Fig.9  The velocity U of different velocity convection schemes

图10  不同速度对流格式结果的中截面速度分布

Fig.10  Streamlines of different velocity convection schemes in the symmetry plane

(a) LU (b) LUV (c) Localblend (d) TVD

图11  不同速度对流格式结果的Q⁃Criterion=5 000 s^-2等值面

Fig.11  The iso surface of Q⁃Criterion=5 000 s^-2 of different velocity convection schemes

(a) LU (b) LUV (c) Localblend (d) TVD

图12  不同速度对流格式结果的表面流线图

Fig.12  The surface streamline visualization of different velocity convection schemes

(a) LU (b) LUV (c) Localblend (d) TVD

从不同速度对流离散格式的对比来看，使用更稳定的对流格式确实在一定程度上减弱了尾迹中的震荡行为，这得益于局部退化到了一阶迎风格式。但是，造成流场震荡的原因是多方面的，不仅仅有格式的影响，还有网格质量的影响，甚至可能是湍流模型和其他因素综合影响的结果。比如文中在使用同一网格和格式进行不同湍流模型计算的时候，就只有RKE模型出现了这种震荡的行为。而关于时均流场计算中会产生流场震荡的原因，还需后续展开进一步地深入研究。

3.4　不同软件平台对比

不同的软件平台之间的对比均采用相同的求解器设置，湍流模型均选择RKE，对流项采用二阶迎风进行离散。由于在OpenFOAM中的RKE模型采用LU二阶迎风发生了尾迹流场震荡，所以也将其使用TVD格式的结果作为对比参考。最终的不同软件平台的C_d值对比如表5所示。由表5可以看出，使用UV计算时，相较于实验的C_d值，OpenFOAM、FLUENT、 STRACCM+的误差分别为 3.6%、 －1.8%、 －3.2%。而在OpenFOAM中使用TVD计算时，其与实验C_d值吻合最好，误差仅为0.7%。从中截面的速度剖面分布（见图13）来看，3个平台均能取得较为合理的结果，但彼此之间也存在一定差异。在斜背处，OpenFOAM和STARCCM+的计算结果较FLUENT更好；而在直背后部，FLUENT则取得了和实验更为一致的结果。

表5  不同软件阻力系数

Tab.5  the drag coefficient of different software

软件	OpenFOAM（LU）	OpenFOAM（TVD）	FLUENT	STARCCM+	实验［17］
Cd	0.289	0.281	0.274	0.270	0.279

图13  不同软件平台计算的流向速度场

Fig.13  The velocity U of different software

图14展示了不同软件平台的湍动能的结果，可以发现定常计算的结果都较实验结果偏小。其中，OpenFOAM计算的湍动能最大，而FLUENT和STARCCM+的计算结果则较低。图15展示了中截面的速度流线图，结果显示定常结果的回流区都较实验的偏长，这与之前的研究结果^［6-7］是一致的。其中，OpenFOAM和STARCCM+的回流区较FLUENT的计算结果更长。另外，FLUENT计算结果的斜背后方分离涡位置更靠下，而且可以清晰观察到下方分离涡的存在，STARCCM+的结果中也可以观察到这一现象，但是这个涡结构较FLUENT的更小。OpenFOAM采用LU计算的流场震荡结果中则没有观察到这一结果，但是使用TVD格式计算的结果中则捕捉到了这一结构。

图14  不同软件结果的湍动能

Fig.14  The TKE of different software

(a) OpenFOAM (b) FLUENT (c) STARCCM+ (d) 实验^[18]

图15  不同软件结果的中截面速度分布对比

Fig.15  The streamlines of different software in the symmetry plane

(a) OpenFOAM (二阶迎风) (b) OpenFOAM (TVD)

(c) FLUENT (d) STARCCM+

图16展示了不同软件平台计算得到的尾部 Q⁃Criterion= 5 000 s^-2的等值面图，由图可以观察到在OpenFOAM的LU格式计算结果中存在类似非定常的涡结构，这是由上文提到的尾迹流场震荡引起，也是造成其阻力系数预测偏高和下方分离涡消失的原因。而采用TVD格式时，OpenFOAM计算得到了比较合理的三维时均尾迹流场结构。对于FLUENT的计算结果，从图中可以看到，尾迹两侧预测出了明显的流向C柱涡，这与之前文献中的仿真结果^［2］有一定差异。之前的实验研究^［17-18，20-21］已经表明，流动在35° Ahmed模型的斜背处完全分离，同时流向C柱涡基本消失或者非常微弱，因此FLUENT的计算明显高估了流向C柱涡的强度。STARCCM+的结果比较合理，和OpenFOAM采用TVD计算的结果较为接近，和文献［2］中定常计算的结果吻合较好。

图16  不同软件结果的Q-Criterion=5 000 s^-2等值面

Fig.16  The iso surface of Q-Criterion=5 000 s^-2 of different software

(a) OpenFOAM (二阶迎风) (b) OpenFOAM (TVD)

(c) FLUENT (d) STARCCM+

图17展示了不同软件平台计算计算的背压和表面流线图，由图可见OpenFOAM采用LU格式计算的斜背上端左右角落的流线几乎平行于斜背上切线，而且背部有异常的压力分布，这正是在涡等值面结果（图16（a））中观察到的类似非定常涡结构造成的，而使用TVD计算的结果则没有观察到这一结果。在FLUENT计算结果中，斜背角落处有较大的负压区，这是由于其预测出的流向C柱涡诱导产生的。

图17  不同软件结果的表面流线图

Fig.17  The surface streamline visualization of different software

(a) OpenFOAM (二阶迎风) (b) OpenFOAM (TVD)

(c) FLUENT (d) STARCCM+

4 结论

本文基于OpenFOAM开源平台，针对35° Ahmed汽车标准模型，研究了RKE、 SA、SST这3种湍流模型，以及LU、 LUV、 Localblend、 TVD这4种速度对流离散格式在Ahmed外流场计算中的表现。同时，还对比了OpenFOAM、 FLUENT、STARCCM+这3个软件平台的差异，得到的主要结论如下：

（1） 在湍流模型方面，对于C_d值，RKE模型预测结果和实验结果最接近，其次是SST模型，而SA模型的误差较大。针对尾迹流场方面，在斜背气流分离处，3个湍流模型的预测相当，均和实验结果较为接近；而在远离模型背部的尾迹区域，3个湍流模型预测和实验均有一定的偏差。

（2） 在速度对流离散格式方面，LUV和Localblend等更稳定的格式可以一定程度上减弱计算不稳定性带来的尾迹振荡问题，但只有TVD格式能够完全的抑制这种震荡行为。

（3） 3种流体仿真软件（OpenFOAM、 FLUENT、 STARCCM+）计算的C_d值均与实验有一定的吻合度，其中OpenFOAM采用TVD格式的误差最小，仅为0.7%。此外，在流场方面，OpenFOAM和STARCCM+可以得到较为合理和相近的三维尾迹流场结构，而FLUENT高估了流向C柱涡的强度。