﻿ 不同冷热壁面对方腔内对流换热影响的格子Boltzmann模拟
 上海理工大学学报  2023, Vol. 45 Issue (5): 448-459 PDF

1. 上海理工大学 能源与动力工程学院，上海 200093;
2. 上海建桥学院，上海 201306

Lattice Boltzmann simulation of convection heat transfer in a cubic cavity with different local cold and hot wall conditions
WANG Jie1, ZHANG Shuanyang1, XU Hongtao1, YANG Mo1,2
1. School of Energy and Power Engineering, University of Shanghai for Science and Technology, Shanghai 200093, China;
2. Shanghai Jian Qiao University, Shanghai 201306, China
Abstract: Based on the lattice Boltzmann method, the multiple-relaxation-time (MRT) model was selected to simulate the three-dimensional natural convection heat transfer in a cubic cavity with local hot and cold walls. The D3Q19 model and the D3Q7 model were used to describe the velocity field and temperature field, respectively. The effects of different Rayleigh numbers Ra (103Ra≤105) and local hot and cold wall position on the natural convection in the three-dimensional cubic cavity were studied. The results reveal that the arrangement of the hot and cold walls has a significant effect on the flow heat transfer, and the larger the Ra, the more obvious the effect. The strongest natural convection heat transfer capacity is obtained by Case 5 with the hot and cold walls in the middle position. As the Ra increases gradually, the natural convection ability increases and the streamline becomes more complicated. Under the same working case, with the increase of Ra, the average Nusselt number Nuav also increases gradually, and the heat transfer capacity is stronger. The Nuav value of Case 5 is the largest at the same Ra, indicating that its natural convection heat transfer capability is the strongest.
Key words: Boltzmann simulation     multiple-relaxation-time model     natural convection     3D cubic cavity

1 计算模型 1.1 物理模型

 图 1 立方腔体中具有不同冷热壁面位置的物理模型 Fig. 1 Physical models with different hot and cold wall positions in a cubic cavity
1.2 数学模型

 $\left\{ \begin{gathered} \frac{{\partial U}}{{\partial X}} + \frac{{\partial V}}{{\partial Y}} + \frac{{\partial W}}{{\partial Z}} = 0 \\ \frac{{\partial U}}{{\partial \tau }} + U\frac{{\partial U}}{{\partial X}} + V\frac{{\partial U}}{{\partial Y}} + W\frac{{\partial U}}{{\partial Z}} = - \frac{{\partial P}}{{\partial X}} + Ma\sqrt {\frac{{Pr}}{{3Ra}}} \left( {\frac{{{\partial ^2}U}}{{\partial {X^2}}} + \frac{{{\partial ^2}U}}{{\partial {Y^2}}} + \frac{{{\partial ^2}U}}{{\partial {Z^2}}}} \right) \\ \frac{{\partial V}}{{\partial \tau }} + U\frac{{\partial V}}{{\partial X}} + V\frac{{\partial V}}{{\partial Y}} + W\frac{{\partial V}}{{\partial Z}} = - \frac{{\partial P}}{{\partial Y}} + Ma\sqrt {\frac{{Pr}}{{3Ra}}} \left( {\frac{{{\partial ^2}V}}{{\partial {X^2}}} + \frac{{{\partial ^2}V}}{{\partial {Y^2}}} + \frac{{{\partial ^2}V}}{{\partial {Z^2}}}} \right) \\ \frac{{\partial W}}{{\partial \tau }} + U\frac{{\partial W}}{{\partial X}} + V\frac{{\partial W}}{{\partial Y}} + W\frac{{\partial W}}{{\partial Z}} = - \frac{{\partial P}}{{\partial Z}} + Ma\sqrt {\frac{{Pr}}{{3Ra}}} \left( {\frac{{{\partial ^2}W}}{{\partial {X^2}}} + \frac{{{\partial ^2}W}}{{\partial {Y^2}}} + \frac{{{\partial ^2}W}}{{\partial {Z^2}}}} \right) + \frac{{M{a^2}}}{3}\varTheta \\ \frac{{\partial \varTheta }}{{\partial \tau }} + U\frac{{\partial \varTheta }}{{\partial X}} + V\frac{{\partial \varTheta }}{{\partial Y}} + W\frac{{\partial \varTheta }}{{\partial Z}} = Ma\sqrt {\frac{1}{{3RaPr}}} \left( {\frac{{{\partial ^2}\varTheta }}{{\partial {X^2}}} + \frac{{{\partial ^2}\varTheta }}{{\partial {Y^2}}} + \frac{{{\partial ^2}\varTheta }}{{\partial {Z^2}}}} \right) \\ \end{gathered} \right.$ (1)

 $\begin{split}& X=\frac{x}{H}\text{，}Y=\frac{y}{H}\text{，}Z=\frac{z}{H}\text{，}{u}_{{\rm{c}}}=\sqrt{g\beta \left({T}_{{\rm{h}}}-{T}_{{\rm{c}}}\right)H}\text{，}\\&Ma=\frac{{u}_{{\rm{c}}}}{{c}_{{\rm{s}}}}\text{，}U=\frac{u}{\sqrt{3}{c}_{{\rm{s}}}}\text{，}V=\frac{v}{\sqrt{3}{c}_{{\rm{s}}}},\; W=\frac{w}{\sqrt{3}{c}_{{\rm{s}}}}\text{，}\\&Pr=\frac{\upsilon }{a}\text{，}P=\frac{p}{3\rho {c}_{{\rm{s}}}^{2}}\text{，}\tau =\frac{t\sqrt{3}{c}_{{\rm{s}}}}{H}\text{，}\\&Ra=\frac{g\beta ({T}_{{\rm{h}}}-{T}_{{\rm{c}}}){H}^{3}Pr}{{\upsilon }^{2}}\text{，}\varTheta =\frac{T-{T}_{{\rm{c}}}}{{T}_{{\rm{h}}}-{T}_{{\rm{c}}}}\\[-12pt] \end{split}$ (2)

 $\begin{split} N u=&\left.\frac{h}{k / H}\right|_{Y=0, X \in[0,1], Z \in[0,1]}=\\&-\left.\frac{\partial \varTheta}{\partial Y}\right|_{Y=0, X \in[0,1], Z \in[0,1]} \end{split}$ (3)

 $N{u_{{\rm{av}}}} = \int_{{Z_1}}^{{Z_2}} {\int_0^1 {Nu{\rm{d}}X{\rm{d}}Z} },\;\; {\text{ }}{{{Z}}_{\text{1}}}{\text{,}}\;{{{Z}}_{\text{2}}} \in \left[ {0,1} \right]$ (4)

2 MRT-LBM模型

 $\begin{split}& {f_i}\left( {{\boldsymbol{r}} + {{\boldsymbol{e}}_i}\Delta t,t + \Delta t} \right) - {f_i}\left( {{\boldsymbol{r}},t} \right) = {{\boldsymbol{\varOmega}} _{{\text{ }}i}} + {F_i}{\text{ }},\\&\qquad \qquad i = 1,2, \cdots ,19 \end{split}$ (5)

 $\boldsymbol{e}_i=c\left[\begin{array}{ccccccccccccccccccc} 0 & 1 & -1 & 0 & 0 & 0 & 0 & 1 & -1 & 1 & -1 & 1 & -1 & 1 & -1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & -1 & 0 & 0 & 1 & 1 & -1 & -1 & 0 & 0 & 0 & 0 & 1 & -1 & 1 & -1 \\ 0 & 0 & 0 & 0 & 0 & 1 & -1 & 0 & 0 & 0 & 0 & 1 & 1 & -1 & -1 & 1 & 1 & -1 & -1 \end{array}\right]$ (6)
 ${F_{{i}}}{\text{ = }} \Delta t{\boldsymbol{G}}\frac{{({{\boldsymbol{e}}_i} - {\boldsymbol{V}})}}{p}f_{{\rm{eq}},i}({\boldsymbol{r}},t)$ (7)
 ${\boldsymbol{G}} = - \beta \left( {T - {T_0}} \right){\boldsymbol{g}}{\text{ }}$ (8)
 $\begin{split}& f_{{\rm{eq}},i}({\boldsymbol{r}},t) = \rho {\omega _i}\left[ {1 + \frac{{{{\boldsymbol{e}}_i}{\boldsymbol{V}}}}{{c_{\rm{s}}^2}} + \frac{{{{({{\boldsymbol{e}}_i}{\boldsymbol{V}})}^2}}}{{2c_{\rm{s}}^2}} - \frac{{{\boldsymbol{VV}}}}{{2c_{\rm{s}}^2}}} \right]{\text{ , }}\\&\qquad i = 1,2, \cdots ,19 \end{split}$ (9)
 ${\omega _i} = \left\{ \begin{gathered} 1/3,{\text{ }}i = 1 \\ 1/18,{\text{ }}i = 2,3, \cdots ,7 \\ 1/36,{\text{ }}i = 8,9, \cdots ,19 \\ \end{gathered} \right.$ (10)
 $\begin{split}& {{\boldsymbol{\varOmega}} _{{\text{ }}i}} = - {{\boldsymbol{M}}^{ - 1}}{\boldsymbol{S}}\left[ {{m_i}\left( {{\boldsymbol{r}},t} \right) - m_{{\rm{eq}},i}\left( {{\boldsymbol{r}},t} \right)} \right]{\text{ ,}}\\&\;\;\qquad \quad i = 1,2, \cdots ,19 \end{split}$ (11)

 $\left\{ \begin{gathered} \rho {\text{ = }}\sum\limits_{i = 1}^{19} {{f_i}\left( {{\boldsymbol{r}},t} \right),{\text{ }}\quad\rho {\boldsymbol{V}} = \sum\limits_{i = 1}^{19} {{f_i}\left( {{\boldsymbol{r}},t} \right)} } \\ {\boldsymbol{m}} = {\boldsymbol{M}}{\left[ {{f_1},{\text{ }}{f_2},{\text{ }} \cdots, {f_{19}}} \right]^{\rm{T}}} \\ \end{gathered} \right.$ (12)

 $\begin{split}& {g_i}\left( {{\boldsymbol{r}} + {{\boldsymbol{u}}_i}\Delta t,t + \Delta t} \right) - {g_i}\left( {{\boldsymbol{r}},t} \right) =\\&\quad - {{\boldsymbol{N}}^{{{ - }}1}}{\boldsymbol{Q}}\left[ {{n_i}\left( {{\boldsymbol{r}},t} \right) - n_{{\rm{eq}},i}\left( {{\boldsymbol{r}},t} \right)} \right]{\text{, }}i = 1,2, \cdots ,7 \end{split}$ (13)

 ${{\boldsymbol{u}}_i} = c\left[ \begin{gathered} 0\quad \;\; 1\quad \;\; - 1\quad \;\; 0\quad \;\; 0 \quad \;\; 0 \quad \;\; 0 \\ 0\quad \;\; 0\quad \;\; 0\quad \;\; 1\quad \;\; - 1\quad \;\; 0\quad \;\; 0 \\ 0\quad \;\; 0\quad \;\; 0\quad \;\; 0\quad \;\; 0\quad \;\;1 \quad \;\; - 1 \\ \end{gathered} \right]$ (14)

 ${\omega _{T,i}} = \left\{ \begin{gathered} 1/4,{\text{ }}i = 1 \\ 1/8,{\text{ }}i = 2,3, \cdots ,7 \\ \end{gathered} \right.$ (15)

 $\left\{ \begin{gathered} T = \sum\limits_{i = 1}^7 {{g_i}\left( {r,t} \right)} \\ {\boldsymbol{n}} ={\boldsymbol{ N}}{\left[ {{g_1},{g_2}, \cdot \cdot \cdot {g_7}} \right]^{\rm{T}}} \\ \end{gathered} \right.$ (16)

Li等[12]通过研究表明，当网格数为50×50×50时，LBM在模拟三维方腔对流问题的研究中具有较高的准确性。因此，本文同样选取50×50×50的网格数来进行模拟研究。

3 程序验证

4 结果与分析 4.1 工况1和工况2的特性分析

 图 2 工况1和工况2的温度分布云图 Fig. 2 Isotherm results of cases 1 and 2

 图 3 工况1和工况2的流线分布（续） Fig. 3 Streamline results of case 1 and 2
4.2 工况3和工况4的特性分析

 图 4 工况3和工况4的温度分布云图 Fig. 4 Isotherm results of cases 3 and 4

 图 5 工况3和工况4的流线分布（续） Fig. 5 Streamline results of case 3 and 4
4.3 工况5的特性分析

 图 6 工况5的温度分布云图 Fig. 6 Isotherm results of case 5

 图 7 工况5的流线分布 Fig. 7 Streamline results of case 5
4.4 不同冷热壁面位置对平均努塞尔数的影响

5 结　论

a. 相同Ra条件下，冷热壁面的相对布置方式对流动换热有重要影响。工况1的对流换热能力最差，工况5的自然对流换热能力最强。在相同工况条件下，随着Ra从103逐渐增大105，自然对流能力增强，流线变得更为复杂。

b. 在相同工况条件下，随着Ra增大，Nuav也逐渐增加，换热能力增强。在相同Ra条件下，工况5的Nuav均最大，表明其自然对流换热能力最强，说明冷热壁面的布置方式减小了上下壁面对流体流动的抑制作用，使浮升力的作用效果更为明显。

 [1] 杨世铭, 陶文铨. 传热学[M]. 4版. 北京: 高等教育出版社, 2006. [2] SHEN R Q, JIAO Z R, PARKER T, et al. Recent application of computational fluid dynamics (CFD) in process safety and loss prevention: a review[J]. Journal of Loss Prevention in the Process Industries, 2020, 67: 104252. DOI:10.1016/j.jlp.2020.104252 [3] 张拴羊, 徐洪涛, 梁天生, 等. Soret和Dufour效应对方腔内双扩散自然对流影响的格子Boltzmann模拟[J]. 空气动力学学报, 2019, 37(2): 207-215. [4] 余端民, 赵明. 圆内开缝圆环形空间自然对流的格子Boltzmann模拟[J]. 上海理工大学学报, 2019, 41(2): 108-115. DOI:10.13255/j.cnki.jusst.2019.02.002 [5] XU H T, LUO Z Q, LOU Q, et al. Lattice Boltzmann simulations of the double-diffusive natural convection and oscillation characteristics in an enclosure with Soret and Dufour effects[J]. International Journal of Thermal Sciences, 2019, 136: 159-171. DOI:10.1016/j.ijthermalsci.2018.10.015 [6] DADVAND A, SARAEI S H, GHOREISHI S, et al. Lattice Boltzmann simulation of natural convection in a square enclosure with discrete heating[J]. Mathematics and Computers in Simulation, 2021, 179: 265-278. DOI:10.1016/j.matcom.2020.07.025 [7] 朱建奇, 侯佳煜, 杲东彦, 等. 基于格子Boltzmann方法的倾斜方腔自然对流模拟[J]. 南京师范大学学报(工程技术版), 2018, 18(4): 19-26. [8] FARKACH Y, DERFOUFI S, AHACHAD M, et al. Numerical investigation of natural convection in concentric cylinder partially heated based on MRT-lattice Boltzmann method[J]. International Communications in Heat and Mass Transfer, 2022, 132: 105856. DOI:10.1016/j.icheatmasstransfer.2021.105856 [9] SHEIKHOLESLAMI M, ELLAHI R. Three dimensional mesoscopic simulation of magnetic field effect on natural convection of nanofluid[J]. International Journal of Heat and Mass Transfer, 2015, 89: 799-808. DOI:10.1016/j.ijheatmasstransfer.2015.05.110 [10] KARATAS H, DERBENTLI T. Natural convection in rectangular cavities with one active vertical wall[J]. International Journal of Heat and Mass Transfer, 2017, 105: 305-315. DOI:10.1016/j.ijheatmasstransfer.2016.09.100 [11] ZHAO Y, WANG L, CHAI Z H, et al. Comparative study of natural convection melting inside a cubic cavity using an improved two-relaxation-time lattice Boltzmann model[J]. International Journal of Heat and Mass Transfer, 2019, 143: 118449. DOI:10.1016/j.ijheatmasstransfer.2019.118449 [12] LI Z, YANG M, ZHANG Y W. Lattice Boltzmann method simulation of 3-D natural convection with double MRT model[J]. International Journal of Heat and Mass Transfer, 2016, 94: 222-238. DOI:10.1016/j.ijheatmasstransfer.2015.11.042 [13] LI L K, MEI R W, KLAUSNER J F. Boundary conditions for thermal lattice Boltzmann equation method[J]. Journal of Computational Physics, 2013, 237: 366-395. DOI:10.1016/j.jcp.2012.11.027 [14] TRIC E, LABROSSE G, BETROUNI M. A first incursion into the 3D structure of natural convection of air in a differentially heated cubic cavity, from accurate numerical solutions[J]. International Journal of Heat and Mass Transfer, 2000, 43(21): 4043-4056. DOI:10.1016/S0017-9310(00)00037-5 [15] PENG Y, SHU C, CHEW Y T. A 3D incompressible thermal lattice Boltzmann model and its application to simulate natural convection in a cubic cavity[J]. Journal of Computational Physics, 2004, 193(1): 260-274. DOI:10.1016/j.jcp.2003.08.008