弱Soret效应混合流体对流系统的分岔与非线性演化*

郑来运 赵秉新 杨建青

1) (宁夏大学机械工程学院, 银川 750021)

2) (宁夏大学数学统计学院, 银川 750021)

3) (宁夏科学工程计算与数据分析重点实验室, 银川 750021)

混合流体Rayleigh-Bénard(RB)对流是研究非平衡耗散系统的自组织斑图及非线性动力学特性的典型模型.本文利用高精度数值方法模拟了底部均匀加热的矩形腔体中混合流体RB对流, 研究了具有极微弱Soret效应(分离比 ψ =−0.02 )的混合流体对流的分岔特性及斑图的形成和演化, 给出了分岔曲线图.获得了行波交替闪动的Blinking状态、局部行波对流和定常对流(SOC)三种稳定状态, 讨论了状态之间的过渡.研究发现从Blinking状态到局部行波对流状态的过渡存在迟滞现象, 过渡时行波频率、对流振幅和对流传热Nusselt数等均有明显的跳跃.在Blinking状态存在的Rayleigh数区间下界附近, 外部施加的不对称初始扰动是形成该状态的诱因.随着Rayleigh数增大, 临界SOC状态经过多次分岔并形成多个具有不同波数的SOC状态后过渡为混沌状态.

1 引 言

非平衡系统在远离热平衡态时就会产生复杂的自组织时空结构, 对其机理和稳定性等的研究是目前涉及非平衡现象科学技术问题的重要课题之一[1].混合流体Rayleigh-Bénard(RB)对流作为研究非平衡耗散系统自组织斑图及非线性动力学特性的一个典型系统, 近年来得到了广泛的关注.由于温度和浓度场之间的Soret耦合效应, 混合流体RB对流中的斑图结构比单组分流体中更为有趣和复杂.特别令人关注的是具有负分离比( ψ <0 )的混合流体中发生在对流临界点附近的对流状态, 这里 ψ 值表征了Soret耦合效应的强度, 其符号决定了混合流体的行为.ψ <0 表示热对浮力的贡献和浓度对浮力的贡献彼此相对, 导致系统经历亚临界Hopf分岔并最终发展成为各种可能的对流状态.早期的实验研究取得了很多重要的成果, 如首次观察到了局部行波对流状态(localized traveling wave, LTW)[2,3]和行波左右交替闪动的Blinking状态[4,5]等, 这些原始发现对之后的研究工作具有重要的意义, 但其机理并不十分清楚.为了深入理解和探讨复杂流动斑图及其非线性特性, 人们开展了基于流体力学基本方程组或其扰动方程组的直接数值摸拟研究, 分析了行波对流状态(traveling wave, TW)的参数依赖性[6], 解释了LTW漂移缓慢的原因[7,8], 讨论了Blinking状态以及“重复过渡”状态(repeated transients)的起源以及长高比的影响[9], 也获得了局部定常对流 (localized stationary overturning convection, LSOC)[10]、 摆动行波对流 (undulation TW, UTW)[11]和含缺陷行波对流[12,13]等有趣的流动斑图.近年来,Mercader等[14]在中等长高比的腔体中, 采用不同的边界条件观察到了多种LSOC对流结构及其多重稳定现象, 并发现当上壁面边界改变时奇宇称LSOC状态整体以一定的速度移动, 且速度大小依赖于Biot数[15].王涛等[16]采用高阶紧致差分方法求解涡量-流函数型的扰动方程, 模拟了矩形长腔内 LTW的形成过程.Watanabe等[17]和 Taraut等[18]模拟了大长高比矩形腔体内LTW状态之间以及与LSOC状态碰撞的情形, 发现碰撞后的对流结构完全依赖于碰撞前LTW状态的特性.赵秉新[19]对水平流作用下具有正的小分离比的混合流体行进波对流进行了模拟研究, 讨论了物性参数对流动的影响.最近, Shevtsova 等[20]介绍了在国际空间站进行的实验结果, 与模拟结果进行了对比,反驳了微重力环境总是影响扩散控制过程的推测,并讨论了实验装置上下壁面温度不均匀对Soret效应的影响[21]; Alonso等[22]模拟了准垂直圆柱体内具有正Soret系数的混合流体对流, 重点讨论了小倾角对斑图形成的影响; Smorodin等[23]通过线性稳定性分析, 研究了高频振动下水平层中混合流体的对流, 分析了振动强度和方向对空间周期行波解的影响; 学者们还模拟研究了流体参数和长宽比[24]、实验装置的倾斜角度[25]对流动的影响, 讨论了从上向下加热的负分离比混合流体对流的开始和非线性状态[26], 以及分离比和长高比对缺陷结构的影响[27,28].

可是, 关于混合流体对流中解的分岔及流动沿着分岔图上部分支的演化过程则研究比较少.文献[24]利用固壁边界条件对矩形腔体内分离比为ψ=−0.1的混合流体进行了模拟研究, 给出了解的分岔图; 文献 [29,30]对 ψ =−0.2 和–0.6 的混合流体进行数值模拟, 获得了几种行波解.对于极微弱分离比下混合流体对流中解的分岔及沿分岔分支演化过程的研究则鲜见报道.综合前人研究成果,我们知道: 在微弱分离比下, 系统中存在Blinking状态和repeated transients等弱非线性对流结构;在 ψ ⩽−0.1 的情形下, 系统沿弱非线性分支从Blinking到LTW对流过渡时对流振幅发生了很大的跳跃, 且存在迟滞现象[11], 但对于如−0.04⩽ψ⩽−0.01的极弱分离比, 人们并没有发现迟滞现象, 认为过渡是光滑的[5].造成这种差异的原因之一是从Blinking分支到LTW分支的跳动可能很小, 实验无法观察到[11]; 另一方面, 考虑到数值方法的低分辨率, 数值模拟中可能将其归结为一种数值误差[8].本文利用文献[19,24]中的高精度高分辨率数值方法求解流体力学基本方程组, 对极弱分离比 ψ =−0.02 下混合流体RB对流进行直接数值模拟, 给出解的分岔图, 并探讨系统沿着非线性分支出现的各种行波对流斑图及其演化过程.

2 数学物理模型与数值方法

2.1 数学物理模型

一薄层双流体混合物(如酒精与水)封入长为L, 高为 d的矩形腔体之中(长高比记为 Γ =L/d ),如图1所示, 考虑二维区域 ( x,z)∈[0,L]×[0,d]上的对流运动.整个装置处在一个均匀的重力场g=−gez中, g 为重力加速度, 方向竖直向下, ez为 z 方向的单位矢量.腔体底部和顶部表面都保持均匀一致的温度 T hot 和 T col , 且 T hot>Tcol.二者之间存在温度差 ∆ T=Thot−Tcol, 它由一个反映浮力与黏性对比的特征量Rayleigh数来表征.这里 α 是热膨胀系数, κ 是热扩散系数, ν是动力黏性系数.随着 Ra 的增大, 浮力效应突出,当 R a 增大到临界值时, 系统中便发生了对流运动.

图1 对流模型示意图Fig.1.Sketch of the two-dimensional convection model.

基于布辛涅斯克(Boussinesq)假设, 选取d,d2/κ, κ /d 和 ρ0κ2/d2分别为长度、时间、速度和压力的特征尺度进行无量纲化, 并定义无量纲温度和无量纲浓度如下

其中下标为0的量是对应物理量的平均值, 则描述该问题的无量纲控制方程组为[24]:

其中 v (u,w),θ,c,p,t 分别为速度场, 温度, 浓度, 压力和时间; 方程组中除 Rayleigh 数 Ra 外, 还包含了其他三个无量纲参数, 即

分别为 Prandtl数, Lewis 数和分离比.其中 ST为Soret系数, β 和D分别是体积膨胀系数和溶质扩散系数.对于纯流体, Rayleigh-Bénard 对流系统初始不稳定的临界Rayleigh数是=1708 , 为方便起见, 通常采用约化的Rayleigh数r=Ra/作为控制参数.在室温下, 酒精与水混合液的Lewis数 Le 在 0.01 附近, Prandtl数介于 5 —20 之间.分离比 ψ 是表征流体非线性特性的参数, 表示温度场对浓度场的Soret耦合效应.对于室温下酒精与水的混合物, ψ 值介于− 0.5—0.2 时, 实验容易控制[6].本文选取与文献[6−8]研究行波对流TW和LTW时相同的参数, 即 P r=10,Le=0.01 ,考虑极微弱负分离比 ψ =−0.02 , 对长高比Γ=12的腔体中对流情况进行直接数值模拟.

求解方程组时必须给出合理的边界条件.腔体所有壁面为无滑移, 浓度不可穿透; 上下壁面恒温,左右壁面绝热.相应边界条件的表达式如下:

其中 ∂ Ω 表示腔体的四个壁面.

2.2 数值方法和特征参数

本文采用具有四阶精度和高分辨率的紧致差分方法数值求解方程组(1)和(2), 该方法已成功用于我们以前的数值计算工作[19,24].该方法采用具有高分辨率的四阶组合紧致迎风格式[19,31]离散非线性对流项, 利用四阶对称Padé紧致格式离散黏性项; 压力泊松方程采用文献[32]中的四阶紧致差分方案进行离散; 使用多重网格方法进行计算以加快迭代过程.文献[24]对数值方法的正确性和网格无关性等进行了验证, 针对本文考虑的极弱分离比情况, 我们也进行了网格无关性检验, 发现241×21网格和精细网格 7 21×61 下特征参数(Nusselt数 N u 、流体混合参数 M和对流振幅|w|max)的最大相对误差均小于3%.在保证计算准确的前提下, 为了减少计算量, 本文采用241×21的网格进行数值模拟.相应地, 我们也进行了时间步长无关性检验, 在 r =1.13 时采用多个不同时间步长 ∆ t 进行了计算.比较发现 ∆ t=4×10−5步长和小时间步长 ∆ t=1×10−7下各个特征参数的最大相对误差均小于1%.于是选择前者作为本次模拟研究的时间步长.

为了定量表征对流特性, 对模拟结果进行深入分析, 我们使用如下特征参数:

1) 对流振幅: 以最大垂向速度| w |max来表征对流振幅的大小.

2) Nusselt 数: 用于度量系统的传热特性, 定义为[6]

其中尖括号 ⟨∗⟩ 表示在垂直方向上求和.Nu表示垂直通过流体层的热通量与热传导状态下的热通量 之 比, 在 热 传 导 状 态 下, 满 足 N u=1 , 所 以Nu−1反映了对流对传热的贡献.特别地，当取值上壁面(冷壁)时, (3)式简化为

3) 流体混合参数 (mixing parameter): 用于表征浓度场相对于传导状态时的变化, 其计算公式为[6]

其中尖括号和上划线分别表示侧向和垂向平均.

4)对于周期为 τ 的非定常流动, 考察上述特征参数在一个周期上的平均值 F ,

其中F代表 | w |max,Nu 和M.

3 结果与分析

3.1 解的分岔

对于本文讨论的负分离比 ψ =−0.02 , 对流发生的临界Rayleigh数为 rc=1.035 , 传导状态的首次不稳定为振荡不稳定, 这种不稳定导致了一个亚临界分岔[1,33](当 r rc时传导状态亦不稳定).图2给出了Nusselt数 N u−1 关于Rayleigh数的分岔曲线, 分岔图中包含了多个定常对流SOC, Blinking和LTW等基本解的分支, 其中点线表示不稳定解; 同时将纵坐标为0的线(点划线)标出, 以便于比较.需要指出的是: 我们重点关注的是稳定解的分支及状态之间的转换, 不稳定解的分支只是示意结果, 并非计算所得.

图2 ψ =−0.02 时 Nusselt数 N u−1 关 于 Rayleigh 数r的分岔曲线.SOC n 表示具有n个涡卷的SOC状态Fig.2.(a) Bifurcation diagram for ψ =−0.02.(b) close-up view of the part of the bifurcation diagram delimited by the straight dashed lines depicted in (a).Where SOC n represents the SOC solutions with n rolls.

从对流临界点 rc出发, 不稳定行波对流状态分支向上弯曲与稳定Blinking状态分支在点处连接.沿着稳定Blinking状态分支,逐渐增大Rayleigh数, Blinking状态的对流强度逐 渐 增 大, 在时过渡到稳定的LTW状态(详细讨论在3.2小节); 继续增大r, 随着对流强度的增大, Nusselt数和对流振幅均逐渐增大, 系统在 r∗=1.022 时过渡为含有12个对流涡卷、波数为3.14的临界SOC状态–SOC12.之后,若沿着 SOC12分支进一步增大 Rayleigh数, 则SOC12分支经数次分岔形成具有不同波数的SOC分支, 这与弱分离比 ψ =−0.1 下有相似的结论[24].若 从 r∗出 发 , 沿 着 SOC12分 支 逐 渐 减 小Rayleigh 数, 则 Nusselt数随之减小, 但系统仍处于稳定的SOC12状态, 直到SOC状态的鞍结分岔点在之下, 随着对流振幅的急剧减小, 系统最终回到了传导状态.

对流发生的临界Rayleigh数 rc、SOC状态的鞍结分岔点、LTW向SOC状态过渡点 r∗等由表1 给出.为便于比较, 我们给出了 ψ =−0.1 下的结果.与 ψ =−0.1 时不同, 极弱分离比− 0.02 下不存在稳定TW状态, 且在SOC12分支下端不存在UTW状态.然而, 此分离比下存在Blinking状态的解, 该状态稳定的 Rayleigh数范围为1.013⩽r<1.0172, 位于对流发生的临界Rayleigh数 rc=1.035 之下的亚临界区域.沿 Blinking状态分支向下, 流动从Blinking状态过渡到传导状态,沿此分支向上则成长为LTW状态.LTW状态位于 对 流 临 界 Rayleigh数 rc之 下, 其 存 在 的Rayleigh 数范围为 ∆ r=0.005 , 较 ψ =−0.10 时更窄.从 LTW状态向下到 Blinking状态和从Blinking向上到LTW状态过渡的Rayleigh数是不相等的, 表明二者之间过渡存在迟滞现象.

表1 分 离 比 ψ =−0.10 和 –0.02 时 , 各 状 态 临 界Rayleigh数的比较Table 1.Comparison of critical Rayleigh numbers for each state, ψ =−0.10 and− 0.02.

3.2 交替闪烁(Blinking)状态

这里讨论极弱分离比下的稳定Blinking状态,重点关注Blinking形成与稳定机理、闪动频率和行波频率随Rayleigh数变化以及向LTW状态过渡.

3.2.1 Blinking 状态的时空结构

在 r =1.015 时, 给温度场施加如下初始微小扰动[24]:

速度和浓度初始值均为零, 从 x0=2Γ/3 的非对称扰动开始, 流场经长时间的演化和发展形成了稳定的Blinking状态, 其流场的时空结构如图3(a)所示.可以看到, 随着时间的发展, 流场中左行波和右行波交替成长和衰落, 斑图结构由左行波和右行波交替控制.当左行波逐渐加强时, 缺陷(左行波和右行波的连接处)的中心向腔体右侧移动, 右行波的成长空间被压缩, 强度逐渐减弱; 当缺陷的中心靠近右端壁时, 在壁面的反射作用下发生转向后向左侧移动; 此后右行波逐渐成长而左行波逐渐减弱, 直到缺陷中心靠近左端壁; 然后在左端壁面的反射作用下再次发生转向, 进入下一个循环, 形成了左右闪动的对流现象.从观测点一 ( 0.13Γ,0.5) 处速度w随时间的变化(图3(b)) 可以看出, 伴随着对流振幅的逐渐增大流场经历了长时间(∆t≈3500) 的对传波 (counterpropagating waves, CPW)状态, 当对流达到一定强度便转变为稳定的Blinking状态.相应的功率谱密度 (图3(c)) 表明Blinking状态具有两个主频, 即闪动频率 ω1和行波频率 ω2.

图4给出了 4 000⩽t⩽4600 时段内两个观测点处垂向速度的变化情况.靠近左侧壁面的观测点一和靠近右侧壁面观测点二 ( 0.87Γ,0.5) 的速度w分别反映了左行波和右行波振幅的变化情况.可以看到, 左行波的振幅增加的时段对应于右行波的振幅衰减的时段, 反之亦然.在左行波振幅较大的时段, 流场主要由左行波控制, 腔体的左半边对流较强; 相应地, 在右行波振幅较大的时段, 流场主要由右行波控制, 腔体的右半边对流较强.

图3 r =1.015 时 (a) 流场时空结构, (b) w (0.13Γ,0.5) 的时间序列和 (c)功率谱密度Fig.3.(a) Spatio-temporal structure, (b) the time series of w (0.13Γ,0.5) and (c) power spectral density for r =1.015.

图4 r =1.015 时, 两个观测点处垂向速度 w 随时间的发展Fig.4.The time series of the vertical velocity w at two monitoring points (a) ( 0.13Γ,0.5) and (b) ( 0.87Γ,0.5) for r =1.015.

3.2.2 Blinking 状态的特性

前人研究结果表明, 两端壁面的反射作用在形成和保持Blinking状态中起到了一定的作用, 那么, 端壁反射作用是不是惟一原因, 有没有其他因素在起作用呢？我们知道, 前述 r =1.015 时Blinking状态是从不对称初始扰动出发计算得到的, 此时我们会想到不对称初始扰动或许是形成Blinking状态的原因之一.为了证实这一点, 在相同的Rayleigh数下, 我们从对称的扰动(即在(6)式中取 x0=Γ/2 )开始计算, 最终获得了与非对称扰动初值时相同的Blinking状态.但是当r=1.013时(为Blinking状态稳定区域下端点),结果却不同: 从非对称扰动出发可获得稳定的Blinking状态, 而对称初始扰动未能发展起来, 流场回到了传导状态.由此可见, 在Blinking状态存在的Rayleigh数范围的边界附近, 非对称的初始条件是诱发Blinking的原因之一, 在该范围的内部, 其诱发作用减弱.在获得稳定的Blinking状态后, 若改变两侧壁边界条件为周期性边界条件后继续计算, 流场便回到了传导状态, 这表明两端壁面反射作用是维持Blinking状态的重要原因[9,11].

图5给出了 r =1.0171 时Blinking和LTW两种状态中, 腔体水平中心线上瞬时垂向速度、温度和浓度波的波形, 以及浓度场斑图结构.可以看出,两种状态浓度波与速度波和温度波之间均存在相位差, 这是导致波动向前传播的直接原因.图5(a)和图5(c)为Blinking状态行波向左壁面传播的情形, 此时处于对流成长的中间阶段, 其对流强度依然明显较弱, 与LTW状态的强非线性对流结构(图5(b)和图5(d), 其浓度波为带尖角的台型结构, 且其振幅远大于Blinking状态的振幅)不同,是一种弱非线性结构.有趣的是, 由于此Rayleigh数位于Blinking状态稳定区域的上界附近, 行波成长起来后最强处可达到LTW状态的水平, 其浓度波形具有与LTW状态类似的台型结构, 但该状态极不稳定, 对流振幅迅速衰减, 之后又再次成长.

在Blinking状态稳定的Rayleigh数范围内,其闪动频率 ω1和行波频率 ω2随Rayleigh数的变化如图6所示.为了比较, 图6(b)中同时给出了LTW状态行波频率随Rayleigh数的变化.随着Rayleigh 数 的 增 大 , 闪 动 频 率 ω1先 减 小 , 在Ra=1734.9(r≈1.01575)处取得极小值, 在此之后逐渐增大并在 R a=1735.8(r≈1.01628) 达到极大值, 最后再次减小.行波频率 ω2随 Rayleigh 数的增大开始缓慢递减, 在 r =0.015 之后迅速减小, 行波传播速度逐渐减缓.当Ra=1737.3(r≈1.0172)时, 行波不再左右闪动, 伴随着行波频率的再次调整, 系统从 Blinking 过渡到 LTW.过渡时, 行波频率的变化是不连续的, 存在明显的跳跃, 突然提升到1.5倍以上的水平 (如图6(b)).

从Blinking到LTW状态过渡时, Nusselt数也有明显的跳跃, 如图7所示.相比于强非线性的LTW状态, 弱非线性Blinking状态的对流振幅相对较小, 对流强度较弱, 反映对流传热贡献的Nusselt数较小.过渡时, 伴随着对流振幅的急剧增长, 对流传热Nusselt数快速增大, 流动达到强非线性对流的水平并维持下来, 便形成了跳跃.通过大量的模拟计算和分析, 我们进一步发现从Blinking到LTW状态的过渡是不光滑的, 存在迟滞 现 象.沿 Blinking 分 支 增 大 Rayleigh 数, 从Blinking状态到LTW状态的过渡点为 r =1.0172 ,而沿LTW分支减小Rayleigh数直到 r =1.0171 ,系统才从LTW过渡到Blinking状态.可以看到,这种滞后非常微小, 仅为 ∆ Ra∼ 10−1, 以致在实验中未能观察到[5].可以预测, 在充分小的分离比下,Blinking和LTW状态之间过渡的迟滞现象就会消失, 从而实现光滑过渡.

图5 r =1.0171 时, Blinking 状态与 LTW 状态流场典型波形和浓度场的比较Fig.5.Comparison of the lateral profiles and concentration fields between the Blinking and LTW states at r =1.0171 : (a) The lateral profile and (c) concentration field of the Blinking state; (b) The lateral profile and (d) concentration field of the LTW state.

图7 Blinking 和 LTW 状态的 N u−1 随 r的变化情况.(b)为 (a)中虚线标注矩形区域的局部放大Fig.7.The variation of N u−1 of the Blinking and LTW states as a function of r.(b) Close-up view of the part delimited by the straight dashed lines depicted in (a).

3.3 LTW状态

本研究获得的LTW状态的流场结构如图8所示, 它是一种行波前缘与壁面接触(称为Wallattached, LTW)[24]的局部行波对流状态, 其对流区域宽度约为 9.0 , 不随Rayleigh数改变而变化.行波传播方向可能是向左也可能向右, 这依赖于流场的初始结构, 但两者Nusselt数 N u 、混合参数M、振幅 | w |max及行波频率等特征参数完全相等, 在分岔图图2中不做区分.事实上, 将左(右)行LTW状态的解沿 Γ /2 做对称投影 x ↔ Γ−x , 即可获得相应的右(左)行LTW状态.由于浓度场和速度场之间存在相位差, 从而导致了一个覆盖整个对流区域的大尺度时均浓度环流, 在前缘和壁面之间形成了一个与大尺度环流反向的二次环流(图8(b)和图8(c)).大尺度浓度环流引起的浓度再分配使得LTW状态稳定下来; 反向二次环流降低了行波对流的速度, 在环形实验装置中(计算中为周期性边界), 它降低了局部对流区域整体漂移的群速度.关于LTW状态的稳定机理等前人已做了大量研究, 我们的结论与前人并无明显不同, 本文不再赘述.

3.4 SOC状态

3.4.1 SOC12状态的形成

图8 LTW 状态的流场结构Fig.8.Structures of the flow field of LTW state: (a) Spatio-temporal structure; (b) a large-scale concentration current; (c) a transient structure of the concentration field.

图9 r =1.13 时流场时空发展和典型时刻的瞬时结构Fig.9.The spatio-temporal development and transient structures of the flow field at typical times for r =1.13.

从零场开始计算, 在 r =1.13 时我们获得了临界SOC状态–SOC12.图9(a)给出了腔体水平中心线上垂向速度的时空结构, 展示了对流从零场开始到发展为SOC12状态的整个过程, 几个关键时刻的流场瞬时结构如图9(b)所示.流场瞬时结构为流线和温度扰动(温度值 θ 减去热传导状态的解,表现了温度偏离传导状态的大小)云图的叠加图,其中实(点)线为正(负)值, 对应涡卷逆(顺)时针旋转, 粗虚线标出了腔体中心 x =Γ/2.从初始零场开始, 随着对流振幅的线性增长, 大约到了t=30时, 腔体中心位置的涡卷率先发展起来并逐渐向两侧扩散, 形成了微弱的CPW状态, 该状态的中心驻波位于腔体的中心; 伴随着对流振幅的持续增长, CPW 状态持续到 t ≈40 ; 之后, 壁面反射波与CPW的左行和右行波在壁面附近产生碰撞和抵消, 使得两侧壁面附近的对流强度降低, 而CPW中心驻波得到了很大的发展( t =46 ); 壁面反射波经过碰撞后进一步向腔体中心传播, 与驻波进行叠加并在中心处湮灭, 使得腔体内对流涡卷数量减少, 在 4 2

对流传热Nusselt数和流体混合参数M随时间的发展如图10(a)所示, 观测点一和观测点三(0.5Γ,0.5)处垂向速度的时间序列由图10(b)给出.在SOC12状态的形成过程中, Nusselt数曲线整体上呈现先振荡增长后振荡衰减之后再快速增长的变化形式, 期间的振荡由反射波与CPW波之间竞争(有抵消也有叠加)所致; 当二者这种竞争结束后, 腔体中对流涡卷快速发展起来, 系统在短时间内达到饱和状态SOC12, 对应于 t >50 后Nusselt数曲线的指数增长.相应地, 对流快速增强使得流体的混掺程度升高, M值迅速减小并最终达到一固定值.

因SOC12是偶宇称的, 该状态的解满足对称性[34]Rx:(u,w,θ,c)(x,z)=(−u,w,θ,c) ( Γ−x,z) ,观测点一和观测点二的垂向速度w相等, 故这里只给出观测点一的时间序列.观测点一靠近壁面(如图10(b)上), 该点处垂向速度的振幅先增长后衰减, 分别对应于CPW波的线性增长和壁面附近反射波与CPW波的碰撞, 振荡则反映了行波经过观测点时相位和振幅的变化.振荡结束时, 在腔体中心区域形成了由两对涡卷构成的驻波结构, 其他区域处于热传导水平, 如图9(b) 所示 t =50 时.之后对流涡卷只在原地滚动, 不再向外移动; 由于没有竞争, 驻波对流振幅快速增长, 这些原地滚动的涡卷逐渐带动其附近流体发生对流, 使其内侧(靠近腔体中心)涡卷得到增强从而达到饱和状态(涡卷宽度为1, 一个波长内包含一对涡卷), 在其外侧则产生旋转方向与之相反的新涡卷; 新涡卷逐渐成长又带动该新涡卷附近流体产生对流, 如此重复直至对流涡卷填满整个腔体.新涡卷从产生到达到最终饱和状态的成长过程中, 由于新涡卷逐渐变大,它与旧涡卷之间的挤压必然导致其自身向外侧移动, 即涡卷在壮大的同时伴随着移动.另外, 观测点一的垂向速度序列显示, 在振荡过后出现了一个单独波动, 这是因为: 大约在 t =65 时刻, 该观测点右侧出现了顺时针旋转的新涡卷, 此时观测点处于该涡卷的左侧边缘, 垂向速度w方向向上(w>0)但其值很小; 随着该涡卷的成长和向左移动,w在 t ≈73 时刻达到最大值; 涡卷进一步向左移动,其中心逐渐移向观测点, 于是w逐渐减小, 并在中心经过观测点时( t ≈80.6 时刻)为 w =0 , 期间该涡卷左侧又产生了逆时针旋转的涡卷; 原顺时针涡卷中心经过观测点后, 观测点处垂向速度方向向下 ( w <0 ), 在 t ≈84 时刻 w 达到极小值; 之后左侧新产生的逆时针涡卷逐渐成长并接触到左侧壁面, 固壁的约束作用使其在成长过程中将原来的顺时针涡卷向右挤压, 原顺时针涡卷改变方向, 向右移动且中心经过观测点后固定下来(观测点最终位于该顺时针涡卷的中心略偏左处), 对应于垂向速度w经过调谐后达到一固定值 w =0.166.观测点三位于几何中心, 该点处垂向速度的变化较为简单, 其振幅先振荡增长后衰减而后再增长, 最后达到一固定水平: 开始的增长体现了CPW中心驻波对流振幅的线性增长; 衰减是由反射波对中心驻波的干扰造成的, 对应于 4 2

图10 r =1.13 时 (a) N u−1 和 M 的变化及 (b) 观测点处垂向速度的时间序列Fig.10.The variation of (a) N u−1 , M, and (b) the vertical velocity at the monitoring points for r =1.13.

3.4.2 SOC12状态的特性

保持 r =1.13 , 到 t =358.56 时刻, 我们得到了稳定的 SOC12状态, 其流场结构如图11所示, 通过对比垂向速度、温度和浓度的波形发现, 三者是同相位的, 不存在相位差; 浓度场的波形呈现有尖角的特殊情形, 沿腔体的长度方向, 尖角呈现向上向下的交替分布.尖角位于相邻对流涡卷的交界面处, 尖角向上对应左侧顺时针(左顺)和右侧逆时针(右逆)的一对涡卷之间 (图11(b)), 左顺右逆的旋转形式将腔体上部的低温高浓度组分向下运输,使该点的浓度高于、温度低于其临近区域.从波形图也可以看出, 当温度场内出现波峰时, 浓度场内即出现波谷, 这一特点与温度高的地方浓度反而低的现象是符合的.从浓度场的分布来看(图11(c)),流体得到了充分混合, 浓度均匀分布.

模拟发现, 沿着SOC12状态分支逐渐增大Rayleigh数的过程中, 系统会经历数次分岔, 并产生多条具有不同涡卷数的SOC状态.这些SOC状态在波数-Rayleigh数空间 ( k,r) 上较大的范围内共存, 其稳定性及向混沌状态的过渡等另文讨论.不同分离比下, SOC12状态 Nusselt数关于Rayleigh数变化如图12所示(注意: 在同一坐标尺度下两条曲线靠得很近, ψ =−0.02 曲线略位于ψ=−0.1的上方.为便于区分, 将 ψ =−0.02 曲线向左平移了0.5个r单位, 其横轴坐标数值标在图上方).可以看到, Nusselt数随 Rayleigh 数指数增长, 当 r ∈[1.008,7.59]时, SOC12状态是稳定的; 在此范围之上, 对流发生二次不稳定并过渡为混沌状态, 不存在与之对应的稳定UTW状态.这与分离比 ψ =−0.1 下的结论不同: 此分离比下的SOC12状态在 1.064⩽r⩽5.69 内稳定, 在该区域的两端均存在稳定的UTW状态.

图11 r =1.13 时 SOC12 状态的流场结构Fig.11.The structure of flow field for the SOC 12 state at r=1.13: (a) The lateral profile on the horizontal centerline of the cavity; (b) the streamlines and the structure of the associated temperature field; (c) the structure of the concentration field.

图12 SOC 12 状态 Nusselt数随 Rayleigh 数的变化Fig.12.The variation of N u with r for the SOC 12 state.

4 结 论

本文利用高精度数值方法求解流体力学基本方程, 模拟研究了极弱分离比 ψ =−0.02 下, 长高比 Γ =12 的腔体内混合流体RB对流解的分岔,Blinking和SOC12等对流结构的形成与特性, 以及各种对流状态之间的过渡.研究发现Blinking状态仅在较小的相对Rayleigh数范围∆r=0.004内稳定; 从Blinking状态向LTW状态的过渡存在迟滞现象, 但这种滞后非常微小, 仅为 ∆ Ra∼10−1,可能实验中未能发现这种现象, 但Kolodner和Surko[5]预测了极微弱分离比下可能存在迟滞现象.本结果从数值计算的角度说明了预测是正确的, 迟滞现象是存在的, 当然这一结论还需实验和理论的进一步验证.系统在 Blinking状态和LTW状态之间过渡时, 行波频率、对流振幅以及Nusselt数等特征参数的变化是不连续的, 有明显的跳跃.在Blinking状态存在的Rayleigh数区间下界附近, 外部施加的不对称初始扰动是形成该状态的诱因, 而在其他区域, 系统的振荡不稳定足以导致Blinking状态的形成.沿对流分岔的定常对流SOC分支, 随着Rayleigh数的增大, 临界SOC12状态经历数次分岔, 形成多条具有不同波数的SOC状态的分支, 最终过渡到混沌状态.与分离比ψ=−0.1时[24]不同, SOC12分支的两端均不存在稳定UTW状态.