雷诺数为22000的二维方柱绕流仿真计算

图 15-1-3 网格划分结果
图 15-1-3 给出了网格划分结果，在划分网格时，约束整个面的网格尺寸为 10mm， 并且定义入口、 出口、 上下边界、 方柱表面的名称集 （Named Selections） ， 以方便在 Fluent 中定义边界条件。
图 15-1-4 局部网格示意
图 15-1-4 给出了局部网格示意图，从图中可以看出，在方柱表面的网格高 度是较高的， 使用尺度化近壁面函数法时并没有划分很细的边界层网格。最后共 得到 48900 个网格，网格最大长宽比为 1.0069，满足仿真计算的网格要求。
（3）对于存在较高导热率的共轭传热问题； （4）对于流场网格长宽比较大的问题。 在本算例中，虽然不存在以上问题，但仍采用双精度 4 核并行计算求解。 在进行计算时， 先进行稳态计算，然后将稳态计算的结果作为瞬态计算的初 场，再进行瞬态计算。 求解器采用压力基求解器，并进行湍流模型的设置。
图 15-1-6 湍流模型的选择
关于雷诺数 22000 方柱绕流 （关于壁面函数法、湍流粘度比的计算等概念可以参考 文献十六）
15.1.1 RNG 模型与尺度化壁面函数在方柱绕流中的应用
图 15-1-1 分析流程
图 15-1-1 给出了二维方柱绕流问题的分析流程， 在 Workbench 中建立模型， 然后使用 Mesh 模块进行分网，最后提交到 Fluent 中进行计算。
图 15-1-15 方柱表面的 Y+值分布范围
仿真后可以知道，近壁面 Y+值的分布范围为 7.07-68.6，在 Fluent 的官方 文献中推荐壁面函数法的 Y+值适用范围为 30-500，对于尺度化壁面函数，由于 在近壁面使用了预防 Y+值过小导致结果恶化的技术，其 Y+值范围下限可以适当 放宽， 但仍应尽量取在 11.25 以上。 过大或过小的 Y+值都可能引起结果的恶化。
图 15-1-5 采用 4 核并行计算
在打开 Fluent 时，按照图 15-1-5 所示，采用单精度、4 核并行计算。双精 度一般用在以下几个方面： （1）存在明显的细长比很大的几何特征，如细长的管道； （2）流体是几个区域，区域之间用小尺寸的管道连接，其中一个区域的流 体压力大大高于整个流域的平均水平；
（在数学上即自由滑移边界条件） ；当流场发展到出口位置时，认为流动已经充 分发展，因此设置出口为自由流 outflow；入口在用速度入口边界条件，入口速 度、湍流强度、湍流长度尺寸均采用算例中所给设定。 在设置参考值时首先选定 compute from->流场名，这样就正确的给出了参 考值中的物理参数（密度、粘度等） ，然后指定 Area 的值为 0.1。 关于参考值的设定是非常有讲究的，根据文献十八 P173 的升阻力系数定义 式： Cd = Fd 1 2 2 ρv S Fl 1 2 2 ρv S
其中，St 为斯特劳哈尔数，f 为漩涡脱落频率，v 为来流速度。根据文献十 七中的实验结果（St=0.13） ，可以反推出漩涡脱落频率应为 3.5s，所以在选择 瞬态步长时就有了选择的依据， 瞬态求解采用二阶隐式求解方式，针对方柱绕流 问题，隐式方法对时间步长的要求不高，很多人取漩涡脱落周期的 1/50 作为瞬 态计算的时间步长，本文取为 0.02s，可以充分地描述出漩涡脱落的瞬态细节。 对于收敛残差的取值， 很多人是取稳态残差稳定值的 1/50， 本文取为 1/100， 即 0.00001。指定每一时间步最多计算 25 个子步，共计算 5000 个时间步，共计 100s。
图 15-1-12 方柱绕流流场的静压云图与动压运图
根据图 15-1-12 中可以看出， 在方柱后方出现了明显的周期性漩涡脱落的现 象。漩涡脱落的现象在总压云图与速度云图中也可以很明显地看出来。 文献十七中认为 RNG 方法在预测平均阻力系数和回流区再附着距离时与实 验及大涡结果比较接近，但是在预测升阻力系数的均方根值时仍具有一定的差 距。采用雷诺平均的 RNG 模型，从模型精度上来说要比大涡模拟至少低一个数 量级。
图 15-1-10 方柱绕流的稳态升阻力迭代历史
从图 15-1-10 中可以看出， 升阻力系数在稳态计算中最终归于稳定，这与残 差收敛至极小值是相互呼应的。 而在有些绕流问题的稳态计算中，也会出现瞬态 特征。如图 15-1-11 所示：
图 15-1-11 某圆柱绕流的稳态计算
图 15-1-11 中给出了某圆柱绕流的稳态计算结果， 在稳态计算结果中出现了 瞬态特征。 在稳态计算中出现瞬态特征的现象，在方柱绕流与圆柱绕流的现象中
15.1 RNG 湍流模型在 22000 雷诺数二维方柱绕流中的应用 【算例内容】 ： 方柱尺寸为 100mm×100mm， 流场区域为 1400mm×3500mm， 方柱距离入口 500mm。流体工质采用水，其密度为 998.2 千 克每立方米，动力粘度为 0.001003 帕秒。 由雷诺数为 22000 可以折算出入口流体的速度应为 0.2210579 米每秒。实验所用来流湍流强度为 2%，取入口湍 流长度尺度为 0.007（入口湍流长度尺度取值参考文献十六 P215） 。 【仿真目标】 ： 对于雷诺数为 22000 的亚临界方柱绕流，合理选择湍流 模型将是关键。 对于分离流动， 不适合采用标准 k-ε 湍流模 型，在高雷诺数湍流模型族中，RNG 湍流模型是值得研究的 选择。而 SST k-ω 模型、雷诺应力模型、分离涡模型和大涡 模拟也是适合的，本节仅讨论 RNG 湍流模型。 要求用两种近壁面处理方式处理方柱表面，第一种方式 是采用尺度化壁面函数，第二种是采用加强型近壁面处理方 法。计算完成后比较两种近壁面处理方法的结果。 文献十七中给出了雷诺数为 22000 的方柱绕流的实验结 果、大涡模拟结果和 RNG 计算结果。可以用来验证本节计算 结果。
需要注意的是，在计算升力与阻力时，升力方向垂直于来流方向（不是垂 直于机翼！ ） ，阻力方向平行于来流方向（不是平行于机翼！ ） 。仅当攻角为 0°或 180°时，升力/阻力才垂直/平行于来流方向。 在 Fluent 二维计算中，机翼长度方向（展向）垂直于几何平面，软件在处 理时机翼长度置为 1 （在流体力学中处理二维问题， 展向长度全部取单位长度） 。 其余参考值需要参考升力系数公式进行计算，并且一定要注意升力/阻力的方向。 对于本二维算例，S 的取值与弦长数值相等，因此本算例 S 取值为 0.1，即 设定参考值时，与尺寸相关的进设置 Area 为 0.1 平方米即可，其余的保持默认 不变；密度等材料的参考值已经通过 Compute from 功能实现了。 之后设置求解器选项， 选择耦合求解器， 对流项的处理采用二阶迎风格式 （对 于之前划分的网格，也可采用更高阶的 QUICK 格式） 。保持默认的收敛残差判定 标准 0.001 不变， 监测了升阻力系数随迭代步数的变化曲线。 以入口条件初始化 整个流场区域，设定最大迭代步数为 1000 步。
图 15-1-2 流场模型
按照算例内容建立流场模型如图 15-1-2 所示，然后进行分网，由于采用的 是尺度化壁面函数法，所以近壁面网格没必要划分地太细。 进入到 Mesh 模块后首先到 mesh 模型树下将 Physics Preference 改为 CFD 与 Fluent，然后设置 Sizing 下的网格相关度、网格光顺等参数。
图 15-1-7 机翼攻角示意图
升力系数：物体（如飞机、导弹）所受到的升力与气流动压和参考面积的
乘积之比，是一个无量纲量。 阻力系数：对于飞行器来说，阻力系数定义为物体（如飞机、导弹）所受 到的阻力与气流动压和参考面积之比，是一个无量纲量。 关于升力系数与阻力系数的计算参考图 15-1-8：
图 15-1-8 计算升阻力系数所用的参考系
图 15-1-9 稳态计算的残差曲线
图 15-1-9 给出了稳态计算残差曲线，从曲线中可以看出，残差曲线基本上 随着计算会继续趋于收敛， 其中连续性方程的残差最大。在计算此类绕流问题的 实际应用过程中， 还会出现当稳态计算的残差位于高位稳定的情况，出现这种情 况时， 应当考虑的一种可能存在的情况就是物理模型中存在瞬态变化的物理效应， 导致稳态迭代时前后迭代的范数残差不能够收敛而位于高位稳定， 这种现象在用 稳态雷诺平均方法求解瞬态问题时经常出现， 此时可以将稳态计算结果作为瞬态 计算的初场，继续计算。所以单看残差收敛情况不一定有太大意义，要结合实际 情况来看。另外发现当启用 Differential Viscosity Model 选项后，本来收敛 至非常小的残差容易转变为高位稳定甚至发散。 本例中升力系数与阻力系数的迭代历史曲线如图 15-1-10 所示：
图 15-1-13 方柱绕流流场的总压云图与速度云图
图 15-1-14 给出了方柱升阻力系数的历史曲线。
图 15-1-14 方柱绕流的升力与阻力系数的历史曲线
从图中可以明显看出， 方柱绕流的升阻力系数随着时间的推移呈周期性的波 动。经过对数据的统计，可以知道，方柱阻力系数的平均值为 2.3（取 80-90s 之间的每个时间步的阻力系数做平均） ，实验值为 2.1，误差为 9.5%，计算出的 漩涡脱落周期为 3s，由此所得的斯特劳哈尔数为 0.15，实验值为 0差曲线将会居于高位，并且升力与阻力系数 的值不一定准确，比如上图中升力系数平均值部位 0。 还有一种现象是升阻力系数趋于平稳，但残差曲线局于高位，造成这种现象 的一种原因是启用了 Differential Viscosity Model 选项，这也说明了 Differential Viscosity Model 选项对计算所造成的不稳定性。 根据文献十八 P119 的关于斯特劳哈尔数的公式： St = fD v
Cl =
以上两式中， Cd 代表阻力系数， Cl 代表升力系数， Fd 为钝体所受到的阻力， Fl 为钝体所受到的升力，ρ 为流体的密度，v 为前方来流的速度，S 代表钝体的 投影面积，在二维问题中 S 的数值与弦长相等（弦长乘以一，二维问题中展向长 度为单位长度一） 。 有关翼型的概念： 翼型：飞机机翼具有独特的剖面，其横断面（横向剖面）的形状称为翼型， 称为翼型 前缘：翼型最前面的一点。 后缘：翼型最后面的一点。 翼弦：前缘与后缘的连线。 弦长：前后缘的距离称为弦长。如果机翼平面形状不是长方形，一般在参 数计算时采用制造商指定位置的弦长或平均弦长 攻角：也称迎角。对于翼形来说，攻角定义为翼弦与来流速度之间的夹角， 抬头为正，低头为负，常用符号 α 表示。其示意图见图 15-1-7：
如图 15-1-6 所示，选择 RNG 湍流模型与尺度化壁面函数处理方法，并且不 启用 Differential Viscosity Model 选项，在 ANSYS Fluent 14.5 的帮助文件 中对该模式有以下论述： The RNG turbulence model in ANSYS FLUENT low-Reynolds-number effects. 可见在考虑到边界层的低雷诺数影响时应该打开该选项，另外，有人在做空 气的涡街仿真时遇到过不勾选该选项就得不到涡街的情况。 （/bbs/viewthread.php?tid=76599） 但是在实际计算中发现，启用该选项后，计算容易变得不稳定，出现了数值 波动甚至是不收敛。因此，在本算例中，不启用该选项。 流线曲率修正是因为使用各向同性涡粘假设的湍流模型对流线曲率不敏感， 因此在流线弯曲时需要进行流线曲率修正，启用该选项。 选择 Fluent 中的 water-liquid 作为流体材料，然后进行边界条件的设置： 设置方柱表面为无滑移壁面条件 wall；流场上下边界为对称边界条件 symmetry has an option of using a differential formula for the effective viscosity to account for the