基于方法 B的有效离散格式
1 引言
在对流体流动、传热及传质 问题 进行数值模拟时,首先要进行区域离散化。视节点
在子区域中位置的不同,可把区域离散化方法分为两类:方法 A(外节点法)及方法 B(
内节点法)。方法 A是指先确定节点的坐标再计算相应的界面,使节点位于子区域角顶的
方法;方法 B是指规定界面位置而后确定节点,使节点位于子区域中心的方法。
当网格不均匀时,方法 B中节点永远位于控制容积的中心,这样方法 B中的节点更能
代表控制容积,界面导数的计算精确度也要高。所求解的区域中物性发生阶跃性的变化时
,采用方法 B时可以较方便地反发生阶跃变化的面作为界面,从而可以避免在同一控制容
积内物性发生突变的情形,若采用方法 A时实现这一点要困难的多。综述之,一般通风空
调工程问题中采用方法 B进行区域离散化。本文的出发点也基于此,希望通过常用的方法
B离散网格建立有效的、强健的差分格式。
2 对流-扩散方程式的离散及基于方法 B的有效差分格式
对流换热问题的数学描述包括质量、动量及能量守恒三种方程式,可用统一张量形式
可表示为:
(1)
上述方程的四项分别是不稳态项,对流项,扩散项以及源项。据[1],无论守恒方程
式的何种离散或差分格式,必须满足四项基本法则,作者提出的守恒型离散格式亦必须满
足这四项原则。现给出一个通用的离散方程式:
(2)
用数值方法求解上述方程(尤其是满足动量守恒的 N-S方程)的主要困难来自非线性
的一阶导数项及压力梯度项的不妥离散方式。到 目前 ,为解决一些相关问题已提出了许
多的差分离散方案:指数、混合、幂函数(乘方格式)、QUICK格式、SHYBRID二阶混合
等方案。这些方案,在具体 应用 中,各有优劣。作者在暖通空调工程数值实验中,基于
B方法生成的交错网格推导出一种有效差分格式,特阐述如下:
控制方程组(1)中包含有两类变量:矢量与标量。{U,V,W}为矢量,扩散系数、温
度、压力、紊流脉动能 k、紊流耗散率ε等为标量。这样交错网格的数据存储可表述为:
标量数据皆存储在节点位置{P、E、W、N、S等};而矢量数据皆存储在截面位置{e、w、n
、s等};图 1中阴影区为标量控制容积区,标量分布均匀,由节点 P等代表。矢量亦均匀
分布于各截面,用 e等代表。下面皆以二维稳态情形为例进行 分析 。
2.1 标量场
首先分析标量场,以φ代表通用标量,由(2)其离散方程的一般形式可写为:
(3)
据幂函数方案 AI(I代表 e,w,n,s等界面)的取法,上式各系数可表述为:
(4)
(4)式中,Peclet数 PI=FI/DI,F代表单位截面上对流流量,而 D代表单元内物理
量的扩散。因为物性参数(标量)亦储存在各个节点上,据前述的四项基本法则,须满足
界面上连续性的要求,故可用调和平均来计算界面参数。以为例可由下式求得:
(5)
其他的系数同理可求。
2.2 矢量场
二维矢量场(即 U、V两个速度矢量的分布),而速度的控制容积与标量的不同,其
示意图如图 2所示。比较图 1与图 2可以看出:矢量控制容积只是在标量控制容积的基础
上沿特定方向发生错移,这也是交错网格的特点之一。以二维稳态离散 U方程为例,由 N-
S方程对上图的控制容积积分可以导出[3]:
图 1 标量控制容积与离散格式图案
图 2 矢量控制容积与离散格式图案
[注 1]:图 2中加"()"符号的节点一致;
[注 1]:图 1及图 2中的 分别等于将要用到的 分别等于将要用到的 。
同(3)所述,亦采用幂函数格式,上式中的系数可以按(4)表达之,仅需调整 I(代
表 E、W、S、N等节点位置)。从图 2可见,Peclet数 P=F/D的求解难在界面交点 N、S等
。界面效点 N处流速 VN须应用加权平均计算,这样,界面流量 FN可表述为:
(6)
扩散系数
(7)
式(7)中系数是标量,而处为界面位置,不存储标量,则可通过周围四个节点上值
的面积加权求解(若三维情形,应为体积加权),这样较容易满足上述的四项基本法则。
所以:
(8)
联立求解,则 DN的计算可进一步简化。同理,可以导出 DS的表达式:
(9)
(10)
至此,已得到整个 U方程系数的求解方法,对于其他矢量方程的系数,同理可求。
2.3 新格式与各种差分格式的特性比较及分析
对于一维稳态无源项的对流扩散问题,其守恒方程为:
(11)
假定边界条件为:
通过计算中间点 x=处Φ值,比较了 理论 解、一阶上风格式、中心差分格式及
本文新的格式:中心差分格式随着 Pe数的增加而产生发散现象;一阶上风格式和本文提
出的差分格式则与理论解有较为接近的值。特别当 Pe数增大时,本文的差分格式的更接
近理论解。
3 算例 应用 、比较及 分析
3.1 单出入口二维通风房间
图 3 具有单个出入口的通风房间
图 4 外掠后台阶的平面图
室内流体流态为稳态紊流,采用 LRN(Low-Reynolds-Number)模型 计算 。压力、速
度通过 SIMPLEC算法耦合求解。入口边界参数给定,出口边界变量的法向梯度为零。以入
口速度 U为基准速度,以入口尺寸 L(1m)为基准尺寸。基于 B 方法 的网格生成为 X×Y
方向网格数为 120×100。之后,作者重新布置网格(80×60),取某垂直截面段(X/L=15)
并比较了两者的水平速度分布。如图 5可见,两者几乎无差别,这可以说明本文格式应用
于工程实际时,对非均匀网格具有较好的适应性。
图 5 两套网格系统截面水平速度的比较
对两套网格系统计算过程的收敛史与差分格式未改进的情形进行比较(图 6),发现
计算精度有较明显的提高,可能的发散趋势被抑制。
图 6 迭代过程收敛史的比较(纵坐标即压力方程的余量)
[注]1-未改进的格式;2-格式改进后的第一套网格;3-格式改进后的第二套网格
3.2 带移动顶盖的方形空腔流
流体流态为稳态层流。以顶盖移动速度 U为基准速度,以空腔进深 D为基准尺寸。基
于 B方法的风格生成为 X×Y方向网格数为 20×20。固体壁面为无滑移条件,因为是层流运
动,故动量方程及连续性方程中的流速边界都是如此。
图 7 截面水平速度的比较(Y/D=)
图 8 截面垂直速度的比较(Y/D=)
取截面(Y/D=)上的速度矢量与 文献 [9]比较,由图 7及图 8可发现两者相差甚
微,再次证明本文改进格式的有效性。笔者在计算过程中还比较了一般程度消耗 CPU时间
一本文所花的时间相差不大,可能是因为物理模型规模较小所致。但是对于规模较大、风
格划分不均匀的情形,本文的离散方程具有不可替代的优势。
3.3 外掠后台阶的流动
如图 4所示,H为后台阶高度,并作为基准尺寸:入口及出口尺寸如图 4所示;未流
速度 U取为 1m/s,并作为基准速度。来流温度 T分布均匀,作为基准温度;壁面为均匀热
流条件。基于 B方法的网格生成为 X×Y方向网格数为 60×50。
流体流态为稳态二维紊流,服从 Bousinesq假设。对此开口流场,出口处采用自由流
速边界条件,压力修正方程界为 Dirichlet条件。壁面附近采用壁面函数法加以修正。中
心轴线上各个变量采用第二类边界条件。
图 9是截面(X/H=10)上的流速分布,与文献[4]结果相吻合。另外,作者还与文献
[5]比较了不同截面上 Nu/NuH数,如图 10所示,两者较接近。
图 9 截面水平流速分布的比较(X/H=15 , Re=×103)
图 10 外掠后台阶局部 Nu数分布的比较(Re=×103)
4 结论及展望
本文提出的差分格式从本质上将讲是对幂函数格式的修正及改进,其中应用了中心差
分、调和平均等思想,从而达到提高计算精度、减少发散的目的。但从一些运算结果来看
,格式并非十分理想,譬如对向后台阶流动(文中未列出),涉及对流换热、壁面边界等
问题 ,结果与文献[2]提供的数据有一定差距;尽管上述运算较原来的算法迭代次数有所
减少,但计算时间并未缩短,并且计算机内存容量要求变有提高,这些问题有待今后改进
。
参考 文献
1 Patankar. . Numerical heat transfer and fluid flow. Hemisphere,
Newyork, 1980.
2 Gooray. A. M. Turbulent Heat Transfer Computations for Rearward-facing
Steps and Sudden pipe Expansions. [J]. Heat transfer, 1985,1:70~76.
3 Larnder. B. E. The numerical computation of turbulent flows. Comp.
Meth. Appl. Mech. Engng 3, 269 (1974)
4 Gartling. D. K. A test problem for outflow boundary conditions -flow
over a backward-facing step. Int Numer Meth Fluids, 1990, 11: 953--967
5 Gooray. A. M and watkins. C. B. Turbulent Heat Transfer Computations for
Rearward-facing Steps and Sudden pipe Expansions. Heat transfer, 1985, 1:70~
76.
6 Chiang, T. P. W. h. Sheu, S. K. Wang. Side wall effects on the structure
of laminar flow over a plane-symmetric sudden expansion. Computers & Fluids,
2000, 29: 467~492
7 T. JONGEN and . Design of an unconditionally stable positive
scheme for the k-εand two-layer turbulence models. Computers & Fluids,
8 G. D. mallinson and G. DE vahl davis. Three-dimensinal natural
convection in a box: a numerical study, J. Fluid Mech. 1977, 83: 1~31
9 Latimer. B. R, Pollard. A. Comparison of pressure-velocity coupling
solution algorithm. Numer Heat Transfer, 1985, 8: 635~652
10 Latimer. B. Rodi and G. Scheuerer, Turbulence models for Near-wall and
low-Reynolds number flows: A Review. AIAA J .23, 1308~1319 (1985)
11 K. ABE and T. KONDOH, A new turbulence model for predicting fluid flow
and heat transfer in separating and reattaching flows---2. Thermal field
calculations, Int. J. Heat Mass Transfer 38, 1467~1481 (1995)
12 . Chinag, Tony W. h .Sheu. Side wall effects on the structure of
laminar flow over a plane-symmetric sudden expansion. Computers & Fluids,
2000, 29: 467~492