好尴尬给哥哥
范文一：fluent自然对流关键设置
fluent自然对流模拟关键点
关于fluent做自然对流的数值模拟，与强制对流的模拟有很大的不同，关键点是自然对流的驱动力是由于温差引起的密度差，进而在重力的作用下，引起流体产生运动。这跟强制对流需要由外界提供动力是完全不同的，所以其设置也是不同的，现把我的一些经验和大家分享。
1， 湍流模型的选择
对于自然对流，湍流模型的选择也是不同的，其主要是要强调壁温的影响和浮升力的影响。具体设置如下:
选择了k-e模型，然后是选择了enhanced wall treatment和full buoyancy effects选项，也就是强调壁温作用和浮升力的作用。
2， 能量方程和重力项都要打开。
3， 材料的设置
(1) 密度采用Boussinesq假设，然后需要设置流体的密度为一个定值，还要设置流体的
热膨胀系数。关于流体的热膨胀系数需要查资料了。如下
(2) 选用不可压缩理想流体假设
就是流体按不可压缩理想流体对待，其含义是，流体的密度是变化的，其变化是由温度变化引起的，而不是由压力变化引起的，如下所示，流体密度项选择incompressible ideal gas选项。
(3) 对上述两种选择的解释
首先，Boussinesq假设比incompressible ideal gas假设，更易收敛。一般情况下选择这个假设即可。
其次，对于Boussinesq假设有其适用范围，因为其假设流体密度是个定值，所以其主要用于流体密度变化小于20%的情况，也就是两壁面温差较小的情况，那么温差值有没有一个经验数据呢,有的，一般如果壁温与流体温度相差在200K以上，Boussinesq假设就不适用了。 第三，incompressible ideal gas假设，其把流体密度看做随温度的变化而变化，所以其适用范围较广，对于Boussinesq假设不能适用的，就用这个假设。
4， solution method设置
对于求解方法的设置，主要是对压力离散方法的选择，要选择PRESTO!或body force weighted
选项。如下所示
以上就是自然对流的关键设置，按照上述设置来进行模拟，完全可以得出浮升力的漩涡流，完毕～
范文二：Fluent实例：辐射与自然对流模拟
Fluent 辐射与自然对流模拟
引言
在这个算例中, 将会解决二维方箱中的辐射与自然对流相结合的问题, 网格采用 四边形单元网格。
在这个算例中将会学到以下知识点:
1. 应用 Fluent 中各种辐射模型 Rosseland ;
2. 使用 Boussinesq model定义密度;
3. 设定辐射与自然对流传热问题的边界条件;
4. 将单一的墙划分为多个墙区域;
5. 对已有的流体物性进行修改;
6. 用隔离求解器求解;
7. 显示速度矢量和流函数等值线,以及温度等值线。
问题描述
将被考虑的问题如图 5.1 所示,一个边长为 L 的正方形箱体,右墙温度为 2000K, 左墙温度为 1000K, 上下墙绝热, 重力向下, 由于热重引起密度梯度所以发展为浮 力流。 箱体中的介质被认为是有吸收性和散射性的, 因此墙壁间的辐射交换因存 在吸收被减弱,同时也因为介质的散射作用而增强了。所有墙壁被认为是黑体, 目的在于应用有效的辐射模型计算箱体中流场和温度场分布,以及墙壁的热流 量,并且对于不同光学深度 aL 比较所表现出的特性。
工质普朗特数大约为 0.71, 基于 L 的雷诺数为 500000, 这说明流动相当于层状流 动,应用 Boussinesq 假设来模拟浮力流动。普朗克数为 0.02,用于考虑传导与辐 射的相对重要性, 其中, T 0 = ( T h + T c )/2。 在这个算例中将有三种 optical thickness 的情况会被考虑到,分别是 aL=0, aL=0.2, and aL=5。注意:物理属性和工作条 件 (重力加速度) 都已经给定以适合于产生的想要的普朗特数, 雷诺数和普朗克
数。如下图所示:
第 1 步: 网格
将网格文件rad/rad.msh 拷至fluent 的工作目录下(就像在指南1 中描述 的一样),并起动fluent 的二维单精度解算器。 1. 读取网格文件rad.msh .
File Read Case...当网格读入的时候,在Fluent 控制窗口会显示相应的信息, 会报告网格有2500 个单元。 2. 检查网格质量。
Grid Check... Fluent 会对网格进行各种各样的检查,并会在控制窗口显示信 息。特别注意最小体积,确保它是正数值。
3. 显示出网格(如图5.2)。
图 5.2:网格显示
注意 :此时所有的墙体为一个整体, wall-4。 下面需要将其分成四个独立的区域,
才能对不同的区域定义不同的边界条件。
4. 将单一的墙分为四部分。
Grid Separate Faces...
(a ) 在 Options . 下选择 Angle 的分离方法;
(b ) 在 Zones 列表下选择 wall-4;
(c ) 定义角度值 89 度。
(d ) 单击 Separate 。
带法向矢量的面在被 89度分隔后会成为单独的区域,既然四个墙的区域是正交 的, wall-4 就被分为四个独立的区域。
5. 重新显示出网格。
(a ) 在 grid display 对话框中选定所有的面后单击 Display 。
注意:现在有四个不同的墙区域而不是一个。
另外:你可以用鼠标右键来验证墙区域名称所对应的墙。 在图形窗口中用鼠标右 键点击任何一个边界, 区域名称, 类型都在 Fluent 控制窗口中显示出来。 当同一 个类型有几个不同区域, 而你又需要迅速的区分开来时, 这个用途的作用就特别 明显。有些时候,可以使内部网格不显示,以便更准确的选择边界。
第 2 步:模型选择
辐射模型采用的是 Rosseland model。
1.
选择默认解算器。
2. 打开Rosseland 辐射模型。
Define Models Radiation...
当你在Radiation Model 面板上单击ok 时,会出现一个信息提示框,告诉你新 的材料物性被添加了,你将在后面设置物性参数,因此现在只需单击ok 确认这 个信息即可。
注意:当你激活辐射模型时, fluent会自动打开能量解算器, 因此不用另外再激 活Energy对话框了。
3.在模型中添加重力影响。
Define Operating Conditions...
(1) 打开 Gravity 。
面板会自动扩展,显示所要额外的输入值。
(2) 设置 y 方向的 Gravitational Acceleration 值为-6.94e-5 m/s 2。
上文中提到过, 重力加速度调整到产生合适的无因次数 (普朗特数, 雷诺数和普 朗克数),见图 5.1 及相关的说明。 (3) 设置工作温度为 1000k.
在下一步要激活的 the Boussinesq model 中将会用到工作温度值。
第 3 步: 材料物性
默认的流体材料为空气, 在此次要解决的问题中工作介质也正是空气。 但却是假 想的流动, 使各个物性都适合于给出特定的无因次数, 因此必须对默认的空气参 数进行修改。应用光学深度 aL=0.2 进行计算。(假设 L=1,吸收系数则应设为 0.2 )。在后面, aL=5 或者无光学介质(aL=0)时,计算得到的结果将会用来 比较不同的光学介质对于的辐射模型的不同。
Define Materials...
1. 在 Density 下 拉 框 中 选 择 boussinesq , 并 设 定 密 度 值 为 1000kg/m3; 关 于 boussinesq 模型的详细内容请看 User's Guide。
2. 设置比热容, Cp=1.103e4 J/kg-K. 3. 设置导热系数为15.309 W/m-K.; 4. 设置粘度为0.001 kg/m-s.; 5. 设置吸收系数为0.2 m -1. ;
提示:通过滚动条使先前面板中不可见的物性显示出来。
6. 在 Scattering Coefficient 和 Scattering Phase Function 中保持默认值, 在要解决的问题中不涉及到散射问题;
7. 设定热扩散系数(用 boussinesq 模型时)为 1e-5K -1。 8. 单击 Change/Create,关闭 Materials 面板。
第 4 步: 边界条件
Define Boundary Conditions...
1. 对底部的墙设置边界条件(wall-4.006) 。
注意:如果不能确定 wall-4.006 对应的是否为底墙, 可以用鼠标右键在图形窗口 中单击底墙,在边界条件面板中 Zone 列表中会自动显示相应区域的名称。 在设定其它墙的边界条件时也可以这样做,一定要保证设定正确的边界条件。
(1)将 Zone Name 改为bottom;
(2)对于绝热墙,保持默认设置热边界条件。(heat flux=0)
注意:Rosseland model 模型不需要设定墙的发射率,对于其它辐射模型,就需 要设定墙的发射率。
2. 对左墙设置边界条件(wall-4)
(1) . 将 Zone Name 改名为 left.
(2)在 Thermal Conditions 下选择 Temperature ,并设置温度值为 1000k; 3. 对右墙设置边界条件(wall-4:007. )
(1)将 Zone Name 改名为 right;
(2)在 Thermal Conditions 下选择 Temperature ,并设置温度值为 2000k; 4. 对顶墙设置边界条件(wall-4:005);
(1)将 Zone Name 改名为 top;
(2)对于绝热墙,保持默认设置的热边界条件(heat flux=0)。
第 5 步: Rosseland 模型求解
1. 为控制解设置参数
Solve Controls Solution...
(1) 在 Equations 和 Under-Relaxation Factors. 下,保持默认值。 (2) 在 Discretization 下 , Pressure 选择 PRESTO! , Momentum and Energy 选 Second Order Upwind。 2. 流场初始化
Solve Initialize Initialize...
(1)将温度设置为 1500k, 单击 Init . 3. 计算时,激活残差曲线显示。
Solve Monitors Residual...
(1) 在 Options 下,选 Plot . (2) 单击 ok 。
注意:Rosseland model 不能解决额外的辐射热传播方程,因此显示不出残差曲 线来。但可以解出能量方程中的导热系数,当使用 the P-1 and DO 模型时,就 可以解决辐射输运方程,能显示出辐射残差曲线。 4. 保存 case 文件。(rad_ross.cas) File Write Case...
5. 开始计算,进行 200 次反复收敛。 Solve Iterate...
大约在 180 次时计算就会收敛。
6. 保存数据文件(rad_ross.dat).
File Write Data...
第 6 步: Rosseland 模型后处理
1. 显示收敛曲线。
Solve Monitors Residual..plot
2. 显示速度矢量 。
Display Vectors...
“CAE/CFD创新工场
” blog.sina.com.cn/caecfd Email:cae_cfd@163.com
图 5.3: Rosseland 模型的速度矢量
3. 显示流量函数等高线(如图 5.4)。
Display Contours...
图 5.4: Rosseland 模型的流函数的等势图
由于箱内的自然对流的存在, 就观察图形中呈现循环模式。 在光学厚度 (0.2)很小 时,辐射不会对流动造成很大影响,得到的流场与没有辐射时相似(如图 5.5)。 但 the Rosseland 模型显示的流场很对称(如图 5.4),跟纯粹的自然对流不一样。 这个差异的产生主要是因为对较小的光学厚度而言 , the Rosseland 模型不是很 合适。
4. 显示温度等势图。 Display Contours...
“CAE/CFD创新工场
” blog.sina.com.cn/caecfd Email:cae_cfd@163.com
图 5.5: Rosseland 模型的温度的等势图
5. 计算每个侧墙的总的壁面热通量。
Report Fluxes...
(1) 在 Options 下选择 Total Heat Transfer Rate。 (2) 在 Boundaries 下选择 right 和 left 。 (3) 单击 Compute 按钮。
热墙和冷墙总的传热率大约为 7.43×105, 两侧墙上总的热通量大致平衡。
“CAE/CFD创新工场 ” blog.sina.com.cn/caecfd Email:cae_cfd@163.com 5. 沿着水平中心线显示 y 向速度图。
在 y=0.5 处,沿着箱体的中心处的水平线创建一个切截面。
Surface Iso-Surface...
1) . 在 Surface of Constant 下拉列表中选 Grid... ,在下面的列表中选 Y-Coordinate 。
2) . 单击 Compute ,看计算域的范围。
3). 在 Iso-Values 中设定值为 0.5,将 New Surface Name 改为 y=0.5.
4) . 单击 Create 来创建切截面。
(2) 在切截面上创建 y 向速度的 XY 散点图。
Plot XY Plot...
在速度剖面图中,右侧热墙呈现上升趋势,左侧冷墙呈现下降趋势。
将绘图数据存储文件。在面板中选择 Write to File 命令,单击 Write... 按钮。在 接下来的 Select File 对话框中,在 XY File 文本输入框中输入文件名,单击 OK 即可。 最后保存 case 和 data 文件(rad_ross.cas 和 rad_ross.dat)。
范文三：圆管自然对流计算和模拟
水平管和竖直管自然对流计算汇总 1. 计算工况表
2. 变化曲线图
圆管自然对流的计算和数值模拟
已知条件如图 1所示:将一圆管分别水平放置和垂直放置在大空间中进行自 然对流换热, 圆管外径 38D mm =,长度 1000L mm =,空气温度 20T C ∞= ,恒壁 温条件 100,150,200,250,300w T C = , 求解自然对流换热系数和换热量以及对流换 热时的空气最大速度。
图 1
一、数值计算
1. 自然对流换热系数和换热量的计算
1) 圆管水平放置计算
以壁温 100w T =℃ 为例,计算过程如下: 特征长度:0.038D m =;
定性温度 ()()21002060m w t t t C ∞=+=+= ;
查空气物性:()0.029W m K λ=?; -62=20.110m ν?; Pr 0.696= 空气的体积膨胀系数:())127360273333v m t K α-=+=+=
格拉晓夫数 Gr :
大空间自然对流的实验关联式为:
()Pr n
Nu C Gr =(1-1)
根据计算的格拉晓夫数 Gr 选择合适的常数 C 和 n (表 1) :
表 1 式(1-1)中的常数 C 和 n
()()3
35
2
62
9.81/333100200.038=3.21020.110v w g t t D Gr αν∞--??-?==??()
由式(1-1)和表 1可得:
故水平圆管换热量:
()()=hA7.9583.140.03811002075.962w t t W ∞Φ-=????-=
按照以上相同的步骤, 在给定恒壁温 100,150,200,250,300℃的情况下,可以 计算出相应的自然对流的换热系数和换热量,计算结果列于表 2中:
表 2水平管计算工况表
h , ()2W m K ?
7.958 9.115
10.045
10.803
11.527
φ, W
75.962
141.388 215.734 296.472 385.128
2) 圆管垂直放置计算
以壁温 100w T =℃ 为例,计算过程如下: 特征长度:1H m =
定性温度 ()()21002060m w t t t C ∞=+=+= ;
查空气物性:()0.029W m K λ=?; -62=20.110m ν?; Pr 0.696= 空气的体积膨胀系数:())127360273333v m t K α-=+=+=
格拉晓夫数 Gr :
()
()1/4
1/45
0.48Pr =0.483.2100.696=10.428
Nu Gr =????()2
3
0.02910.4287.9583810
Nu h W m K D λ-??===??()()3
3
9
2
62
9.81/333100201=5.831020.110v w g t t L Gr αν∞--??-?==??()
对于竖圆柱按照竖壁同用一个关联式必须满足:
经验算,并不满足情况,应该按照文献【杨世铭 . 细长竖圆柱外及竖圆管内 自然对流传热】中的关联式进行计算。
表 3竖圆柱自然对流关联式
由表 1可得: 先计算 ()()1/41/41/490.038=Pr 5.83100.6969.591
D D
Ra Gr H H ?
?=???= ()()1
1/41/4
1/4
-1
90.590.52 =0.59+0.529.595.83100.696 =162.594
D Nu Ra Ra H -??
??=+???? ???????
????????
故水平圆管换热量:
()()=hA4.7153.140.03811002045.008w t t W ∞Φ-=????-=
按照以上相同的步骤, 在给定恒壁温 100,150,200,250,300℃的情况下,可以 计算出相应的自然对流的换热系数和换热量,计算结果列于表 2中:
表 4水平管计算工况表
1/435H
d H Gr ≥ )20.029162.594
4.715H
1
Nu
h W m K λ??=
==?
h , ()2W m K ?
4.715 5.369 5.899 6.335 6.754
φ, W
45.008
83.390 126.703 173.860 225.649
2. 水平管(H )和竖直管(V )自然对流换热系数和换热量的对比图形
图 2换热系数
图 3换热量
3. 计算结果分析
由图 2和图 3可知:
1) 水平放置的圆管自然对流的换热系数和换热量都明显高于竖直放置的圆
管;
2) 随着温度的增加,两者换热系数和换热量都逐渐呈线性增长; 3) 水平圆管自然对流换热系数相对增加较多。
二、数值模拟
1. 水平圆管的数值模拟
1) 物理模型
如图 4所示, 本文采用的物理模型为大空间自然对流, 外边界设置为压力出 口边界,与大气相通,内边界为高温管道壁面,圆管直径按照实际尺寸设计。用 ICEM-CED 建立的模型长为 380mm ,宽为 380mm ,圆管直径 38mm
,位于中心位
置。
图 4
2) 网格划分
本次模拟的网格为结构化网格, ICEM 网格划分需要对物理模型进行分块处 理(block ) ,块的划分采用 O-block , O — block 易于对内边界做网格加密处理, 块的划分和网格的生成如图 5和图 6所示。
图
5
恒温壁面 压力出口边界
20℃空气
图 6
4) 网格质量和网格无关性验证
经网格无关性验证后,网格质量符合要求,网格划分合理。 5) 计算结果与分析
自然对流是由于空气温度差引起的密度差, 从而产生浮升力推动空气运动的 现象,实质属于可压缩流动。在 Fluent 中气体模型采用 Boussinesq 可以得到比 较好的模拟结果。 Boussinesq 近似是将动量方程中密度定义为时间的函数,而能 量方程中的密度视为常量。 在 Fluent 中设置好参数和边界条件后, 计算结果如下: ①壁温 100℃模拟结果
图 7 温度云图
图 8 压力云图
图 9速度云图
图 10旋涡
最大速度可在云图中直接读出:0.476m/s。 ②壁温 150℃模拟结果
图 11 温度云图
图 12 压力云图
图 13 速度云图
图 14 旋涡
最大速度可在云图中直接读出:0.537m/s ②壁温 200℃模拟结果
图 15温度云图
图 16压力云图
图 17速度云图
图 18旋涡
最大速度可在速度云图中直接读出:0.585m/s ②壁温 250℃模拟结果
图 19温度云图
图 20压力云图
图 21速度云图
图 22旋涡
最大速度可在云图中直接读出:0.697m/s ②壁温 300℃模拟结果
图 23温度云图 图 24压力云图
图 25速度云图
图 26旋涡
最大速度可在云图中直接读出:0.736m/s 速度随着温度变化的汇总表:
表 5 水平管最大速度计算工况表
max u , m/s
0.476
0.537
0.585
0.697
0.736
结论分析:
1) 自然对流换热强弱取决于高温壁面温度与周围流体温度差的大小,温差
越大,换热发展越迅速,流动越强烈;
2) 随着壁面温度的增加,最大空气流速也在随之增加;
3) 在温差的驱动下形成上升流,并在压差作用下上升流两侧形成漩涡。 2. 竖直圆管的数值模拟
1) 物理模型
如图 4所示, 本文采用的物理模型为大空间自然对流, 由于物理模型左右对 称, 故只需模拟其中的一侧即可, 同样外边界设置为压力出口边界, 与大气相通, 内边界为高温管道壁面,圆管直径按照实际尺寸设计。用 ICEM-CED 建立的模型 长为 2000mm ,宽为 570mm ,圆管直径 38mm ,位于中心位置。
图 27
2) 网格划分
本次模拟的网格为结构化网格, ICEM 网格划分需要对物理模型进行分块处 理(block ) ,内边界的网格加密处理,块的划分和网格的生成如图 28和图 29所 示。
图
28
恒温壁面
20℃空气
绝热
绝热
网格加密
图 29
6) 网格质量和网格无关性验证
经网格无关性验证后,网格质量符合要求,网格划分合理。
7) 计算结果与分析
自然对流是由于空气温度差引起的密度差, 从而产生浮升力推动空气运动的 现象,实质属于可压缩流动。在 Fluent 中气体模型采用 Boussinesq 可以得到比 较好的模拟结果。 Boussinesq 近似是将动量方程中密度定义为时间的函数,而能 量方程中的密度视为常量。 在 Fluent 中设置好参数和边界条件后, 计算结果如下:①壁温 100℃模拟结果
图 30温度云图
图 31压力云图
图 32速度云图
图 33旋涡
最大速度可在云图中直接读出:0.831m/s。 ②壁温 150℃模拟结果
图 34温度云图
图 35压力云图
图 36速度云图
图 37旋涡
最大速度可在云图中直接读出:1.05m/s ②壁温 200℃模拟结果
图 38温度云图
图 39压力云图
图 40速度云图
图 41旋涡
最大速度可在速度云图中直接读出:1.18m/s ②壁温 250℃模拟结果
图 42温度云图
图 43压力云图
图 44速度云图
图 45旋涡
最大速度可在云图中直接读出:1.29m/s ②壁温 300℃模拟结果
图 46温度云图
图 47压力云图
图 48速度云图
图 49旋涡
最大速度可在云图中直接读出:1.39m/s 速度随着温度变化的汇总表:
表 6竖直管最大速度计算工况表
max u , m/s
0.84
1.05
1.18
1.29
1.39
结论分析:
1) 自然对流换热强弱取决于高温壁面温度与周围流体温度差的大小,温差 越大,换热发展越迅速,流动越强烈;
2) 随着壁面温度的增加,最大空气流速也在随之增加;
3) 在温差的驱动下形成上升流,并在压差作用下上升流两侧形成漩涡。 3. 水平圆管和竖直圆管自然对流的最大速度对比
1) 现将模拟的最大速度汇总,如表格 7:
表 7最大速度对比表格
水平管
max
u , m/s0.476 0.537 0.585 0.697
0.736
垂直管
max
u , m/s0.84 1.05 1.18 1.29 1.39
2) 水平圆管和竖直圆管自然对流的最大速度曲线图,如图 50:
图 50
结论:①竖直管的自然对流最大速度明显高于水平管自然对流的最大的速度; ②随着温度的增加, 两者的最大速度都逐渐呈线性增加, 且增加的幅 度越来越小;
范文四：有限体积法与LBM分区耦合模拟方腔自然对流
第 45 卷 第 5 期 Vol . 45 No . 5 西 安交 通大学 学报
2011 年 5 月Ma y 2011 J O U RN AL O F XI′A N J IA O TO N G U N IV ER SI T Y
有限体积法与 L BM 分区耦合模拟方腔自然对流
栾辉宝 , 徐辉 , 陈黎 , 陶文铨
()西安交通大学能源与动力工程学院 , 710049 , 西安
摘要 : 在已有的密度分布函数重构算子的基础上 , 推导出了温度分布函数的重构算子 ,解决了格
() 子 Boltzma n n 方法 LB M与有限体积法耦合计算传热问题的关键难题. 选二维方腔自然对流对
3 6 耦合方法进行了考核. 在瑞利数 R a = 10,10范围内 ,耦合结果同商业软件 FL U EN T 结果符合 得很好 ,并且各物理量在耦合界面处连续且光滑过渡. 通过残差曲线可以看出 , 耦合模型在密网格 以及大瑞利数情况下 ,数值稳定性要好于单一 LB M .
关键词 : 格子 Boltzma n n 方法 ;有限体积法 ;多尺度 ;耦合
() 中图分类号 : TB13 文献标志码 : A 文章编号 :0253Ο987 X 201105Ο0078Ο06
Coupl ing of FVM an d L BM f or Natural Convect ion in a Square Cavity
L U A N H ui bao , XU H ui , C H EN L i , TA O We nqua n
( )School of Ener gy a nd Po wer Engi neeri ng , Xi′a n J iao to ng U ni ver sit y , Xi′an 710049 , Chi na Abstract : O n t he ba si s of t he e xi sti ng de n sit y di st ri butio n f unctio n reco n st r uctio n op e rato r , t he t e mp erat ure di st ri b utio n f unctio n reco n st r uctio n op erato r wa s de rived to calculat e heat t ra n sf er by
( ) ( ) co up li ng of t he lat tice Boltzma n n met ho d LB M a nd t he fi nit e vol ume met ho d FV M . The p re se nt co up li ng met ho d wa s vali dat e d by t he 2D nat ural co nvectio n f lo w i n a squa re cavit y wit h
3 6 ( ) va rio u s Ra ylei gh n umber s R af ro m 10to 10. The re sult s f ro m t he co up li ng met ho d a gree well wit h t ho se by co mmercial sof t wa re FL U EN T , a nd all t he p hysical qua ntitie s cro ss t he co up led i n2 t e rf ace smoo t hl y . Acco r di ng to re si dual hi sto r y c ur ve s it i s li kel y t hat t he n umerical st a bilit y of t he p re se nt met ho d a re bet t e r t ha n t ho se of t he p ure LB M at fi ne gri d numbe r s a nd hi gh R a. Key words : lat tice Boltzma nn met ho d ; fi nit e vol ume met ho d ; multi scale ; co up li ng
自然界与工程领域中的许多物理现象或过程常. 对于本身包含多种尺度的物理过 应的控制方程
常发生在涵盖几个数量级的空间及时间范围内 , 被 程 , 试图只用一种数值模拟方法来研究 , 即使是可 [ 1Ο2 ] 称之为多尺度物理学现象. 近年来 , 多尺度数值 能的话 , 其有效性、经济性和对过程本质揭示的深 模拟为研究多尺度现象提供了一条很好的途径. 常 度都可能受到限制. 如仅采用宏观方法 , 则会存在见的数值模拟方法根据其适用范围通常可以分为 3 [ 3Ο4 ] 一些不足 , 诸如无法预知细微部分的过程以及引入 ( ) 个层次 : 宏观层次, 如有限差分法 FDM、有限
( ) ( ) 体积法 FV M、有限元法 F EM等 ; 介观层次 , 如 一些经验性很强的参数 ; 如果仅采用介观或者微观
() 格子 Boltzma nn 方法 LB M、直接模拟蒙特卡洛法 模型 , 需要的计算机资源往往是难以承受的. 例如
() D SM C等 ; 微观层次 , 例如分子动力学模拟法 在纳微米流动中 , 分子动力学受限于计算机资源的
( ) ( ) MD S、量子力学 Q M S等. 各种数值模拟方法都 需求 , 不能用于整个流动区域的模拟 , 而 Na vie r2
受限于一定的努森数 K n 取值范围 , 并且有与之对 Sto ke s 方程不能描述纳微米流动中连续介质假设
不成立的区域.
() () 收稿日期 : 2010Ο09Ο07 . 作者简介 : 栾辉宝 1983 - , 男 ,博士生 ; 陶文铨 联系人,男 ,教授 ,博士生导师 ,中国科学院院士.
() 基金项目 : 国家自然科学基金重点资助项目 50636050.
网络出版地址 : ht tp :/ / www . cnki . net/ kcms/ detail/ 61 . 1069 . T. 20110117 . 1044 . 001 . ht ml 网络出版时间 : 2011Ο01Ο17
ht tp : ?www . jdxb. cn
第 5 期 栾辉宝 ,等 : 有限体积法与 LBM 分区耦合模拟方腔自然对流 79
鉴于单一方法的不足 , 多尺度耦合模型就应运 Pat a n ka r 和 Sp al di ng 提出 , 际上首次由著名学者
( 而生了. 通过建立起几种层次上方法的耦合模型 ,并被命名为求解压力耦合方程的半隐算法 SIM2
) PL E 算法. 近年来 , 作者所在研究团队对 SIM2有望将宏观方法的高效性与介观、微观方法的精确 [ 1Ο2 ] 性结合起来 , 相互补充 , 扬长避短. 在宏观方法 PL E 算 法 进 行 了 改 进 , 先 后 提 出 了 CL EA R 和
中 , FV M 推导出的离散方程可以保证具有守恒特ID EAL 算法. 数值实验证明 , 新的算法能在一定程 性 , 而且离散方程系数的物理意义明确 , 是目前流 度上克服 SIM PL E 算法的不足 , 并且大大提高了计
[ 3Ο4 ] 动与传热的数值计算中应用最广泛的一种方法. 算的收敛速度和算法的健壮性.
宏观部分的模拟 FV M 部分采用基于同位网格同时 ,LB M 是近 10 年来发展最迅速的介观模拟方
[ 10 ] 的 ID EAL 程序.法 , 具有物理图像清晰、边界条件易于实施、并行
[ 5Ο7 ] 11 3 密度分布函数重构算子效率高等突出优点. 发展 FV M 与 LB M 的宏观2
介观耦合模型是非常有应用前景的. 已有的关于 FV M 与 LB M 耦合的难点在于 : 如何将耦合界 FV M 与 LB M 耦合计算的模型均是基于等温流动 面上的宏观信息量准确地传递给 LB M 中需要的分 [ 8Ο9 ] 考虑到许多实际流动问题都包含有温度的. 的布函数. 以图 1 所示耦合模型为例 , 计算区域被划 变化或受到热效应的影响 , 本文构造了用于计算非 分为 2 块 , 左侧采用 FV M 求解 , 右侧采用 LB M 求
等温流动问题的 FV M 与 LB M 耦合模型.解 , 中间是重合区. M N 是 FV M 区域的右边界 ,
A B 是 LB M 区域的左边界. M N 边界上的值可以由 1 FV MΟL B M 耦合方法 LB M 区域对应点上的宏观值得到 , 然而 A B 边界
上的粒子分布函数就不能直接由 FV M 区域对应点 11 1 L BM 模型
上的宏观量直接得到 , 需要根据重构算子进行重() LB M 采用标准的二维九速 D2Q9模型. 求解[ 8Ο9 ] ( ) 构. 作者先前的工作 已经推导出了式 5所示的 [ 5Ο6 ] 温度场采用双分布函数 ( DD F) 格子模型, 其中
密度分布函数重构算子 , 并成功地应用到模拟二维 ( ) 密度分布函数 f x , t用于模拟速度场 , 而温度分 i
等温流动问题中. 在此为了节省篇幅 , 仅给出密度 ( ) 布函数 gx , t用于模拟温度场.i
分布函数重构算子的表达式 在模拟方腔自然对流问题时 , 为便于处理由温eq - 2 τΔ( f i = f i [ 1 - ftU i β cs U iα5 xuβ +[ 4 ] α 差引起的浮升力项 , 常常采用 Bo u ssi ne sq 假设. 2- 1 ρ ν ρ) νxuβ +Sαβ5 x ]( )5 5 α α 重力项中的密度可以表示为
式中 : U α = cα - uα ; Sαβ = 5 uα + 5 uβ.iixx β α ( )ρ ρβ( ) 1 = [ 1 - T - T]0 0
式中 : ρ和 T是冷面的密度和温度;β为体胀系数. 0 0
重力项可以表示为
ρρβ( )( )G = g - gT - T2 0 0 0
式中 : g 为重力加速度.
Bo u ssi ne sq 方程使用 LB M 求解可以通过在演
[ 5 ] ( ) 化方程右侧加入外力项来实现, 式 2改写为
( ΔΔ) f x + ct , t +t=i i
1 eq ( ) ( ) f x , t- [ f x , t- x , t ] + F3( ) ( ) i i f i i τ f
Fi =
图 1 FV M 与 LBM 分区耦合示意图 2 3 9 3 1 2 ( c?F) ( c?F) + - F i i ωρ1 -i 2 4 2 τ2cccf
( )4 1 . 4 温度分布函数重构算子
) : F = - gβ( T - 式中 T0 . 对于模拟传热问题 , 耦合边界上仅有密度分布
11 2 FVM 模型函数重构算子是不够的 , 还需要通过宏观物理量重
FV M 是在流动与传热问题求解中应用最广泛 构出温度分布函数 , 即要得到关于温度分布函数的 的一种数值方法. 在 FV M 中 , 压力修正算法是求 重构算子.
解不可压流场的主导方法. 该算法于 1972 年在国 使用 Chap ma n2En sko g 方法 , 对时间导数、空
ht tp : ?www . jdxb. cn
西 安 交 通 大 学 学 报 第 45 卷 80
g的二阶项可以表示为 间导数和温度分布函数进行多尺度展开i ( ) ( )12 2 ( )εε6 1 5= 5+5t t t ( )( )( )( )( )2 2 0 1 1 = - τΔt 5t + 1 - gi g D g = giii( )τ1 2g ( ) ε75 = 5 xx α α ( )( ) ( )1 0 12 2 ( ) ( )( )( )( )20 1 1 0 ( )εε 8 g= g+g +gτ - Δ i i i i - τΔ 5g t t g i - tD i [ Di g i ] = g 2 [ 8Ο9 ] 0 1 2 类似于 f 的推导过程, g在ε、ε和ε上的表 i i ( ( )( ( )) ) - 2 20 12 0 ( )- τΔt{ 5g - αc[ D ] g } 19 gt i T s i i达式为 ( ) 式 19中最后一个花括号内的第 2 项较第 1 项 ( )eq 0 0 ( )ε= gi g 9 : i 是个小量 , 可以忽略. 上式可以简化为( )( )1 0 = - Δt τ D g i + g i ( )( )( )( ) ( )2 20 02 i τΔ= - t5= - τΔg t[ 5 T i t T + gg g 5i g t i 2 (Δ) O[ t]( )10 ( ) ( )( )( ) 02 0 12 τΔ5g 5= - t{ a 5 g [ 5 uβ ] ] T + ui tgT T i x β α ( ) ( )( ) ( )11 20 = - Δtτ + 5t g i g [ D i g i ] - i ( )( ) - 2 - 1 0 12 - 1 (ρuβ) } = - τΔt ωρgc βc T5a T g[ 5 ] T +s i it T i x 2α 2 (Δt) τ g( ) ( )12 0 3 ( )+ O[ (Δt) ] [ D i ] gi 11 2 ( ) ( )( )- 2 - 1 11 1 ωρν(ρ ( 5 ) ) uβ + 5 x uα ] =c βc T[ 5 xi ix s α α β ( ) ( ) 式中 : D= 5+ cα 5 . 对式 10、11求零阶速度i t ix α ( )( ) 0 12 - 1 2 )- τΔt{ a T g( 5 T + gT i x α εε矩 , 可以得到 t=t 和 t =t 时间尺度上的宏观 1 2 方( ) ( )( )11 1 - 2 - 1 ) (ωρνρ5 x uβ + 5 x uα+ 5 c βc T[i is x α α β 程 ( )( )( )1 1 1 ( )( )1 1 ( ) ρ5uβ + 5 uα5 ]} =xxx β α α ( )( )T + 5 uα T = 0 12 5 t x α ( )( ) - 1 - 2 - 1 0 12 ( )( )2 2 1 ρ )- τΔt{ a T g( 5 T + ωcβcTν? gT i x i i s α ( )T = a [ 5 ] T5 13 T xt α ( ) ( )( ) 12 11( )ρ( ) ρ20 [5uβ + S β α 5 ]} x x αα1 2Δt . 式中 : a = c τ T - s g 2 最终 , 可以得到温度分布函数 g的重构算子 i
( ) ( )( )12 2 0 eq 根据链导法则 , 可得 = g- τΔt ? εε g gi g= g+g ++i i i i eqeq eq - 2 - 1 eq ( )5g = 5 g 5g + 5g 5uβ14 t i T t u t i i ωT [ U αg5 T + U αT cβc5 uβ -i i xi i i s x β α α eq eq eq - 1- 1 eq2 ( )5 g = 5 g 5 T + 5g 5 uβ15 xi T i xui x ( )α α β α ρωT cβ5 ρ] - τΔt{ a T g 5 T + i i x gT xi αα eq ( )- 2 - 1 2 1 根据 g的表达式 , 对其求导得 i )ωρνρ( ρuβ + Sβα cβc T[5 5 ]}= s i i x x αα eq - 2 - 2 - 1 eq 2 ω( ) ω5g = 5 [T 1 + ccγuγ]= cT c βui ui s i s i i τΔ () ggi α 5 x T + a T 5 i { 1 - T ]} + t T [ U β β x αα eq τω( ) ( )f T cβf -i 16 gi i i - 1 +τΔtω ρρ ( )cβT5 21g i i x eq β ( )τ- 2 - 1 eq 0 U βf f i i ω( ) gT = 5 g = 5 [T 1 + ccγuγ]i T i T i s i
( ) 式中 : f i 可以通过密度重构算子式 5求得. ( )17
从而可以构造出 g的一阶项 i 耦合计算实例2 ( )( )( )1 1 0 = - τΔtD gi gig i = 为了验证温度分布函数重构算子能否正确应用 ( )( )( )( )1 0 1 0 ) 于模拟传热流动问题 , 选取具有基准解的方腔自然 - τΔt ( 5 g i + cα5 g i = gti x α
对流来进行研究. 方腔自然对流是研究传热流动的 ( )( )( )( )0 1 0 1 5t uβ + - τΔt[ 5 g 5 T + 5 g gT i tui β
[ 11Ο12 ] ( ) ( )( ) ( )01 01 经典问题 , 已有很多学者对此进行了研究 . 图 1 ( ) 5g 5 uβ] = T + cα5 g 5 ui xi T i x β αα
( ) ( )( ) 01 左壁面为高温壁面 Th = 1 给出了计算模型 , K, τΔ- t[ -g ( ) 5 g 5 uα T -T i x α ( ) 右壁面为低温壁面 T= 0 K, 上下壁面绝热. 方 c - 1 0( )( )1 ρ) 5(ρδg i u 5 uαuβ + pαβ +x α β 腔被分为左右 2 个区域 , 左侧采用 FV M 进行求解 , ( )( )( )( )0 1 0 1 ( ) ci α 5 T g i 5 T + 5g 5 uβ]= xuix α β α 右侧采用 LB M 进行求解 , 两者在耦合界面上进行 ( )( )( )( )0 1 0 1 - τΔt[ - uα5 g 5 T - uα5g 5 uβ - gT ixuix αβ α信息传递. 区域尺寸为 H = 1 m , 网格数为 201 × ( )( )- 1 0 1 201 . 流体的普朗特数 P r =ν/ a 固定为 01 71 , 密度 ρ5 ( )5 x δg u pαβ + i α β 3 - 5 2 ( ) ( )( )( )ρν 01 0 1 = 11 0 k g/ m, 运动黏度= 11 3 ×10 m/ s , 重力cα5 g 5 T + cα5g 5 uβ ] = i T i xi uix α β α2 () 加速度 g = 91 8 m/ s 方向垂直向下. 在瑞利数 ( )( )( ) ( )0 1 01 - τΔt[ U α5 g 5 T + U α5 g 5 uβ - gi T ix i ui x β αα 33 6 βΔνRa = gH T/ a为 10, 10的范围内进行模拟 , 并( ) ( )( )- 1 01 - 1 eq 1 ρ= - τΔt [ U αT g5 T + 5g i 5 x p ] gi i x uα β β 将 FV MΟLB M 结果同商业软件 FL U EN T 的结果 ( )( )- 1 1 1 - 2 ρ ωρ( ) ωuβ - T c β5 ]18 U αT c βc 5 i ix i i is x α β 进行了对比.
ht tp : ?www . jdxb. cn
第 5 期 栾辉宝 ,等 : 有限体积法与 LBM 分区耦合模拟方腔自然对流 81
3 6 : H 为方腔几何高度. 式中 图 2,7 给出了 R a = 10和 R a = 10时的流线 、
表 1 给出了耦合模型的计算结果同文献[ 12 ]结 等温线、u 速度等值线、v 速度等值线 、压力等值线
果的对比 , 其中 N u代表高温竖壁上局部 N u 的 ma x 以及涡量等值线图. 可以看出 , 在流线、等温线以及 最大值 , N u代表高温竖壁上局部 N u 的最小值. min 速度线图中 , 耦合方法同 FL U EN T 结果非常一致 , 从表中可以看出 , 耦合模型计算结果同文献值符合
并且这些物理量在耦合界面上连续且光滑过渡. 在 得很好.
压力等值线以及涡量等值线中 , 耦合模型在界面两 本文还进一步比较了单一方法与耦合方法的收
侧物理量存在偏差. 其原因为 FVM 采用不可压缩模 敛特性 , 收敛判据定义如下
(型 , 标准格子模型具有弱可压性 , A B 边界 图 1 所
) 示上的密度无法由 FVM 直接传递给 LBM , 而需要 ( u ( i , j , t +Δt) ( ) u i , j , t| + - | ??根据 LBM 区域内部邻近边界的值进行外推得到. ( Δ)v i , j , t +t ( ) ) - v i , j , t| / | ??最后 , 统计了局部 N u 及其平均值 N u, 定义 m
( ) ( ) ) ( )( u i , j , t| +| i , j , t| ?e 24 v |如下 ??
3 5 T 图 8a 给出了 R a = 10 情况下单一热格子模型 H ( )22 N u = - ΔT 5 x w 随网格数目变化时的收敛曲线. 可以看出 , 随着网 H 格的加密 , 最终收敛的残差值增大. 其原因为过分 ( )N um 23 = N u d y0 ?
36( ) ( ) aR a = 10 bR a = 10
图 2 不同 R a 时的流线图
36( ) ( ) aR a = 10 bR a = 10
图 3 不同 R a 时的等温线图
36( ) ( ) aR a = 10 bR a = 10
图 4 不同 R a 时的 u 速度等值线图
ht tp : ?www . jdxb. cn
西 安 交 通 大 学 学 报 第 45 卷 82
36( ) ( ) aR a = 10 bR a = 10
不同 R a 时的 v 速度等值线图 图 5
36( ) ( ) aR a = 10 bR a = 10
不同 R a 时的压力等值线图 图 6
36( ) ( ) aR a = 10 bR a = 10
不同 R a 时的涡量等值线图 图 7
N u 及其最大值和最小值位置与文献值的对比 表 1
3456 R a = 10 R a = 10 R a = 10 R a = 10 参数 本文 本文 本文 本文 文献[ 12 ] 文献[ 12 ] 文献[ 12 ] 文献[ 12 ]
N um 11 110 1 114 1 230 1 245 1 507 1 510 1 807 1 806 1224488
1133771717N u1 561 1 581 1 539 1 539 1 738 1 637 1 710 1 442 max
( ) y/ H01 099 01 099 01 145 01 143 01 083 01 085 01 375 01 368 max
N u01 690 01 670 01 584 01 583 01 746 01 773 01 978 11 001 min ( ) y/ H01 995 01 994 01 995 01 994 01 997 01 999 01 997 01 999 min 细密的网格要大大增加计算机的运算次数 , 舍入误 LB M 的收敛曲线. 可以看出 , 下的耦合方法与单一
[ 4 ] 差会增加. 从图 8 b 可以看出 , 单一热格子模型 耦合方法收敛曲线的振荡较 LB M 方法明显减弱 , 随着 R a 的增加 , 收敛曲线的振荡越加剧烈 , 数值 并且收敛时的残差值会有所降低. 其原因是耦合方 不稳定性越加明显 , 收敛时的残差值会随之增大. 法引入了宏观 FV M 方法求解部分区域 , 这部分区 这种不稳定性是热格子模型本身固有的不足. 域可以避除 LB M 方法固有的弱可压缩性及数值不
6 图 9 比较了 R a = 10、网格数为 400 ×400 情况稳定性 , 进而改善了计算结果.
ht tp : ?www . jdxb. cn
第 5 期 栾辉宝 ,等 : 有限体积法与 LBM 分区耦合模拟方腔自然对流 83
参考文献 :
E Weinan , EN GQ U IS T B , L I Xia ntao , et al . Hetero2 [ 1 ]
geneo us multi scale met ho ds : a review [ J ] . Co mmun Co mp ut Phys , 2007 , 44 : 367Ο450 .
[ 2 ] 陶文铨. 传热与流动问题的多尺度数值模拟 : 方法与
应用[ M ] . 北京 :科学出版社 , 2009 .
陶文铨. 计算传热学的近代进展 [ M ] . 北京 : 科学出 [ 3 ] 版社 , 2005 .
陶文铨. 数值传热学[ M ] . 2 版. 西安 : 西安交通大学 [ 4 ] 出版社 , 2000 .
郭照立 , 郑楚光 , 李青 , 等. 流体动力学的格子 Boltz2 3( ) aR a = 10 ma nn 方法[ M ] . 武汉 : 湖北科学技术出版社 , 2002 . 何[ 5 ]
雅玲 , 王勇 , 李庆. 格子 Boltzmann 方法的理论及
应用[ M ] . 北京 : 科学出版社 ,2008 . [ 6 ] H E Yaling , L I Qing , WA N G Yo ng , et al . L at tice Bo2
ltzma nn met ho d a nd it s applicatio n in engineering t her2 [ 7 ] mop hysics [ J ] . Chi ne se Science Bulletin , 2009 , 54
() 22: 4177Ο4134 .
徐辉 , 栾辉宝 , 陶文铨. LBM 与宏观数值方法界面信
息耦合的重构算子 [J ] . 西安交通大学学报 , 2009 , 43 [ 8 ] () 11: 6Ο10 .
XU H ui , L U A N H uibao , TAO Wenqua n. A reco n2
st r uctio n operato r fo r interf ace co upli ng bet ween LBM ( ) b网格数为 400 ×400 a nd macro2numerical met ho d of finite2f a mily[J ] . J o ur2 图 8 热格子模型随着 R a 以及网格数变化的收敛曲线
() nal of Xi′a n J iao to ng U niver sit y , 2009 , 43 11: 6Ο10 .
L U A N H B , XU H , C H EN L , et al . N umerical ill us2
t ratio ns of t he co upling bet ween lat tice Boltzmann [ 9 ] met ho d and finite2t ype macro2numerical met ho ds [ J ] .
N umer Heat Transf er : B , 2010 , 57 : 147Ο170 .
SU N D L , Q U Z G , H E Y L , et al . A n efficient seg2
regated algo rit hm fo r inco mp re ssible f l uid f lo w and [ 10 ] ( heat t ra nsf er p ro blems2ID EAL Inner Do ubly Iterative
) Efficient Algo rit hm fo r L inked Equatio ns : p art I
Mat hematical fo r mulatio n a nd sol utio n p rocedure [ J ] .
N umer Heat Transf er : B , 2008 , 53 : 1Ο17 . 6 图 9 R a = 10、网格数为 400 ×400 时耦合方法与
LBM 方法收敛曲线对比
[ 11 ] 童长青 , 何雅玲 , 王勇 , 等. 封闭空腔自然对流的格
子2Boltzmann 方法动态模拟 [J ] . 西安交通大学学报 ,
() 2007 , 41 1: 32Ο36 . 结论3
TON G Cha ngqing , H E Yaling , WA N G Yo ng , et al . 在已有的密度分布函数重构算子基础上 , 推导 Simulatio n of t ransient nat ural co nvectio n in squa re
出了温度分布函数的重构算子 , 解决了 FV M 与双 cavit y wit h inco mp re ssible t her mal lat tice2Boltzmann
分布函数 LB M 耦合计算非等温流动的难题. 以二 met ho d [ J ] . J o ur nal of Xi′a n J iao to ng U niver sit y ,
维方腔自然对流为例 , 对非等温耦合模型进行了验 () 2007 , 41 1: 32Ο36 .
证. 耦合模型的结果同文献以及 FL U EN T 的结果 [ 12 ] BA RA KO S G , M I TSO UL IS E. Nat ural co nvectio n
f lo w in a squa re cavit y revi sited : lamina r a nd t ur bulent 在定性定量上都符合得很好 , 证明了密度以及温度
mo del s wit h wall f unctio n [ J ] . Int J N umer Met ho d 重构算子的正确性以及耦合方法的可行性. 通过比
Fl uids , 1994 , 18 : 695Ο719 . 较大 R a 以及密网格情况下的收敛曲线 , 发现耦合
(编辑 荆树蓉) 方法可以改善热格子模型的数值不稳定性.
ht tp : ?www . jdxb. cn
范文五：封闭腔内水自然对流换热数值模拟
封闭腔内水自然对流换热数值模拟 第58卷第11期
2007年11月
化工
JournalofChemicalIndustryandEngineering(China) Vo1.58NO.11
November2007
封闭腔内水自然对流换热数值模拟
苏燕兵,陆军,白博峰
(西安交通大学动力多相流国家重点实验室,陕西西安710049)
摘要:为了揭示封闭腔内非Boussinesq流体在浮力驱动下所特有的流动换热现象和形成机理,采用CFD软件
Fluent对封闭腔内水的自然对流进行数值模拟,得到矩形封闭腔高宽比,Rayleigh数,倾斜角度,壁面温度差对
流动和传热的影响规律.研究结果表明:由于水的密度在3.98~C达到最大,两竖壁面温度跨越这一点时会引起
流动图像反转;具有流动反转的双涡结构降低了对流换热平均Nusselt数;相同Rayleigh数下,高宽比为1对应
对流换热平均Nusselt数最大值;倾斜角度对平均Nusselt数影响与Rayleigh数和温度边界条件有关.
关键词:自然对流;流动反转;高宽比;倾斜角
中图分类号:TK124文献标识码:A文章编号:0438一u57(2007)11—2715—06 Numericalsimulationofnaturalconvectionand heattransferofwaterincavities
SUYanbing,LUJun,BAIBofeng
(StateKeyLaboratoryofMultiphaseFlowinPowerEngineering,Xi'an
JiaotongUniversity,Xi'an710049,Shaanxi,China)
Abstract:Torevealthefeaturesofflowstructureandheattransferandthemechanismofthenon
Boussinesqliquidflowdrivenbythermalbuoyancyincavities,thenaturalconvectionofwaterinsquare
andrectangularenclosureswasnumericallysimulatedwithCFDsoftwareofFluent.Theeffectsofaspect
ratioofthecavity,Rayleighnumber,inclinationangleandtemperaturedifferencebetweenthetWOwallsof
thecavityontheflowandheattransferwereinvestigated.Theresultsshowthattheflowpatterninverses
ifthetWOwallstemperatureofthecavitywasgreaterandlessthan3.98~C,respectively,atwhichthe
waterdensityismaximum.Theflowpatterninversionhasadoublevortexstructureanddecreasesthe
averageNusseltnumberofnaturalconvectiveheattransfer.AtafixedRayleighnumber,theaverage
Nusseltnumberreachesamaximuminthesquarecavityattheaspectratioof1.Theinclinationofthe
squarecavityhasmorecomplexinfluenceontheaverageNusseltnumberwhichdependsnotonlyonthe
Rayleighnumberbutalsoonthethermalboundaryconditions.
Keywords:naturalconvection;flowpatterninversion;aspectratio;inclinationangle 自然对流换热因其在工程中的重要应用,如电
2006—12—04收到初稿,2007—08—07收到修改稿.
联系人:白博峰.第一作者:苏燕兵(198O一)男,硕士研
究生.
基金项目:国家自然科学基金项目(50536040);新世纪优秀
人才支持计划项目(NCET一04—0923).
子元件冷却,换热器,凝固融化,晶体生长等而被
广泛研究.通常,在数值模拟浮力驱动流动时要 引入Boussinesq假设来简化模型,然而对于液态 Receiveddate:2006—12—04.
Correspondingauthor:Prof.BAIBofeng.E—mail:bfbai@
mail.Xjtu.edu.en
Foundationitem:supportedbytheNationalNaturalSeienee
FoundationofChina(50536040)andProgramforNewCentury
ExeellentTalentsinUniversity(NCET—O4O923). ,111n
?
27j6?化工第58卷
的水,镓,铋,液氦,HgCdTe等在特定温度
时密度有最大值的流体,如果采用基于线性关系的 Boussinesq假设,在计算包含此温度的自然对流 时,则不能获得真实的流动图像].
矩形腔自然对流换热方面已经有大量文献报 道,但大多数是针对Boussinesq流体].非
Boussinesq流体一般以水为研究对象.Hideo等] 做了方腔水对流倾斜实验并进行简单数值模拟. Lin等[6建立了比较完整的数学模型研究方腔内水 在跨越密度最大时的自然对流.Braga等]对矩 形腔内的水自然对流进行了模拟并且用实验验证了 温度分布和流动图像.Wei等[9在数值模拟方面 做了比较多的工作,重点考察了高宽比对流动的影 响.Braga等[8和李恒等[1分别用实验研究了矩形 腔和圆柱腔内水凝固过程的自然对流.上述这些工 作已经重点揭示了水的流动反转流场,温度场变 化,但是对换热的分析较少,对物理机理解释也不
够深入..
本文用Fluent商业软件,数值模拟研究矩形 腔内在水密度最大值附近的自然对流换热,重点分 析流动反转及其物理机理.模拟得到热壁面温度, 高宽比,倾斜角度等对流动和对流换热的影响 规律.
1物理模型和数值方法
1.1控制方程和边界条件
,左壁面 选取的物理模型和坐标系如图1所示
为高温Th,右壁面为低温T,其他壁面绝热.对 较小的Rayleigh数(Ra一10.,10)分别进行二 维和三维的稳态层流模拟,计算结果表明,在此 Rayleigh数范围,三维流动效果微弱.取三维结果 中如图1所示的中截面上的温度,速度分布与二维 结果进行比较,非常接近,因此在本文中为了节约 计算时间均采用二维模型进行模拟.
水的密度在常压下3.98~C达到最大,选用抛 adiabaticOT/OY=0
I
lkL
//////////,…,,,,,,0/
aT/ay=0
图1矩形腔物理模型
Fig.1Schematicdiagramofphysicalmodel
物线形式表达式
lD一[1一y(T—To)](1)
式中T0是对应水密度最大值时水的温度,热 膨胀系数y一8.0×10K一.
定义Rayleigh数和Prandtl数
Ra—gY(Th—)L./(Aav) Pr—v7a
用如下的量纲1参量
X一百217,Y一青,U一U, V一,P一,一器
对稳态的Navier-Stokes方程及能量方程进行量纲
1化可得控制方程组
oU
十一
aV
XoY—O(2)a.…
u孺au十au—Pr(+)一aP.RaPr2COS~(3)
uav
十y
av
—
Pr(+)一+RaPrO2sin(4) ua0
十y
a0
一+(5)
物理模型对应边界条件为
o(o,y)一Oh,(1,y)一Oo
OY!一aY!一ou
U—V一0,在4个边界壁面上 1.2局部和平均Nusseit数
Nusselt数是对流换热研究中的一个重要参 量m].对于本问题,冷,热壁面上的局部Nusselt
数计算方法为
Nu—
I一?
平均Nusselt数用式(7)计算
Nuav一q(7)
1.3数值方法
用Fluent6.0软件对上述模型进行模拟求解. 计算中采用四边形均分网格,Fluent中设定稳态 隐式有限容积SIMPLE算法,计算精度选双精度 格式.压力项用普通格式,对流项用QUICK格 式,能量项用二阶迎风格式,并选用合适的亚松弛 因子.因为模型通过改变几何尺寸来改变Ra大 小,所以对不同模型都要进行网格无关计算验证. 2模拟结果与分析
2.1密度变化引起的流动反转
冷壁面温度恒设定为273K,图2给出了不同
第11期苏燕兵等:封闭腔内水自然对流换热数值模拟?27J7?
热壁面温度时的温度场和流场计算结果.从图上可 以看出:在热壁面温度为277K时,热壁面上部和 冷壁面下部的温度梯度较大,流场只有一个逆时针 方向涡结构.在热壁面温度为281K时,水密度最 大对应的等温线在垂直中线附近,其两侧的密度变 化趋势不一样,形成一对流动方向相反的涡,密度 反转面就像是互不相溶的液液界面.在热壁面温度 为285K时,密度反转面开始向右下角倾斜,逆时 针涡被顺时针大涡挤压到右下角区域.上述物理模 型在Bossinesq假设下,流场只会有一个涡;对出 现这种流动反转的双涡结构,此时Bossinesq假设 不再适用.
(b)Th一281K
(c)』h一285K
图2Ra一10热壁面不同温度时的温度场和流场 Fig.2Temperaturefieldandstreamlineatdifferent
hotwalltemperaturewithRa一10
图3是Ra===10时水平中截面上量纲1速度 ,温度在量纲1横坐标x方向的分布曲线. X曲线表明:在T一277K时冷热壁面附近水的 流动方向与Th一281K,285K时相反,并且Th一 277K,281K时两壁面附近绝对值相等,Th一 285K时的绝对值在热壁面附近大于冷壁面附近 的.Th===277K与Th一285K时有两个极值点; 而T一281K速度有3个极值点,且在竖直中截 面达到绝对值最大,曲线对于X一0.5轴对称. 因为y一0.5不穿过T一285K时右下角逆时针 涡,所以X曲线右侧只能看到小涡存在,导致 (b)temperature
图3y一0.5截面上量纲1速度V和温度0分布 Fig.3Distributionsofnon—dimensionalvelocity
Vandtemperature0atsectiony一0.5
变化率减慢.
X曲线可以用一X曲线解释:水在3.98? 时密度最大,在重力作用下必然向下流动,大的 x曲线变化率对应大的速度值.从x曲线可以 看出:T一277K时,热壁面附近水密度大向下流 动;T===281K时,水密度最大值在X一0.5附 近;T===285K时,在冷壁面下半部分附近出现密 度反转.T一277K时的X曲线在冷热壁面附近 变化率大小相等;T一285K时的x曲线在热壁
面变化率大于冷壁面;T一281K时x曲线有3 个梯度区,并且x一0.5附近变化率大于冷,热壁 附近变化率.
图4给出的是冷热壁面上局部Nusselt数对量 纲1纵坐标的Nu—Y曲线.Th一277K时Nu最大 值在冷壁面y一0.02,热壁面y一0.97处,即对流 最强烈的腔体左上角和右下角;T一281K时的 冷,热壁面Nu最大值都在Y一0.04,即腔体左,
K时热壁面Nu最大值在 右两个下角;T一285
?27】8?化工第58卷
y一0.03,冷壁面附近由于小涡存在Nu不再是单 调变化,最大值在Y一0.96,最小值在Y一0.4 (两个涡的分离点).在方腔角点区域,因为壁面黏 滞作用使得水流动极其缓慢,所以N"减小很快. I
l
li
?Th=285K矗A
?,oooooooo....II.--.."AA^?^^^ooo00.20.40.60.81.0
y
(b)coldwall
图4方腔热壁面和冷壁面局部Nusselt数分布 Fig.4DistributionoflocalNusseltnumberon
hotwaUandcoldwallforRa一10
2.2高宽比影响
引入高宽比A=H/L,保持L不变,通过改 变H来改变A.图5是热壁面温度为12?,Ra为 10.情况下不同高宽比(A一0.25,0.5,1,2,5,
10)时模拟计算的流线图.从图中可以看出:当
A?1时,随着A的增大,腔体右下部分逆时针方
向的涡的影响区域也逐渐增大,顺时针方向的大涡
被此逆时针涡挤得越来越贴近热壁面.当A<1 时,逆时针方向涡经历先减弱后增强的过程;A一
0.5时逆时针涡只占右下角很小部分,流动也微
弱;A—O.25时逆时针涡则占据整个腔体高度,涡
流强度也增大.
图6给出不同热壁面温度时的平均Nusselt数
(Nu.)一高宽比(A)曲线(A一0.25,0.5,1, 2,5,10),从图中可以看到:Th相同时,Ra越
大Nu越大,即对流越强换热越快.Ra相同时
图5Th一285K,Ra=10时不同高宽比对应的流场
Fig.5Streamlineatdifferentaspectratiowith
Th一285K.Ra一10
(a)丁:277K
(c)Th一285K
图6不同热壁面温度条件下时高宽比
对平均Nusselt数的影响
Fig.6EffectofaspectratioonaverageNusselt
numberwithdifferenthotwalltemperature
加84
第11期苏燕兵等:封闭腔内水自然对流换热数值模拟?27j9? 謦令=
30.~--6o.=90.Ca=120.妒150.妒180.
图7方腔在Ra一10,11h一285K时不同倾角下流场分布 Fig.7StreamlineatdifferentleananglewithRa一10,11h一285Kinsquarecavity
A一1附近N‰最大;A>1时Nu随着A增大而
减小,并且减小速度逐渐缓慢;A<1时N随着 A的减小而减小,并且AdO.5时迅速减小.相同 Ra,相同高宽比,一281K对应的Nu最小. 分析流场图和Nu一A曲线,可以得到结论: 双涡结构越明显,换热性能越差.这是因为涡的交 界面上流体流动方向一致,互不掺混,相当于多了 一
层热阻隔断了冷热壁面之间的直接联系,这样大 大减小了热传递速率,而且两个涡强度相差越小换 热效果越差.从结果上看,高宽比对换热影响的规 律比较一致.
2.3倾斜角度影响
如果方腔倾斜,重力在两直角坐标轴上都有分 量,必然会对流动产生很大的影响,图7给出了在 T一285K,Ra一10.时,倾斜角度?分别为30., 60.,90.,120.,150.(热壁面在下为0.,热壁面 在上为180.)时的流场.从图中可以看出:?< 60.,以单涡结构为主,随倾斜角度的增大过渡到 越来越明显的双涡结构.
方腔倾斜影响流场,也必然会影响到换热效 果,图8给出不同热壁面温度下的N"一?曲线, 可以看出Ra仍然是影响换热的最主要因素.Ra一 10.时Nu基本不随倾角变化;Rn>10.时,丁h一 277K,N随倾斜角度增加先增后减,?一120.附 近最大,?:0.时最小,即热壁面在上比热壁面在 下换热效果好;T一281K,Nu对于?一90.基本 呈轴对称分布,且在?一90.达到最小;T一285K 时,Nu先增后减,约在?一6O.附近达最大值, ?一180.时最小.在相同Rayleigh数下,T一281
K时N"受?的影响比Th一277K,285K小. 形成上述情况的原因主要有两点:一是形成了 流动微弱的"死水区",换热以热传导为主,热传 导换热效率远小于热对流;另一方面是形成强度近 似的双涡,阻断了冷,热壁面的直接联系,增加了 热阻.因此,为了达到良好的换热效果,选择合适 的倾斜角度非常必要.
妒/(.)
(a)丁h=277K
妒/(.)
(b)丁h一281K
(c)丁h一285K
图8方腔不同热壁面温度下平均Nusselt 数随倾斜角度的变化
Fig.8AverageNusseltnumberunderdifferentlean
angleofsquarecavitywithdifferenthotwalltemperature
?
2720?化工第58卷
3结论
本文用数值模拟方法研究了封闭腔内水的自然 对流,得到非Bossinesq流体流动,换热特性,并 且得到了不同物理,几何条件对流动,换热影响. 主要结论如下:(1)当封闭腔体内存在水温比较缓 慢跨越T0区域时,则描述水自然对流的 Boussinesq假设不再适用;(2)流动反转的双涡结 构阻断了冷,热壁面的直接联系,减弱了对流换热 效果;(3)Rayleigh数对流动,传热具有决定性的 影响,高宽比A,倾斜角度?对流动和换热有重要
影响;(4)在一定Ra和温差?丁下,一定的高宽 比或者倾斜角度对应着最佳的传热特性. 符号说明
A——矩形腔高宽比
H,L——分别为腔体高度,宽度,121 ^——对流传热系数,W?121?K q——热流量,W?121
?T_一冷热壁面温差,K
U,——量纲1速度
x,y——量纲1坐标
口——热扩散系数,m2?s 倾斜角度,(.)
A——水的热导率,W?121?K 卜热膨胀率,?_1
运动黏度,Pa?S
量纲1温度
下角标
av——平均值
c——冷壁面
h——热壁面
References
[1]J~rgenZierep,HerbetOertalJr.ConvectiveTransportand
[2]
[3]
[4]
[5]
[6]
[7]
[8]
[9]
[10]
[11]
[12]
InstabilityPhenomena.Karlsruhe:Braun,1982 GrayDD,GiorginA,ThevalidityoftheBoussinesq approximationforliquidsandgases.Int.J,HeatMass Transfer,1976,19:545—55l
GdeVahlDavis.Naturalconvectionofairinsquarecavity: abenchmarknumericalsolution.Int.J.Numer.Meth. Fluids,1983,3:249—264
DengBin(邓斌),TaoWenquan(陶文铨).Three-
dimensionalnumericalsimulationofturbulentflowandheat transfercharacteristicssideofshell—and-tubeheat
exchangers.JournalofChemicalIndustryand Engineering(China)(化工),2004,55(7):
1053—1059
HideoInaba.TakeyukiFukuda.Naturalconvectioninan inclinedsquarecavityinregionsofdensityinversionof water.J.FluidMech.,1984,142:363—381
LinDS.NansteelMW,Naturalconvectionheattransfer inasquareenclosurecontainingwaternearitsdensity maximum.Int.J.HeatMassTransfer,1987,3O: 2319—2329
BragaSL,ViskantaR.Transientnaturalconvectionof waternearitsdensityextremuminarectangularcavity. Int.J,HeatMassTransfer,1992,35:86I-875 BragaSL,ViskantaR.Effectofdensityextremumonthe solidificationofwateronaverticalwallofarectangular cavity.Experiment.Therm.FluidSci.,1992,5:
703,713
WeiTong,JeanNKoster.Naturalconvectionofwaterina rectangularcavityincludingdensityinversion.Int.J.Heat FluidFlow,1993,14:366—375
WeiTong.Aspectratioeffectonnaturalconvectioninwater nearitsdensitymaximumtemperature.Int.J.HeatFluid Flow,1999,20:624—633
LiHeng(李恒),BaiBofeng(白博峰),LuJun(陆军),
GuoLiejin(郭烈锦).Experimentalstudyofthermal
convectionduringthewatersolidificationprocessincylinder cavity.JournalofEngineeringThermophysics(32程热物
理),2006,27(6):977—980
TaoWenquan(陶文铨).NumericalHeatTransfer(传热
学).2nded.Xi'an:Xi'anJiaotongUniversity
Press,2001
转载请注明出处范文大全网 » fluent自然对流关键设置
