2. 武汉科技大学 省部共建耐火材料与冶金国家重点实验室,湖北 武汉 430081;
3. 大连理工大学 能源与动力工程学院,辽宁 大连 116024
2. The State Key Laboratory of Refractories and Metallurgy, Wuhan University of Science and Technology, Wuhan 430081, China;
3. School of Energy & Power Engineering, Dalian University of Technology, Dalian 116024, China
多孔介质内的自然流动现象广泛存在于自然界和工程技术领域,如:建筑物隔热、地下污染物迁移、谷物干燥、太阳能收集、填充床内催化反应放射性废物容器冷却等.二维多孔介质方腔内的自然对流模型是研究多孔介质热对流的基本模型.其特点是方腔左右壁面分别维持恒定的一高一低的温度,上下壁面为绝热壁面.由于此模型边界条件简单,因此常用来验证不同数值解法的优劣性.
近几十年来,不同的学者运用了不同的方法对多孔介质内自然对流问题进行了研究.文献[1-2]使用交替隐式有限差分方法研究了多孔介质倾斜方腔和梯形腔内自然对流,并且应用有限容积法研究了局部非热平衡下多孔介质方腔内的自然流动[3];文献[4]使用Galerkin有限元法研究了多孔介质波纹曲面腔内的自然流动;文献[5]应用稳定的有限元二阶投影法研究了随机构造的多孔介质内的自然流动;文献[6]开发了一种模拟多孔介质内自然流动的热格子-波尔兹曼模型等.不同的研究者采用不同的数值方法对物理问题进行研究,而据作者了解,目前并没有人采用谱方法对多孔介质内自然对流进行模拟.
谱方法的主要优点是精度高,是与有限差分、有限元并驾齐驱的求解偏微分方程的三大数值方法之一.其对光滑函数具有指数性逼近的谱精度, 可以以较少的网格点得到较高的精度, 具有谱解析性和全域性.因此,谱方法在流体流动[7]、辐射换热[8]、磁流体[9]、量子力学[10]等方面有了广泛的应用.本文主要应用Chebyshev配置点谱方法来模拟局部热平衡状态下多孔介质方腔内的自然流动,验证谱方法在该问题上的适用性及其精度,并在此基础上分析无量纲参数Ra对流动和温度分布的影响.
1 物理问题与数学模型本文所研究的物理模型与坐标系统如图 1所示,该模型中,二维方腔的左侧垂直壁面为维持恒定温度Th的热壁面,右侧垂直壁面为维持恒定温度Tc的冷壁面,且Th>Tc,上下壁面为绝热壁面.同时假设:流动为稳定的、不可压缩的层流流动,用Darcy模型描述多孔介质内动量方程;多孔介质内固体与流体处于局部热平衡状态;多孔介质结构均匀,具有各向同性渗透率;除密度以外,多孔介质各物性参数认为是常数.
根据以上假设,描述图 1所示物理问题的无量纲控制方程、初始条件和边界条件如下:
(1) |
(2) |
(3) |
式中:
其中:Ra为瑞利数;Ψ为无量纲的流函数;ν为运动黏度,m2·s-1;u为x方向速度,m·s-1;v为y方向速度, m·s-1;x, y为直角坐标系坐标,m;K为多孔介质渗透率,m2;g为重力加速度,m·s-2;β为热膨胀系数,K-1;T为温度,K;t为时间,s;α为热扩散率,m2·s-1;Γ为热扩散率比;X,Y分别为无量纲的x方向和y方向坐标.
2 Chebyshev配置点谱方法为了使数值求解过程更为稳定,在应用Chebyshev配置点谱方法求解控制方程及边界条件 (1)~(3) 时,在方程 (1) 上增加非稳态项 (∂Ψ/∂τ)[11].采用Chebyshev-Gauss-Lobatto (CGL) 配置点对无量纲化的控制方程进行空间上的离散.
2.1 数值求解过程对控制方程中非稳态项采用一阶有限差分格式进行离散,时间步长为Δτ.方程中的Laplacian项在时间上采用隐式 (n+1时层),其他项采用显式 (n时层) 进行处理.在Chebyshev配置点谱方法的理论中,无量纲流函数和温度可用配置点上的值由Lagrange插值多项式进行近似[12],将插值式代入控制方程,并引入初始条件与边界条件可得
(4) |
(5) |
这里,系数矩阵的各个元素的计算式分别为
(6a) |
(6b) |
(6c) |
(6d) |
(6e) |
(6f) |
式中:e1=Dy(1)(1, 1);e2=Dy(1)(Ny+1, 1);e3=Dy(1)(1, Ny+1);e4=Dy(1)(Ny+1, Ny+1);e=e1e4-e2e3; 矩阵D(1)和D(2)分别为一阶和二阶导数矩阵[13].
由于A, B, P和Q矩阵的特征值和特征向量均为实数,可采用高效的矩阵对角化方法[13]对方程组进行求解.求解过程中,时间步长Δτ=10-5,计算终止判据为前后时层待求参数之差小于10-8.
2.2 网格独立解测试对于二维多孔介质方腔,竖直壁面的平均努塞尔数为
(7) |
为了获得平均努塞尔数的网格无关解,选取了三套网格进行了测试,分别为24×24,30×30和36×36, 结果见表 1.可知,当网格数增大至36×36时,平均努塞尔数趋于稳定,最大的误差小于0.86%.以下的计算采用的节点数均为36×36.
表 2对当前数值算法所求出的竖直壁面努塞尔数与其他研究者所求出的结果进行了比较,从表中看出,当前数值解法与其他解法所得结果误差均小于3%.
为了测试本数值方法的精确性,构造了一个精确解对其进行了误差测试.
精确解的形式如式 (8) 所示,精确解对应的控制方程和边界条件详见作者前期工作[16].
(8) |
采用标准差E2进行测试:
式中:Nx为x方向上总的配置点数;Ny为y方向上总的配置点数;下标i, j分别为x和y方向配置点节点标识.
图 2和图 3分别为当X等于0.5时,不同的Y对应节点及Y等于0.5,不同的X对应节点上数值计算解与精确解之间的误差.从图中可以看出,其流函数最大误差为1.76×10-12,而温度最大误差为2.79×10-11.由此表明,采用谱方法计算二维多孔介质方腔的自然流动问题具有很高的计算精度.
不同Ra下,流场和温度场的分布如图 4所示,从图中可以看出,随着Ra的增大,流场内中心涡的形状发生变化,且流动逐渐加强;温度场内等温线发生偏移,等温线上部向冷壁面集聚,等温线下部则向热壁面集聚,致使等温线聚集的部位传热加强.
不同Ra下,竖直壁面平均努塞尔数的变化趋势如图 5所示.从图中可以看出,随着Ra的增大,平均努塞尔数增大,且二者在测试范围内 (10≤Ra≤1 500) 呈现线性变化的规律.
1) 本文构造了一个精确解对误差进行了测试,结果表明本方法具有很高的计算精度.
2) 在验证数值方法精度的基础上,分析了无量纲参数Ra对多孔介质内流场、温度场及竖直壁面平均努塞尔数的影响.发现,随着Ra的增加,腔体内流动得到加强;垂直壁面上的平均努塞尔数随着Ra的增大而增大,在10≤Ra≤1 500区间内,二者呈现线性变化的规律.
[1] | Baytas A C, Pop I. Free convection in oblique enclosures filled with a porous medium[J]. International Journal of Heat and Mass Transfer, 1999, 42(6): 1047–1057. DOI:10.1016/S0017-9310(98)00208-7 |
[2] | Baytas A C, Pop I. Natural convection in a trapezoidal enclosure filled with a porous medium[J]. International Journal of Engineering Science, 2001, 39(1): 125–134. |
[3] | Baytas A C, Pop I. Free convection in a square porous cavity using a thermal nonequilibrium model[J]. International Journal of Thermal Sciences, 2002, 41(6): 861–870. |
[4] | Misirlioglu A, Baytas A C, Pop I. Free convection in a wavy cavity filled with a porous medium[J]. International Journal of Heat and Mass Transfer, 2005, 48(9): 1840–1850. DOI:10.1016/j.ijheatmasstransfer.2004.12.005 |
[5] | Ma X, Zabaras N. A stabilized stochastic finite element second-order projection method for modeling natural convection in random porous media[J]. Journal of Computational Physics, 2008, 227(8): 8448–8471. |
[6] | Gao D Y, Chen Z Q, Chen L H. A thermal lattice Boltzmann model for natural convection in porous media under local thermal non-equilibrium conditions[J]. International Journal of Heat and Mass Transfer, 2014, 70: 979–989. DOI:10.1016/j.ijheatmasstransfer.2013.11.050 |
[7] | Subich C J, Lamb K G, Stastna M. Simulation of the Naviere-Stokes equations in three dimensions with a spectral collocation method[J]. International Journal for Numerical Methods in Fluids, 2013, 73(1): 103–129. |
[8] | Li B W, Sun Y S, Yu Y. Iterative and direct Chebyshev collocation spectral methods for one-dimensional radiative heat transfer[J]. International Journal of Heat and Mass Transfer, 2008, 51(5): 5887–5894. |
[9] | Zhang J K, Li B W, Chen Y Y. Hall effects on natural convection of participating MHD with thermal radiation in a cavity[J]. International Journal of Heat and Mass Transfer, 2013, 66: 838–843. DOI:10.1016/j.ijheatmasstransfer.2013.07.090 |
[10] | Li Z C, Chen S Y, Chen C S, et al. A spectral collocation method for a rotating Bose-Einstein condensation in optical lattices[J]. Computer Physics Communications, 2011, 182(6): 1215–1234. DOI:10.1016/j.cpc.2011.02.002 |
[11] | Khandelwal M K, Bera P, Chakrabarti A. Influence of periodicity of sinusoidal bottom boundary condition on natural convection in porous enclosure[J]. International Journal of Heat and Mass Transfer, 2012, 55(11): 2889–2900. |
[12] | Canuto C, Hussaini M Y, Quarteroni A, et al. Spectral methods:fundamentals in single domains[M]. Berlin: Springer, 2006. |
[13] | Peyret R. Spectral method for incompressible viscous flow[M]. New York: Springer, 2002. |
[14] | Manole D M, Lage J L. Numerical benchmark results for natural convection in a porous medium cavity[J]. ASME Publications Heat Transfer Division, 1993, 216: 55–55. |
[15] | Moya S L, Ramos E, Sen M. Numerical study of natural convection in a tilted rectangular porous material[J]. International Journal of Heat and Mass Transfer, 1987, 30(4): 741–756. DOI:10.1016/0017-9310(87)90204-3 |
[16] | Chen Y Y, Li B W, Zhang J K. Spectral collocation method for natural convection in a square porous cavity with local thermal equilibrium and non-equilibrium models[J]. International Journal of Heat and Mass Transfer, 2016, 96: 84–96. DOI:10.1016/j.ijheatmasstransfer.2016.01.007 |