沿程损失阻力系数的FLUENT数值模拟(计算流体力学作业)

计算流体力学课程作业

作业题目：

学生姓名：

沿程损失阻力系数的 FLUENT数值模拟

易鹏

动力工程及工程热物理12级

机械与运载工程学院

学生学号： 专业年级： 学院名称：

2012年5月2日

沿程损失阻力系数的 FLUENT数值模拟

一、 引言

沿程损失（pipeline friction loss）是指管道内径不变的情况下，管内流体流过一段距离后的水头损失。其中边界对水流的阻力是产生水头损失的外因，液体的粘滞性是产生水头损失的内因，也是根

lv2本原因。沿程能量损失的计算公式是：hf=λ。其中：l为管长，d2gv2为沿程损失系数，d为管道内径，为单位重力流体的动压头（速度

2g水头），v为流体的运动粘度系数。粘性流体在管道中流动时，呈现出两种流动状态，管道中的流速vvcr(vcr为层流向湍流转变的临界流速)为层流，此时整个流场呈一簇互相平行的流线。则vvcr时为湍流，流场中的流体质点作复杂的无规则的运动。沿程损失与流动状态有关，故计算各种流体通道的沿程损失，必须首先判别流体的流动状态。

沿程损失能量损失的计算公式由带粘性的伯努利方程

v12p1v22p2++z1=++z2+hf推出，可知，h=P1-P2 其中：

f2gρg2gρgρgv2 ——单位质量流体的动能（速度水头）。流体静止时为0。 2g——单位质量流体的势能（位置水头）。 zpρg——单位质量流体的压力能（压强水头）。

1

ΔpL=λ ，计算出达西摩擦因子又由量纲分析的定理，得出 1dρV222LVμVd2Δpdh=λλ=,由于Re=，ν=，则λ=f(Red)。 2， 则fρνD2gLρV关于沿程损失最著名的是尼古拉茨在1932～ 1933年问所做的实验(右图为实验装置图)。其测得曲线如图1，从此得出了几个重要结论：

1．层流区 Re＜2320为层流区。在该区域内，管壁的相对粗糙度对沿程损失系数没有影响。

2．过渡区 2320＜Re＜4000为由层流向湍流的转换区，可能是层流，也可能是湍流，实验数据分散，无一定规律。

3．湍流光滑管区 4000＜Re＜26．98（d/ε）8/7，为湍流光滑管区。勃拉修斯（p.Blasius）1911年用解析方法证明了该区沿程损失系数与相对粗糙度无关，只与雷诺数有关，并借助量纲分析得出了4×10e3＜Re＜10e5范围内的勃拉休斯的计算公式为

0.31640.25

Re 2

湍流光滑管的沿程损失系数也可按卡门一普朗特（Karmn-Prandtl）公式

11/2进行计算。

2lg(Re1/2)0.8

当105＜Re＜3×106时，尼古拉兹的计算公式为

0.00320.221Re0.237

4．湍流粗糙管过渡区 26．98（d/ε）8/7＜Re＜2308（d/ε）0.85为湍流粗糙管过渡区。该区域的沿程损失系数与按洛巴耶夫（Б.H.Лo6aeв）的公式进行计算，即

2dV1.42lgRe1.42lg1.273qv 5．湍流粗糙管平方阻力区 2308(d/ε)0.8511/22lgd1.74 2进行计算。

3

图1尼古拉茨曲线

三、数值模拟

1、前处理

因为层流有精确解所以在此不做讨论，而湍流状态下如果考虑圆管的粗糙度也是十分复杂，而且在粗糙表面的流动很难模拟，所以在此我们重点研究湍流水力光滑区的达西摩擦因子与Re的关系。FLUENT中通过改变流速或者粘度系数来控制Re，并进行数值模拟，计算出管中试验段两端的压力的差值，即可得到沿程损失阻力系数（达西摩擦因子），再将所得的值与上图水力光滑区曲线或布拉休斯公式对比，判断其是否正确。

1．模型 建立一个半径r=21mm，长l=3m的圆截面直管，其中前2m是前置段，用来让湍流充分发展，后1m为实验段。假设

4

其材料是光滑的，没有摩擦，内部流体为水。设水的ρ为1000kg/m3，粘度系数μ为0.001kg/(m*s)。下图就是使用K-epsilon湍流模式试算的velocity inlet后端vy云图，说明在试验段之前设置前置段还是十分有必要的。

使用gambit可以很容易的建模，直接使用cylinder命令建立方案中的模型。但是在此未使用这种方法，由于液体的粘性力作用，在壁面附近有比较大的速度梯度，而且在入口端是湍流发展段，所以需要端面使用边界层网格加密，轴向在入口处加密。具体步骤是：

1．做半径为0.021的圆。

2．做出x=0.021，y=0，z=5的点，并连接圆上与其对应的两点。

3．为该线mesh，选择ratio 1.05，让线网格在入口处加密。在此同时将将入口端面的圆分成50等分线网格（数目自定，但是这样已经足够）

4．使用sweep命令，选上with mesh选项，让直线绕圆周旋

5

转成圆柱面，并且将网格自动画好。如右图。

5.端面上创建边界层网格，first percentage(第一层边界层网格的高度关于宽度的百分比)在这里取了15,rows取5层，Growth factor取1.1。(注意一个问题，就是在画边界层网格时有个方向选择问题，打开edge的list里面，每个edge其实可以点多次，具体多少次看该edge属于多少个face，通过试验，就可以看到边界层具体会向哪个方向生成)。具体设置如右图。

6.为端面直接画面网格，由于之前 端面的圆已经分好了网格和边界层网格，不用设定参数gambit自动画网格，完成后如下图。

7．在生成体的选项中选择sweep，

勾选with mesh选项，让圆端面沿管轴线方向扫过，即可完成体网格的绘制。

8．最后选择求解器（solver）Fluent 5/6，设置z=0处端面为Velocity in、圆柱面为wall和z=3处端面为outflow。

6

9．Export mesh。注意：不要选择2D模型输出的选项。

下图是网格完成后的模型。一共生成了50800个体网格。

四、数值模拟及数据处理

由于是光滑圆管（或水力光滑），则达西摩擦因子λ只是Re的函数。而在Re<2000时，圆管中的流动属于层流，泊肃叶也做过此范围内流动的大量实验，得出经验公式λ=dV64，Red定义为，在这里

νRedV为距入口10m（即试验段的起始端）的截面平均流速，湍流时的Re也如此定义。又因为在圆管流动中雷诺数Re>2000才进入湍流状态，并且在20007

使用ANSYS12.0中的Fluent作为流场模拟的软件，在这里圆管属于细长结构中的流动用双精度（Double Precision）模式模拟较精确。准备使用k-epsilon，增强壁面函数的k-epsilon和S-A湍流模式分别计算。而且由于流动是湍流，并且网格在内部并不是和流速垂直的，所以使用二阶迎风格式，来提高精度，并且设置残差到10e-5，以提高精度。通过Fluent的Report菜单中的surface Integrals命令可以获得入口和出口的压力和速度的平均值。下图为Surface Integrals的窗口，其中Inexp是实验段的入口截面，out就是出截面：

Rein为入口雷诺数，v为入口速度，p1是试验段起始端的压力，p2是试验段结束端的压力。Red为实验段起始处雷诺数。1和2分别为模拟算出的达西摩擦因子和用布拉休斯公式算出的达西摩擦因子。

8

再通过此表数据作出拟合曲线与布拉休斯公式的解对比，分析误差。 K-epsilon湍流模式计算结果

Rein 3500 4000 4500 5000 6000 7000 9000 12000 15000 20000 ν 0.0833 0.0952 0.1071 0.1190 0.1428 0.1666 0.2142 0.2856 0.3570 0.4760 p1 -22.77689 -27.525547 -32.735207 -38.257824 -50.02042 -62.615845 -90.401306 -136.95033 -188.45871 -282.94147 p2 1 2 误差 -48.133842 0.1023 0.0411 148.74% 0.0386 136.74% 0.0376 131.67% 0.0359 122.82% 0.0346 115.23% 0.0325 102.65% 0.0302 0.0286 0.0266 87.95% 76.35% 61.26% -58.737934 0.09643 0.0398 142.38% -70.20034 0.0915 -82.344009 0.0872 -108.3573 -136.4162 0.0801 0.0745 -198.27414 0.0658 -302.46423 0.0568 -417.95313 0.0504 -630.13373 0.0429 可见误差相当之大，究其原因，应该是标准k-epsilon在壁面区使用了不够精确的近壁函数的半经验公式，以及工况中流场为层流向湍流的过度区。在FLUENT中对K-epsilon做如下修改：

增

9

强壁面函数的K-epsilon湍流模式计算结果

Rein 3500 4000 4500 5000 6000 7000 9000 12000 15000 20000 ν 0.0833 0.0952 0.1071 0.119 0.1428 0.1666 0.2142 0.2856 0.357 0.476 p1 p2 1 2 误差 24.77% 22.15% 20.11% 18.51% 16.27% 15.06% 14.49% 15.31% 17.55% 13.79% -10.265804 -14.505714 0.051327 0.041136 -12.627173 -17.870281 0.048595 0.039785 -15.308632 -21.644625 0.0464 0.038631 0.037627 0.03595 -18.074289 -25.591486 0.04459 -24.361353 -34.508709 0.0418 -31.505663 -44.656143 0.039799 0.034591 -48.221485 -68.535301 0.03719 -79.798485 -113.64851 0.03486 -119.07154 -170.06151 0.033607 0.032484 0.03023 0.02859 -190.21982 -271.87988 0.030274 0.026606 对比标准k-epsilon的精度高多了，但是仍然不够精确。如果将网格划分得更精细些，将更好的控制误差。 S-A湍流模式计算

Rein 3500 4000 4500 5000 6000 7000 9000 12000 15000 20000 ν 0.0833 0.0952 0.1071 0.119 0.1428 0.1666 0.2142 0.2856 0.357 0.476 p1 p2 1 2 误差 5.00% 3.73% 2.84% 2.20% 1.37% 0.90% 0.97% 0.22% 0.17% 2.13% -8.285567 -18.98961 0.043193 0.041136 -10.27976 -23.63747 0.041268 0.039785 -12.46471 -28.73887 0.039726 0.038631 -14.83306 -34.28072 0.038453 0.037627 -20.10822 -46.64936 0.036444 0.03595 -26.07841 -60.67434 0.034901 0.034591 -40.20676 -93.95073 0.032798 0.032484 -65.64599 -153.9034 0.030297 -96.69733 -227.0509 0.028638 0.03023 0.02859 -162.0646 -381.9362 0.027171 0.026606 10

以上是两种湍流模式模拟的曲线图，很明显S-A模式模拟出的结果优于k-epslion模式的。

五、总结

虽然模拟的结果和经验公式还是有误差，但是经验公式本身也是不精确的，而且由于湍流模式和数值计算中都有不可避免的误差，然而能做到5%以内已经满足工程上的需要了。而为什么S-A这种适用于低雷诺数的湍流模式，反而在高雷诺数情况下的计算结果更精于k-ε，通过查找文献，有以下的解释，一个湍流模型要想精确地求解出流动阻力，必须考虑近壁区低雷诺数的影响，特别是必须能很好地模拟出近壁区的时均速度轮廓。Spalart—Allmaras湍流模型可很好地满足上述要求，要知道Spalart—Allmaras湍流模式最早是用于解决飞行器阻力问题，而且结果也证明Spalart-Allmaras湍流模型在模拟流动阻力方面的优势，但是S-A模式的计算成本也较高。标准的k-ε

11

模型在近壁区采用壁面函数的半经验果公式，其误差较大。如使用修正的k-ε模式，即在近壁区的壁面函数做了修正，通过对比修正k-ε的确可以提高精度，如果使用了更精细网格会有更不错的结果。

12