搜索

x

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

二维湍流热对流最大速度Re数特性及流态突变特征Re

何建超 方明卫 包芸

引用本文:
Citation:

二维湍流热对流最大速度Re数特性及流态突变特征Re

何建超, 方明卫, 包芸

Scaling of Reynolds number based on maximum velocity and characteristic Reynolds number in two-dimensional thermal turbulence convection

He Jian-Chao, Fang Ming-Wei, Bao Yun
PDF
HTML
导出引用
计量
  • 文章访问数:  569
  • PDF下载量:  13
  • 被引次数: 0
出版历程
  • 收稿日期:  2022-02-28
  • 修回日期:  2022-06-09
  • 上网日期:  2022-09-16
  • 刊出日期:  2022-10-05

二维湍流热对流最大速度Re数特性及流态突变特征Re

  • 1. 北京航空航天大学, 航空科学与工程学院, 北京 100191
  • 2. 中山大学, 航空航天学院, 深圳 518107
  • 通信作者: 包芸, stsby@mail.sysu.edu.cn
    基金项目: 国家自然科学基金(批准号: 11772362)资助的课题.

摘要: 本文计算系列二维湍流热对流, Prandtl(Pr)数和Rayleigh(Ra)数范围分别为0.25—100和1×107—1×1012, 研究Reynolds(Re)数的变化规律. 以最大速度计算的Re数与Ra数存在标度律关系, 但中间出现间断. 研究表明, 大尺度环流形态由椭圆形到圆形的突变引起流动失稳, 导致最大速度值间断下降, 影响Re数变化趋势的连续性. 所有Pr数对应的流态突变特征Re数为常值, Rec约为1.4 × 104, 即当Re数达到特征Rec时, 大尺度环流形态会发生从椭圆形到圆形的突变. 间断点对应的RacPr数之间存在标度关系Rac-Pr1.5. 对Ra数进行补偿平移, 所有Pr数的ReRaPr–1.5的变化曲线重合, 不同Pr数有相同的间断临界点位置, RacPr–1.5 = 109.

English Abstract

    • 热对流现象广泛存在于自然界和工业设计中, 研究热对流特性有重要的意义. Rayleigh-Bénard(RB)热对流是从众多热对流过程中抽象出来的典型物理模型之一, 是在一个封闭的空间内下底板加热上底板冷却产生热对流运动和热输运的系统[1]. RB热对流系统存在丰富而复杂的流动和热输运现象, 一直受到国内外学者的关注和研究.

      在RB热对流系统中, Rayleigh数(Ra)、Prandtl数(Pr)和宽高比($\varGamma$)是控制系统的3个重要无量纲参数, 而Nusselt数(Nu)和Reynold数(Re)分别是反映系统传热效率以及湍流强度的无量纲参数. 多年来, 两个响应参数与控制参数之间的关系是RB热对流的研究重点[1-3], 而其中影响最广泛的是Grossmann和Lohse提出的GL理论[4-8]. GL理论在较大范围内预测了Nu (Ra, Pr)和Re(Ra, Pr)的变化趋势, 至今也得到了比较多数据的验证[8-10]. 在GL理论中, Ra-Pr相图被分为多个区域[4,8], 在不同区域中耗散率与Ra, Pr的关系有所差异. 也即, 在同一Pr数下, 随着Ra数增大, 耗散率的发展趋势会出现变化, 其他相关的物理量的行为也会有所变化[11]. 早期关于RB热对流的研究中, 发现到达临界RaRac (1707)之后, 流体运动会从一种定常(time-independent)状态转换变为随时间变化的非定常状态[11-15], 而这个转变的Ra数与Pr数存在依赖关系[12]. “芝加哥对流实验”[16-18]发现随着Ra数的增大, 系统的流态会从对流状态变为湍流状态, 同时系统的Nu数和Ra数关系也发生了转变. 其中, 当RB系统转变为湍流状态后, 系统中的主要结构, 包括羽流、大尺度环流等, 会随着Ra增大而发生突变[19-20]. 当大尺度环流由椭圆形转变为圆形时, 系统的Nu数也会偏离GL理论预测曲线[10], 但在更大的Ra后, Nu数的变化趋势再次与GL理论预测的一致[19-20].

      有数值模拟研究[21]表明不同Pr数下, RePr的标度律会出现变化, 而ReRa的标度律均为0.53, 与之前[18,22]获得的0.46不同. 最近, 在准二维的实验中Li等[23]Pr = 11.7—145.7的范围内发现Re的标度律在0.53—0.60之间, 与Chen等[24]Pr = 4—7之间发现的0.55十分接近. 在二维DNS中, Xu等[25]发现在Pr = 0.025时ReRa的标度律为0.50, 而Werne等[26]Pr = 7时的标度律为0.54. 这些研究表明, Re(Ra)的关系也会随Pr数变化而略有变化, 这变化值得去进一步深入研究.

      本文采用高效并行直接求解方法PDM-DNS[27], 在“天河二号”超级计算机上进行了多组Pr数和Ra数的二维湍流热对流DNS模拟, Pr数从0.25到100, Ra数范围为1 × 107—1 × 1012, 跨度为5个量级, 总共133个计算算例. 本文对多组Pr数和Ra数的二维湍流热对流中反映湍流特性的Re数变化规律进行研究, 并探讨热对流流态突变时对应的典型特征Ra数和Re数的特性及其变化规律.

    • 在RB热对流的数值模拟中, 通常引入Oberbeck-Boussinesq(OB)近似, 对方程进行简化. 在OB近似下, 无量纲化的热对流方程为

      $ \begin{split} &\nabla \cdot {\boldsymbol{u}}=0,\\ &\dfrac{D\boldsymbol{u}}{Dt}=-\nabla p+\dfrac{1}{\sqrt{Ra/Pr}}{\nabla }^{2}{\boldsymbol{u}}+\theta {\boldsymbol{k}},\\ &\dfrac{D\theta }{Dt}=\dfrac{1}{\sqrt{Ra\cdot Pr}}{\nabla }^{2}\theta ,\end{split} $

      其中, $\boldsymbol{u}$为无量纲速度矢量, $ \theta $为无量纲温度, $ p $为压力, $\boldsymbol{k}$为单位垂向矢量, 方向与重力方向相反. 数值计算中边界条件为壁面速度均采用无滑移条件, 温度为侧壁采用绝热条件, 上下底板采用恒温条件, 上底板恒温冷却$ {\theta }_{\mathrm{t}\mathrm{o}\mathrm{p}}=-0.5 $, 下底板恒温加热$ {\theta }_{\mathrm{b}\mathrm{o}\mathrm{t}}=0.5 $.

      大规模湍流热对流的DNS计算由于计算量巨大必须通过并行计算进行, 其中压力泊松方程的并行求解是实现高效并行计算的关键. 利用高效并行直接求解压力泊松方程的PDD算法[28], 建立了热对流DNS的并行直接求解方法(parallel direct method of DNS, PDM-DNS), 并展现了很好的并行效率[27]. 在数值求解过程中, 采用投影法进行计算, 时间和空间均是二阶精度; 压力泊松方程采用快速傅立叶变换(FFT)进行解耦, 然后直接求解三对角方程组, 而动量方程和温度方程采用显式格式推进.

      本文使用高效的PDM-DNS方法, 在“天河二号”超级计算机上进行了多组Pr数和Ra数的2D湍流热对流DNS模拟. 图1给出了计算算例的Ra-Pr相图. 本文中研究的Pr数从0.25到100, Ra数从1 × 107到1 × 1012, 跨度为5个量级, 总共133个计算算例. 本文的二维湍流热对流计算数据丰富, 相图中呈现的结果是目前二维RB系统较完整的数据.

      图  1  二维热对流计算算例Ra-Pr相图

      Figure 1.  The Ra-Pr Phase diagram explored in the present study.

      为了充分识别系统的流动, 本文采用了上下边界加密的网格进行数值模拟, 计算网格大小满足Kolmogorov尺度和Batchelor尺度, 即${\eta }_{K} = {HP{r}^{1/2}}/ {{\left[Ra\left(Nu-1\right)\right]}^{1/4}}$$ {\eta }_{B}={\eta }_{K}{Pr}^{1/2} $. 并且, 边界附近的网格满足Shishkina等[29]提出的要求, 即:

      $ \begin{split} &{N}_{\theta ,\mathrm{B}\mathrm{L}}\geqslant \sqrt{2}Nu{Pr}^{-0.5355+0.033\mathrm{l}\mathrm{o}\mathrm{g}Pr} , \\ &{N}_{v,\mathrm{B}\mathrm{L}}\geqslant \sqrt{2}N{u}^{1/2}{Pr}^{-0.1785+0.011\mathrm{l}\mathrm{o}\mathrm{g}Pr}, \end{split}$

      其中, a ≈ 0.482, $ {N}_{\theta , \mathrm{B}\mathrm{L}} $$ {N}_{v, \mathrm{B}\mathrm{L}} $分别表示温度边界层和速度边界层内的最少网格数. 时间步长Δt小于Kolmogorov时间尺度$ {\tau }_{K}=\sqrt{Pr/\left(Nu-1\right)} $的1/1000. 本文所有算例均在流场充分发展后进行统计, 统计时间均大于200个无量纲时间, 并且在不同统计时间段求出的Nu数误差约为1%[20].

      在本文的算例中, 有少数算例会出现大尺度环流翻转的情况, 具体包括Pr = 1.2, Ra = 2 × 107, Pr = 2.0, Ra = 1 × 107—1 × 108, Pr = 4.3, Ra = 5 × 107—2 × 108, Pr = 100, Ra = 2 × 108, 共9个算例. 当仅统计顺时针旋转或逆时针旋转的Re数和Nu数时, 两者的结果是一致的. 因此, 可以将顺时针旋转的数据进行水平翻转, 叠加到逆时针的旋转的数据上, 得到时间更长的时间平均场. 如果统计数据包括大尺度环流的翻转过程, 并且不作方向统一的处理, 那么将对数据统计造成一定的影响, 比如Re数会减小约50%, Nu数波动约为3%[30]. 因此, 关于翻转算例的统计, 仅统计顺时针和逆时针的流动, 并作方向统一处理, 不考虑发生大尺度环流翻转的过程.

    • Re数是反映热对流系统的湍流流动特性的特征参数. 本文研究系列Pr数及Ra数下的二维RB系统湍流Re数与Ra数和Pr数之间的变化规律及特性.

      图2给出了一个典型流态的平均速度场云图. 图2可见, 该速度分布成圆环状, 方腔四个角落中的速度都较小, 速度分布在中心速度很小, 随着半径的变化速度快速增大, 到接近一半的位置速度达到最大, 而后速度随半径的增大速度减小.

      图  2  平均速度分布及最大速度点

      Figure 2.  The time-averaged velocity filed and the maximum velocity point.

      经过计算全场的速度大小, 找到绝对速度的最大值, 即图2中白点处的速度, 作为计算湍流热对流Re数的特征速度, 其无量纲的计算公式如下:

      $ Re=\sqrt{Ra/Pr}\cdot {\left(\sqrt{{{u}}^{2}+{v}^{2}}\right)}_{\mathrm{max}}. $

    • 首先, 以Pr = 0.7的系列Ra数热对流为典型计算结果, 探讨热对流的Re数随Ra数的变化规律.

      图3给出了Pr = 0.7的系列Ra数热对流的Re数随Ra数的变化情况, 图中为双对数坐标. 可以看到, Re数随Ra数的增加逐渐增大, 表明系统的湍流强度随Ra数增加逐渐变强. 特别地, Re数的变化在Ra = 109处出现了间断平移, 分成了较明显的两段.

      图  3  Re数与Ra数关系

      Figure 3.  Re as a function of Ra

      对于Re数与Ra数的变化关系中出现的间断平移, 是由二维湍流热对流的流态突变造成的, 将在下节内容中详细讨论.

    • 在低Ra数时, 二维湍流热对流的流态是倾斜的椭圆大尺度环流加两个角涡的流态. 随Ra数增高, 角涡会发生脱离导致流动失稳, 系统的流态突变成圆形的大尺度环流流态[19].

      图4给出了Pr = 0.7时流态从椭圆突变成圆形前后典型流态的平均速度场云图, 图中明显可见两种完全不同的速度分布. 由于不同Ra数时速度值得变化范围较大, 为了清楚地反映速度分布情况, 每个Ra数对应的速度云图色标范围是不一致的, 均设置为各个算例的最小速度到最大速度. 在Ra数较低的图4(a)(b)中, 速度分布形态基本一致, 均为倾斜的椭圆形和两个角涡, 图中白点为最大速度点, 都出现在椭圆和角涡的相切处. 大尺度环流与角涡相切的位置速度都会比较大, 而接近另外两个角落的速度较小. 随着Ra数增高, $ Ra= 1\times 1{0}^{9} $时流态发生突变, 大尺度环流形态突变为圆形, 如图4(c)(d)所示. 突变后高速区域呈圆环状, 整体速度相近, 最大速度随机产生在圆环上, 并且最大速度相较突变前略有减小, 这一点在图4的色标范围上可以看出.

      图  4  平均速度场变化及流态突变(Pr = 0.7) (a) Ra = 2 × 108; (b) Ra = 5 × 108; (c) Ra = 1 × 109; (d) Ra = 2 × 109

      Figure 4.  The time-averaged velocity fieldsc and the sudden change of flow pattern (Pr = 0.7): (a) Ra = 2 × 108; (b) Ra = 5 × 108; (c) Ra = 1 × 109; (d) Ra = 2 × 109.

      计算Pr = 0.7时所有Ra数情况下的最大速度值, 探讨最大速度值随Ra数的变化规律.

      图5给出了Pr = 0.7时不同Ra数的平均速度场中最大速度值UmaxRa数的变化情况. 可以看到, 最大速度的变化明显分为两段. 当$Ra\leqslant 5\times 1{0}^{8}$, 大尺度环流流态为椭圆, 随着Ra数的提高平均速度场中的最大值Umax逐渐变大. 在$ Ra= 1\times 1{0}^{9} $时, 大尺度环流流态突变为圆形, 最大速度Umax间断式突然下降. 之后随着Ra数增大, 平均速度场中的最大值Umax再次逐渐变大. 当$ Ra\geqslant 5\times 1{0}^{10} $后, 最大速度Umax出现波动, 但速度值变化不大.

      图  5  最大速度随Ra数变化

      Figure 5.  The maximum velocity as a function of Ra

      二维湍流热对流在较低Ra数时呈椭圆形大尺度环流形态, 由于羽流基本沿椭圆路径运动, 流动相对集中稳定, 使平均场速度较大. 随着Ra数的提高, 速度增大, 带动角涡脱落造成流态失稳混乱, 导致流态从椭圆到圆形的突变[31]. 当流态失稳突变后, 羽流运动呈现较为混乱的绕行, 会使平均后的速度场整体减小, 所以导致图5中最大速度间断式减小.

      热对流在流态突变时平均场最大速度Umax的间断式减小, 是造成图3Re数随Ra数的变化规律出现间断的原因.

    • 本文计算了10个不同的Pr数系列的算例, 数据较多. 为清晰展示结果, 在讨论二维湍流热对流Re数特性的Pr数影响时, 首先选取3个典型的Pr数0.7, 4.3和10进行分析.

      图6给出了3个Pr数系列的Re数的变化特性, 蓝色六边形、红色方形和绿色实心圆点分别为${Pr}={0.7}$, 4.3和10, 黑色线表示$ 1.4\times 1{0}^{4} $. 由图6中可见, 3个Pr数系列的Re数都分别随Ra数的增大而增大, 同时Re数随Ra数分布有两个阶段, 阶段之间存在明显的间断. 将间断处的Re数定义为流态突变特征ReRec, Ra数定义为流态突变特征RaRac. 有意思的是, 不同Pr数时流态突变对应的特征ReRec的值基本一致($ R{e}_{\mathrm{c}}\approx 1.4\times 1{0}^{4} $). 也即, 无论Pr数为多少, 当Re数达到Rec时, 平均场的大尺度环流形态一定会从椭圆形突变为圆形.

      图  6  Re数与Ra的关系图

      Figure 6.  Re as a function of Ra

      这一现象的原因与上文讨论的Pr = 0.7的一样, 是由于最大速度值在流态突变处出现骤减, 导致Re数变化出现间断平移. 另外, 随着Pr增大, 发生间断平移的位置会向右移, 表明随Pr增大, Rac会逐步增大.

      计算所有Pr数下二维湍流热对流的Re数, Pr数范围为0.25—100, 讨论这个范围内不同Pr数情况下, 热对流的Re数与Rac的变化规律.

      图7给出相应的Pr从0.25到100的系列Ra数下Re数分布. 图中发现不同Pr数下湍流Re数分布特征与之前3个典型Pr数情况相似. Re数从一百左右到几十万, 跨度4个量级. 除了Pr = 50和100时Re数数值不够大, 其他的结果都在流态突变特征Re数处存在间断, 间断点对应的Ra数依次往后移动. 然而间断处的Re数, 也就是流态突变特征Re数的值基本相同, 是一个常数, $ {Re}_{\mathrm{c}}\approx 1.4\times 1{0}^{4} $. 也即, 对于不同系列Pr数和Ra数, 当由最大速度定义的Re数达到特征Rec时, 流态就会发生由椭圆型大尺度环流到圆形大尺度环流的流态突变, 流体突变特征ReRecPr数和Ra数变化无关.

      图  7  不同Pr数下Re数与Ra数的关系

      Figure 7.  Re as a function of Ra at different Pr

      在不同Pr数下, Rac的位置逐渐右移, 表明RacPr增大而增大. 在二维RB热对流中, 流态突变Ra数与Pr数有Ra-Pr1.5的标度律关系[30]. 为了探讨不同Pr数下Re数的转变共同特征, 对横坐标Ra数进行补偿平移, 研究ReRaPr–1.5的关系.

      图8给出Re数随RaPr –1.5的分布. 图中可见, 10个Pr数的Re(RaPr –1.5)分布完全重合在一起, 表明在不同Pr数下, Re数的变化规律基本一致, 呈现很好的自似性.

      图  8  不同Pr数下Re数与RaPr–1.5的关系

      Figure 8.  Re as a function of RaPr–1.5 at different Pr

      图8中黑色虚线表示$ {Re}_{\mathrm{c}}\approx 1.4\times 1{0}^{4} $, 红色虚线表示$ R{a}_{\mathrm{c}}{Pr}^{-1.5}={10}^{9} $. 两条虚线的交点对应不同Pr数下的流态突变. 这表明不同Pr数都有基本相同的流态突变特征点位置, $ R{a}_{\mathrm{c}}{Pr}^{-1.5}={10}^{9} $, 同时$ Re={Re}_{\mathrm{c}}\approx 1.4\times 1{0}^{4} $. 这是两个常数, 可以用于预测不同Pr数和Ra数下流态发生突变的时机以及系统Re数的变化趋势, 对区分不同参数下的流态以及流态突变特性的探讨有重要意义. 当$ Ra{Pr}^{-1.5} < {10}^{9} $时, 系统的流态为椭圆形; 当$ Ra{Pr}^{-1.5}= {10}^{9} $时, 系统的流态会发生突变, 系统的Re数达到$ {Re}_{\mathrm{c}}\approx 1.4\times 1{0}^{4} $; 随后$ Ra{Pr}^{-1.5} > {10}^{9} $, 系统的流态为圆形.

      另外, 对图8中所有数据进行Re-(RaPr–1.5)γ的拟合, 得出: 流态突变前, Re数的标度律为Re-Ra0.61Pr–0.92, 流态突变后为Re-Ra0.55Pr–0.83.

      流态突变的临界Re数为常数, 反映出二维湍流热对流特殊的流动形态变化特征, 也将为二维湍流热对流的流动特性研究以及对应的传热特性变化特征等问题的研究提供新的思路. 更多的价值还需要进一步深入的研究.

    • 本文采用高效并行直接求解计算方法PDM-DNS完成了系列Pr数和Ra数的二维湍流热对流的DNS模拟, Pr数从0.25到100, Ra数从1 × 107到1 × 1012, 跨度为5个量级, 总共133个计算算例, 所呈现的结果是目前二维湍流热对流系统相当完整的数据. 本文研究以平均场最大速度为特征速度的Re数特性以及最大速度的变化规律, 发现大尺度环流形态由椭圆形变为圆形的突变对Re数特性的影响. 研究结论如下:

      1)典型算例Pr = 0.7时的Re数特性结果表明, ReRa数的变化存在标度律关系, 但中间出现明显的间断现象. 从基本流态特征出发进行研究发现, 大尺度环流形态由椭圆形突变为圆形会引起最大速度值的间断式突然下降, 导致Re数特性的间断现象出现.

      2)不同Pr数的Re数特性研究表明, Re数随Ra数的变化特性有明显的自相似性. 并且发现, 流态突变对应的特征ReRec为常数, 与Pr数和Ra数变化无关, $ {Re}_{\mathrm{c}}\approx 1.4\times {10}^{4} $. 当Re数达到特征Rec时, 大尺度环流形态会发生突变, 从椭圆形变为圆形. 这为RB热对流的流态区分提供新的方法, 并且进一步根据Re数的自相似性进一步预测Re数的变化规律.

      3)流态突变特征RaRacPr数之间存在标度关系Rac-Pr1.5. 对横坐标Ra数进行补偿平移, ReRaPr–1.5的变化曲线完全重合, 表明不同Pr数时Re数随Ra数的变化规律基本一致, 流态突变前后分别有Re-Ra0.61Pr–0.92Re-Ra0.55Pr–0.83的标度律关系. 不同Pr数有基本相同的间断特征点位置, $ R{a}_{\mathrm{c}}{Pr}^{-1.5}={10}^{9} $.

参考文献 (31)

目录

    /

    返回文章
    返回