EI、Scopus 收录
中文核心期刊
快速检索 年期检索 高级检索

高径比对GaAs熔体液桥热毛细对流失稳的影响

周游, 曾忠, 刘浩, 张良奇

周游, 曾忠, 刘浩, 张良奇. 高径比对GaAs熔体液桥热毛细对流失稳的影响. 力学学报, 2022, 54(2): 301-315. DOI: 10.6052/0459-1879-21-227
引用本文: 周游, 曾忠, 刘浩, 张良奇. 高径比对GaAs熔体液桥热毛细对流失稳的影响. 力学学报, 2022, 54(2): 301-315. DOI: 10.6052/0459-1879-21-227
shu
Zhou You, Zeng Zhong, Liu Hao, Zhang Liangqi. Effect of aspect ratio on thermocapillary convection instability of GaAs melt liquid bridge. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(2): 301-315. DOI: 10.6052/0459-1879-21-227
Citation: Zhou You, Zeng Zhong, Liu Hao, Zhang Liangqi. Effect of aspect ratio on thermocapillary convection instability of GaAs melt liquid bridge. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(2): 301-315. DOI: 10.6052/0459-1879-21-227
shu
周游, 曾忠, 刘浩, 张良奇. 高径比对GaAs熔体液桥热毛细对流失稳的影响. 力学学报, 2022, 54(2): 301-315. CSTR: 32045.14.0459-1879-21-227
引用本文: 周游, 曾忠, 刘浩, 张良奇. 高径比对GaAs熔体液桥热毛细对流失稳的影响. 力学学报, 2022, 54(2): 301-315. CSTR: 32045.14.0459-1879-21-227
shu
Zhou You, Zeng Zhong, Liu Hao, Zhang Liangqi. Effect of aspect ratio on thermocapillary convection instability of GaAs melt liquid bridge. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(2): 301-315. CSTR: 32045.14.0459-1879-21-227
Citation: Zhou You, Zeng Zhong, Liu Hao, Zhang Liangqi. Effect of aspect ratio on thermocapillary convection instability of GaAs melt liquid bridge. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(2): 301-315. CSTR: 32045.14.0459-1879-21-227
shu

高径比对GaAs熔体液桥热毛细对流失稳的影响

  • *.

    重庆大学航空航天学院, 重庆 400044

  • †.

    重庆交通大学西南水运工程科学研究所, 重庆 400016

基金项目: 国家自然科学基金(12172070, 12102071)、中央高校基本科研业务费(2021 CDJQY-055)和重庆市教委科学技术研究(KJQN202100706)资助项目
详细信息
    作者简介:

    曾忠, 教授, 主要研究方向: 流体力学、热毛细对流、多相流. E-mail: zzeng@cqu.edu.cn

  • 中图分类号: O35

EFFECT OF ASPECT RATIO ON THERMOCAPILLARY CONVECTION INSTABILITY OF GaAs MELT LIQUID BRIDGE

  • *.

    College of Aerospace Engineering, Chongqing University, Chongqing 400044, China

  • †.

    Chongqing Southwest Research Institute for Water Transport Engineering, Chongqing Jiaotong University, Chongqing 400016, China

HTML全文

    摘要
  • 摘要: 采用基于谱元法线性稳定性分析方法, 研究了高径比对GaAs熔体(Pr = 0.068)液桥热毛细对流失稳的影响, 同时结合能量分析揭示了热毛细对流的失稳机制. 研究结果表明: 与典型低普朗特数(例如Pr = 0.011)熔体静态失稳模式和典型高普朗特数(例如Pr>1)熔体振荡失稳模式不同, GaAs熔体热毛细对流失稳模式依赖于液桥高径比(As). 随高径比的变化, GaAs熔体热毛细对流存在两种失稳模式. 高径比As 在0.4≤As≤1.18范围内, 热毛细对流失稳是从二维轴对称定常对流转变为三维周期性振荡对流(振荡失稳); 高径比在1.20≤As≤2.5范围内, 热毛细对流失稳是从二维轴对称定常流动转变为三维定常流动(静态失稳). 典型的高普朗特数熔体液桥热毛细对流失稳机制是热毛细机制; 典型的低普朗特数液桥热毛细对流失稳机制是水动力学惯性机制. 本文基于扰动能量分析的结果表明: GaAs熔体热毛细对流失稳同时包括水动力学惯性失稳机制和热毛细失稳机制的贡献, 其中水动力学惯性失稳机制占主导作用, 两种机制对热毛细对流失稳能量贡献的占比随高径比的变化而变化.
    Abstract: In this paper, we explore the effect of aspect ratio on the instability of thermocapillary convection in GaAs melt (Pr = 0.068) liquid bridge by using the linear stability analysis in the context of spectral element method. Besides, we provide physical insight on the underlying instability mechanism via energy analysis. Differing from the cases of typical low Prandtl number (such as Pr = 0.011) and typical high Prandtl number (such as Pr > 1), which correspond to stationary instability and oscillatory instability respectively, the instability of the thermocapillary convection of GaAs melt (Pr = 0.068) is of note due to its noticeable dependence on the aspect ratio (As). In particular, we observe two instability modes for the flow considered here with the variation of the aspect ratio. When the aspect ratio As ranges from 0.4 to 1.18, thermocapillary flow transits from two-dimensional axisymmetric steady convection to three-dimensional periodic oscillatory convection (oscillatory instability). While for 1.20 ≤ As ≤ 2.5, the stationary instability appears and the two-dimensional axisymmetric steady flow transits to three-dimensional steady flow. As for the instability mechanism of the thermocapillary convection, the liquid bridge of high Prandtl number is characterized by thermocapillary mechanism, while the case of low Prandtl number features the hydrodynamic inertia mechanism. Based on disturbance energy analysis, it is shown that the instability of the present thermocapillary convection arises from the combined action of the hydrodynamic inertial instability and thermocapillary instability, in which the hydrodynamic inertial instability mechanism is dominant, and the specific proportion of these two contributions varies with the aspect ratio.
    • HTML全文

  • 引言

    热毛细对流是沿流体−流体界面非均匀温度导致的表面张力梯度驱动的对流, 其广泛存在于晶体生长[1-2]、液滴迁移[3-4]、3D打印[5-6]等工业应用中, 其中晶体生长是其重要的研究背景之一. 浮区法是一种无坩埚高品质单晶生长技术, 随着航空航天技术的发展, 空间微重力环境为浮区法生长更大尺寸和更高品质的单晶提供了可能. 在空间微重力环境中, 热浮力对流极度衰减, 热毛细对流成为浮区法晶体生长中熔体物质与热输运的主要方式. 由两块不同温度同轴圆盘以及圆盘之间的液柱组成液桥模型, 其可当成源于浮区法晶体生长技术的半浮区简化模型. 近40年来, 液桥模型已被广泛运用于热毛细对流研究.

    熔体普朗特数对液桥热毛细对流失稳机制和失稳模式有重要影响. Smith和Davis[7]基于线性稳定性分析研究了不同普朗特数无限大平面液层中热毛细对流的不稳定性, 发现了2种对流失稳类型. Levenstam等[8-9]通过三维数值模拟和实验对低普朗特数熔体液桥热毛细对流研究时发现, 热毛细对流第一次失稳是从二维轴对称定常流动转变为三维定常对流. Chen等[10]采用微扰动法研究了低普朗特数液桥热毛细对流失稳, 研究表明: 当Pr≤0.1时, 液桥热毛细对流的第一次失稳为静态失稳. Ichiro等[11]通过实验对高普朗特数熔体液桥热毛细对流进行了研究, 研究表明随着温差的增大, 热毛细对流第一次失稳是从二维轴对称定常流动转变为三维振荡对流. Wanschura等[12]通过线性稳定性分析研究了不同普朗特数对液桥热毛细对流的影响, 结果表明热毛细对流失稳存在两种失稳类型:当Pr≤0.05时, 液桥热毛细流动失稳从二维轴对称定常对流转变为三维定常对流, 当0.5≤Pr≤4.8时, 液桥热毛细对流失稳从二维轴对称定常对流转变为三维振荡对流, 两种失稳模式具有不同失稳机制. Chen等[13]采用与Wanschura等[12]同样的液桥模型, 通过线性稳定性分析研究了10−10≤Pr≤8.0 的热毛细对流失稳, 结果表明: 当Pr≤0.06时, 热毛细对流失稳是由水动力学惯性机制引起, 并且其失稳模式为静态失稳; 当Pr≥0.1时, 失稳是由热毛细机制引起, 并且其失稳模式为振荡失稳. 2001年, Levenstam等[14]结合数值模拟和线性稳定性分析研究了0.001≤Pr≤7范围内液桥热毛细流的失稳机制, 他们的研究表明: 当0.001≤Pr≤0.05时, 热毛细对流失稳是由水动力学惯性机制引起, 其失稳模式为静态失稳, 当0.057≤Pr≤0.068, 热毛细对流失稳主要是由水动力学惯性机制引起, 但其失稳模式转变为振荡失稳. 当0.07≤Pr≤1.80时, 热毛细对流失稳是由水动力学惯性机制和热毛细机制共同作用引起, 其失稳模式为振荡失稳, 当0.183≤Pr≤7时, 热毛细对流失稳是主要是由热毛细机制引起, 其失稳模式为振荡失稳. 上述文献结果表明: 液体普朗特数对液桥热毛细对流的失稳机制和失稳模式具有重要影响. 对高普朗特数液体, 热毛细对流第一次失稳是从二维轴对称定常对流转变为三维周期性振荡对流, 其失稳机制是热毛细失稳机制; 对低普朗特数液体, 热毛细对流第一次失稳是从二维轴对称定常对流转变为三维定常对流, 失稳机制为水动力学惯性失稳机制.

    针对锡熔体(Pr = 0.009)的热毛细对流, Li等[15]数值模拟研究了液桥高径比对热毛细对流失稳的影响, 结果显示: 当高径比在0.6~2.2之间时, 热毛细对流失稳临界Marangoni数(Mac)随着高径比的减小而单调增加, 其热毛细对流失稳模式均为静态失稳. Rybicki和Floryan[16]数值模拟研究了液桥的热毛细对流, 结果表明高径比对液桥热毛细流动有重要影响. Velten等[17]通过地面实验研究了高径比对液桥热毛细对流的影响, 结果表明随着高径比的增加Mac逐渐减小. Chen和Hu[18]通过线性稳定性分析研究了液桥的高径比对半浮区液桥热毛细对流的影响, 随着高径比的增加Mac数和临界频率随之减小. Nishino等[19]结合空间实验和线性稳定性分析方法, 研究了高径比对高普朗特数液桥热毛细对流不稳定性的影响, 其研究结果表明: 当Pr = 67时, 在高径比As = 1.25左右, 临界频率和Mac出现突降现象, 其失稳模式为振荡失稳. 吴勇强等[20]通过地面实验研究了液桥高径比对液桥起振温差的影响, 研究表明随着高径比的增大液桥内的起振温差逐渐呈减小趋势. 王佳等[21]对大尺寸液桥的浮力−热毛细对流进行了地面实验, 实验发现在高普朗特数情况下, 在Mac附近, 流场内会有行波现象出现, 随着高径比的变化,其流动模式发生改变. Kang等[22]在天宫二号上进行了高普朗特数液体桥热毛细对流的实验研究,以液桥高径比和体积比为特征, 研究了几何参数对液桥热毛细对流失稳临界条件的影响, 首次完整地得到了微重力条件下高径比−体积比参数空间内热毛细对流不稳定性的临界条件和振荡特性. Wang等[23]基于地面实验, 研究了存在重力场作用时, 几何参数(高径比和体积比)对大普朗特数大尺寸液桥热毛细对流失稳临界条件的影响, 在不同的几何参数下出现了6种不同的失稳状态. Liu等[24]通过线性稳定性分析研究了高径比对环形液池热毛细对流(Pr = 0.011和Pr = 1.4)的影响, 研究发现: 当高径比较小时, 失稳模式为振荡失稳; 当高径比较大时, 失稳模式为静态失稳, 这表示高径比对环形液池热毛细对流的失稳模式转变有影响, 但是该现象尚未在液桥热毛细对流研究的文献中报道.

    目前, 虽然文献已有很多高径比对液桥热毛细流失稳影响的研究, 但研究主要聚焦于高普朗特数熔体和低普朗特数熔体的热毛细对流失稳, 而较少关注普朗特数在0.057≤Pr≤0.8[14]范围液桥的热毛细流失稳, 特别缺乏这一范围内高径比对其热毛细失稳影响的研究. 本文针对GaAs熔体(Pr = 0.068), 采用谱元法对基态解进行求解, 然后基于线性稳定性分析确定失稳临界参数Mac, 最后结合能量分析方法研究对流失稳的物理机制, 研究高径比对液桥热毛细对流失稳临界参数、失稳模式和失稳机制的影响.

    1.   数学模型和数值方法

    1.1   物理模型与基本控制方程

    采用图1所示的液桥模型, 液桥由中心轴线重合的圆盘和圆形液柱组成, 圆柱形液桥高度为H, 半径为R, 高径比定义为As = H/R. 液桥上端为低温固壁(φc), 下端为高温固壁(φh), 熔体当成不可压缩牛顿流体. 假设表面张力为温度的线性函数, 圆柱自由表面绝热并且忽略自由表面变形.

    图  1  半浮区液桥模型
    Figure  1.  Floating half-zone model
    下载: 全尺寸图片 幻灯片

    分别用H, $\upsilon /H$, ${H^2}/\upsilon $和$\rho {\upsilon ^2}/{H^2}$作为长度、速度、时间和压力的特征物理量, 对流体的控制方程进行无量纲化, 其中$\upsilon $是运动学黏性系数, $\;\rho $是密度, 无量纲温度定义为$T{\text{ = (}}\varphi - {\varphi _c})/({\varphi _h} - {\varphi _c})$, 其中${\varphi _c}$为液桥低温圆盘固壁温度, ${\varphi _h}$为高温圆盘固壁高温, $\varphi $表示表示液桥内任意一点温度. 热毛细对流失稳前是二维轴对称定常流, 失稳后变为三维对流. 采用柱坐标($r,\theta ,z$), r为径向坐标, θ为周向坐标, z为轴向坐标. 无量纲控制方程为

    $$ \left. \begin{array}{l} \nabla \cdot {\boldsymbol{u}} = 0 \hfill \\ \dfrac{{\partial {\boldsymbol{u}}}}{{\partial t}} + \left( {{\boldsymbol{u}} \cdot \nabla } \right){\boldsymbol{u}} = - \nabla p + {\nabla ^2}{\boldsymbol{u}} \hfill \\ \dfrac{{\partial T}}{{\partial t}} + \left( {{\boldsymbol{u}} \cdot \nabla } \right){ T} = \dfrac{1}{{Pr}}{\nabla ^2}T \end{array} \right\} $$ (1)

    无量纲边界条件可以表示为:

    (1) 自由表面(r = 1/As)

    $$ u = 0,\;\;\frac{{\partial w}}{{\partial r}} = - \frac{{Ma}}{{Pr}}\frac{{\partial T}}{{\partial z}},\;\;\frac{{\partial T}}{{\partial r}} = 0 $$ (2)

    (2) 热壁(z = 0)

    $$ u = 0,\;\;w = 0,\;\;T = 1 $$ (3)

    (3) 冷壁(z = 1)

    $$ u = 0,\;\;w = 0,\;\;T = 0 $$ (4)

    其中u表示径向速度, w表示轴向速度, 无量纲参数Pr以及Ma分别定义为

    $$ Pr = \frac{\upsilon }{\kappa } $$ (5)
    $$ Ma = \frac{{{\gamma _T}\Delta \varphi H}}{{\rho \upsilon \kappa }} $$ (6)

    其中$\kappa $为热扩散系数, ${\gamma _T}$为流体表面张力系数, $\Delta \varphi = $$ {\varphi _h} - {\varphi _c}$.

    1.2   线性稳定性分析

    为了得到高精度的失稳临界条件, 采用基于时间分裂法的谱元法求解方程式(1)~式(4), 得到定常轴对称的基态解, 然后在所得二维基态解上施加一个三维小扰动[25]

    $$ ({\boldsymbol{u}},p,T) = ({{\boldsymbol{u}}_0},{p_{0,}}{T_0}) + ({\boldsymbol{\hat u}},\hat p,\hat T) $$ (7)

    其中${{\boldsymbol{u}}_0}$, ${p_0}$和${T_0}$表示速度、压力和温度的基态解, ${\boldsymbol{\hat u}}$, $\hat p$和$\hat T$表示扰动速度、扰动压力和扰动温度. 将式(7)代入到方程式(1)~式(4)中, 并忽略高阶扰动项可导出扰动控制方程

    $$ \left. \begin{array}{l} \nabla \cdot {\boldsymbol{\hat u}} = 0 \hfill \\ \dfrac{{\partial {\boldsymbol{\hat u}}}}{{\partial t}} + \left( {{{\boldsymbol{u}}_0} \cdot \nabla } \right){\boldsymbol{\hat u}} + \left( {{\boldsymbol{\hat u}} \cdot \nabla } \right){{\boldsymbol{u}}_0} = - \nabla \hat p + {\nabla ^2}{\boldsymbol{\hat u}} \hfill \\ \dfrac{{\partial \hat T}}{{\partial t}} + \left( {{{\boldsymbol{u}}_0} \cdot \nabla } \right)\hat { T} = \dfrac{1}{{Pr}}{\nabla ^2}\hat T \hfill \\ \end{array} \right\} $$ (8)

    扰动边界条件为:

    (1) 自由表面(r = 1/As)

    $$ \left.\begin{split} & \hat u = 0\\ &\frac{{\partial \hat v}}{{\partial r}} - \frac{{\hat v}}{r} + \frac{{Ma}}{{Pr}}\frac{1}{r}\frac{{\partial \hat T}}{{\partial \theta }} = 0\\ &\frac{{\partial \hat w}}{{\partial r}} + \frac{{Ma}}{{Pr}}\frac{{\partial \hat T}}{{\partial z}} = 0\\ &\frac{{\partial \hat T}}{{\partial r}} = 0 \end{split}\right\}$$ (9)

    (2) 热壁(z = 0)

    $$ \hat u = \hat v = \hat w = 0,\;\;\hat T = 0 $$ (10)

    (3) 冷壁(z = 1)

    $$ \hat u = \hat v = \hat w = 0,\;\;\hat T = 0 $$ (11)

    将扰动量表示为正则模形式

    $$ \hat \varphi = \tilde \varphi (r,z){{\rm{e}}^{\lambda t + {\rm{i}}k\theta }} + c.c. $$ (12)

    其中φ表示u, v, w, T和p, c.c.为复数虚部项的共轭. 复数λ可分为实部λr和虚部λi, 其中实部λr表示为线性增长率, 虚部λi表示为振荡频率ƒ, k表示为失稳波数. 其中λr<0表示流动稳定, λr>0表示流动不稳定, 而λr = 0则代表流动失稳的临界状态, 线性稳定性分析正是要捕捉λr = 0时的临界参数, 在临界条件下, λi = 0表示失稳是静态失稳, λi≠0表示失稳是振荡失稳. 在将表达式(12)代入式(8)~式(11), 采用谱元法对扰动方程进行离散, 得到一个广义特征值问题, 即

    $$ {\boldsymbol{A\tilde X}} = \lambda {\boldsymbol{B\tilde X}} $$ (13)

    其中, A, B为系数矩阵, ${\boldsymbol{\tilde X}}$表示扰动场变量的特征向量, 利用ARPACK程序包[26-27]求解式(13)对应的特征值问题. 基于隐式重启Arnoldi方法[28]的ARPACK程序包是一个开源代码库, 可以高效地求解大型非对称稀疏矩阵的特征值. 通过式(13)确定失稳临界波数kc后, 参考Boppana和Gajjar[29]的处理过程, 可以类似地确定失稳临界参数

    $$ M{a_c} = M{a_{j - 1}} + \frac{{(M{a_j} - M{a_{j - 1}}){\lambda _{r,j - 1}}}}{{{\lambda _{r,j - 1}} - {\lambda _{r,j}}}} $$ (14)
    $$ {f_c} = {\lambda _{i,j - 1}} + \frac{{M{a_c} - M{a_{j - 1}}}}{{M{a_j} - M{a_{j - 1}}}}({\lambda _{i,j}} - {\lambda _{i,j - 1}}) $$ (15)

    1.3   扰动能量分析

    确定流动失稳临界参数后, 可通过扰动能量分析的方法来研究流动失稳机制. 对扰动方程(8)中的动量方程乘以$\hat u$并对求整个求解域进行积分, 可得到扰动动能变化率表达式[12,30-31]

    $$ \frac{1}{{{D_k}}}\frac{{\partial {E_{{\rm{kin}}}}}}{{\partial t}} = \frac{1}{{{D_k}}}\int_V {\frac{{\partial {e_{{\rm{kin}}}}}}{{\partial t}}{\rm{d}}V = } {M_z} + {M_\varphi } + {I_V} - 1 $$ (16)

    其中$\partial {e_{{\rm{kin}}}}/\partial {t}$表示扰动动能的分布, Dk表示黏性耗散率, Mz和Mφ分别表示由扰动引起的热毛细力沿液桥轴向和周向做功的功率, IV项表示从基态流场到扰动流场的动能传递率. 则各个能量项的具体表达式为

    $$ {D_k} = \int_V {{{(\nabla \times {\boldsymbol{\hat u}})}^2}{\rm{d}}V - 2\int_S {\frac{{{{\hat v}^2}}}{r}{\rm{d}}S{n_r}} } $$ (17)
    $$ {M_z} = \frac{1}{{{D_k}}}\int_S {\hat w\frac{{\partial \hat w}}{{\partial r}}{\rm{d}}S} $$ (18)
    $$ {M_\varphi } = \frac{1}{{{D_k}}}\int_S {\left(\hat v\frac{{\partial \hat v}}{{\partial r}} - \frac{{{{\hat v}^2}}}{r}\right){\rm{d}}S} $$ (19)
    $$ \begin{split} & {I_V} = - \frac{1}{{{D_k}}}\int_V {{\boldsymbol{\hat u}}({\boldsymbol{\hat u}} \cdot \nabla ){{\boldsymbol{u}}_0}} {\rm{d}}V =\\ &\qquad\frac{1}{{{D_k}}}\int_V \left( - \hat u\hat u\frac{{\partial {u_0}}}{{\partial r}} - \hat u\hat w\frac{{\partial {u_0}}}{{\partial z}} - \frac{{\hat v\hat v{u_0}}}{r} -\right.\\& \qquad \left. \hat u\hat w\frac{{\partial {w_0}}}{{\partial r}} - \hat w\hat w\frac{{\partial {w_0}}}{{\partial z}} \right) {\rm{d}}V = \\ &{I_{V1}} + {I_{V2}} + {I_{V3}} + {I_{V4}} + {I_{V5}} = \sum\limits_{i = 1}^5 {{I_{Vi}}} \end{split} $$ (20)

    其中

    $$ {I_{V1}} = - \frac{1}{{{D_k}}}\int_V {\hat u\hat u\frac{{\partial {u_0}}}{{\partial r}}{\rm{d}}V} $$ (21)
    $$ {I_{V2}} = - \frac{1}{{{D_k}}}\int_V {\hat u\hat w\frac{{\partial {u_0}}}{{\partial z}}{\rm{d}}V} $$ (22)
    $$ {I_{V3}} = - \frac{1}{{{D_k}}}\int_V {\frac{{\hat v\hat v{u_0}}}{r}} {\rm{d}}V $$ (23)
    $$ {I_{V4}} = - \frac{1}{{{D_k}}}\int_V {\hat u\hat w\frac{{\partial {w_0}}}{{\partial r}}} {\rm{d}}V $$ (24)
    $$ {I_{V5}} = - \frac{1}{{{D_k}}}\int_V {\hat w\hat w\frac{{\partial {w_0}}}{{\partial z}}{\rm{d}}V} $$ (25)

    其中IV1, IV2和IV3表示基态径向速度向扰动场的能量传递, 而IV4和IV5则表示基态轴向速度向扰动场的能量传递. V表示液桥的体积, S表示自由表面.

    2.   结果与讨论

    2.1   网格依赖性和程序有效性验证

    选取Pr = 0.068, 液桥高径比As = 1,计算网格采用8r × 9z单元划分, 单元内分别采用5r × 5z, 6r × 6z, 7r × 7z阶谱离散, 对应网格节点为41 × 46, 49 × 55, 57 × 64. 将不同网格的计算结果进行对比, 如表1所示, 不同网格计算的波数吻合, Mac相对误差较小. 综合考虑计算精度和计算效率, 本文以下计算采用6r × 6z阶谱离散.

    表  1  基于不同网格分辨率得到的热毛细对流(Pr = 0.068)失稳临界Marangoni数(Mac)和波数kc
    Table  1.  The critical Marangoni number (Mac) and wave number kc for the instability of thermocapillary flow with Pr = 0.068 under different mesh resolution conditions
    Mesh resolution41 × 4649 × 5557 × 64
    Mac / kc1179.7/31182.4/31182.7/3
    下载: 导出CSV 
    | 显示表格

    选取液桥高径比As = 1的热毛细对流进行程序有效性验证, 计算网格采用8r × 9z单元划分, 6r × 6z阶谱离散, 对应网格节点为49 × 55. 将计算Mac转化为Rec ($R{e_c} = {{M{a_c}} \mathord{\left/ {\vphantom {{M{a_c}} {Pr}}} \right. } {Pr}}$)与Levestam等[14]的结果进行对比, 如表2所示.

    针对Pr = 0.009熔体热毛细对流的失稳, 基于我们自主开发谱元法程序计算结果与Li等[15] 的临界参数(失稳临界Marangoni数(Mac)和波数kc)进行对比, 如表3所示.

    以上与Levestam等[14]和Li等[15]的结果对比表明: 本文计算热毛细对流失稳临界参数与他们结果吻合, 我们程序的有效性得到了验证.

    表  2  计算失稳临界参数(Rec, fc和kc)与文献[14]结果对比
    Table  2.  The critical values (Rec, fc, kc) from the present computations against Ref. [14]
    PrRec/ fc/kc (present)Rec/ fc/kc [14]Relative error of Rec/fc
    0.06817379/205.3/317229/206/30.87% / 0.36%
    0.18314021/409.1/314276/415/31.78%/1.42%
    1.02544/64.9/22551/65/20.27%/0.15%
    下载: 导出CSV 
    | 显示表格
    表  3  计算失稳临界参数(Mac和kc)与文献[15]结果对比
    Table  3.  The critical values (Mac and kc) from the present computations against Ref. [15]
    AsMac/kc (present)Mac/kc[15] (LSA)Mac/kc[15] (DNS)
    0.632.27/332.21/334.73/3
    0.823.12/223.18/224.71/2
    1.016.76/217.00/219.77/2
    下载: 导出CSV 
    | 显示表格

    2.2   稳定性曲线与失稳模式

    GaAs熔体液桥热毛细对流的基态是二维轴对称定常对流. 固定液桥的高度为1, 液桥半径的变化范围为0.4≤R≤2.5, 则对应液桥的高径比变化范围为2.5≥As≥0.4. 离散采用6r × 6z阶谱离散, 计算网格采用(8r-14r) × 9z单元划分, 所以液桥轴向有54个网格, 径向根据高径比的不同有48~84个网格.

    如表4所示, 失稳波数随着高径比的增加而减小, 失稳波数的变化符合Zeng等[32]提出的高径比和失稳波数之间的关系1.57≤kcAs≤3.14. 如图2所示, 临界Mac数随着高径比增加, 整体呈下降趋势, Mac数在高径比为0.67和高径比为1.18处出现峰值点. 表4中的线性稳定性分析结果表明, 在高径比0.40≤As≤2.50区间内, 存在两种失稳模式, 当0.4≤As≤1.18时, 热毛细对流从二维轴对称定常对流转捩为三维振荡对流; 当1.20≤As≤2.5时, 热毛细对流从二维轴对称定常对流转捩为三维定常对流.

    表  4  Pr = 0.068时不同高径比下的失稳临界参数
    Table  4.  The critical values at different aspect ratios for Pr = 0.068
    AsMacfckcType
    0.402033.1±490.76



    0.431920.0±493.65
    0.481936.3±507.45
    0.531770.8±504.74
    0.591838.4±536.14
    0.671988.8±233.04type I (oscillatory instability)
    0.711278.0±200.34
    0.771083.9±201.74
    0.91995.3±155.23
    1.001182.4±205.33
    1.051246.8±63.62

    1.111273.6±62.72
    1.181393.6±36.72
    1.201348.40.02
    1.251058.60.02type II (stationary instability)
    1.33691.30.02
    1.43531.90.02
    1.54483.60.02
    1.67350.50.01
    1.82262.80.01
    2.00234.60.01
    2.22248.40.01
    2.50326.00.01
    下载: 导出CSV 
    | 显示表格
    图  2  Mac随高径比(As)的变化(Pr = 0.068)
    Figure  2.  The variation of Mac versus the aspect ratio As for Pr = 0.068
    下载: 全尺寸图片 幻灯片

    在图3中 IV项表示IV1项到IV5项的总和, M项表示Mz和Mφ的总和. 由图3可知, 失稳主要是由水动力学惯性机制主导. 结合表5扰动能量平衡分析的结果可知, IV4和IV5是造成流动失稳的主要原因, 且Mz和Mφ也占据一定比例, 表明失稳是由水动力学惯性机制占主导, 同时热毛细机制对失稳有所贡献. 在高径比As = 1.82处, 由水动力学惯性机制主导的IV4 + IV5的占比最大, 而热毛细机制引起的Mz + Mφ的占比最小, Levestam等[14]针对高径比As = 1熔体液桥热毛细对流失稳研究表明: 当Pr = 0.068时, 液桥热毛细对流是振荡失稳, 其失稳机制是水动力学的惯性失稳机制, 这与本文的分析结果一致. 下面具体介绍随高径比变化的两种失稳模式.

    图  3  扰动能量随高径比As变化
    Figure  3.  The variation of the kinetic energy budgets with respect to the aspect ratio
    下载: 全尺寸图片 幻灯片
    图  3  扰动能量随高径比As变化(续)
    Figure  3.  The variation of the kinetic energy budgets with respect to the aspect ratio (continued)
    下载: 全尺寸图片 幻灯片
    表  5  Pr = 0.068时不同高径比能量分析结果
    Table  5.  The energy analysis results at different aspect ratios for Pr = 0.068
    AsIV1/%IV2/%IV3/%IV4/%IV5/%(IV4 + IV5)/%Mz/%Mφ/%(Mz + Mφ)/%
    0.402.45−7.90−1.4451.6947.9299.613.543.677.21
    0.432.17−12.31−1.4559.6344.46104.094.043.347.38
    0.484.39−6.12−1.8148.9848.3897.363.112.996.10
    0.534.14−11.26−1.8458.6944.07102.763.552.536.08
    0.596.69−3.11−2.5344.1949.8194.002.592.224.81
    0.670.639.55−1.1118.1963.8782.062.576.178.74
    0.71−0.1810.44−2.6720.4062.9583.352.496.498.98
    0.770.8512.58−4.1619.4163.6583.061.865.757.61
    0.910.339.57−4.2922.3959.0681.453.366.509.86
    1.007.1016.23−5.7514.4361.5375.961.485.046.52
    1.05−1.95−23.35−1.0571.5926.6598.2418.958.4027.35
    1.112.49−16.152.0065.2123.2788.4815.517.0622.57
    1.189.54−6.576.3455.9019.2275.1210.954.3615.31
    1.2010.29−2.227.5350.6020.4171.019.733.6613.39
    1.254.911.527.6644.2827.3371.619.144.9314.07
    1.33−2.441.758.7843.1229.8672.9810.728.1818.90
    1.43−4.711.349.2243.8028.8172.6110.8810.6421.52
    1.54−4.842.519.1843.2727.7170.9810.0412.1422.18
    1.675.33−23.21−13.15117.748.12125.8616.78−11.675.11
    1.825.42−19.76−13.34118.697.33126.0215.32−13.651.67
    2.004.22−14.91−12.64108.4611.71120.1716.04−12.883.16
    2.222.39−9.98−11.0393.1817.04110.2217.92−9.518.41
    2.50−0.11−5.68−8.5777.1319.5796.7020.53−2.8717.66
    下载: 导出CSV 
    | 显示表格

    2.3   失稳机制

    2.3.1   失稳模式I (振荡失稳)

    由表4可知, 当0.40≤As≤1.18时, 失稳临界频率不等于零, 即热毛细对流从定常二维轴对称流动失稳变为三维周期性振荡对流. 随着高径比的增加, 失稳的临界波数从kc = 6逐渐减少为kc = 2. 由表5和图3可知, 在失稳模式I中, 失稳主要是由水动力学惯性机制作用, 其中IV4和IV5是造成水动力学惯性失稳的主要原因; 在这一失稳模式中, 热毛细机制对失稳也有一定贡献. 通过具体分析, 失稳模式I在不同高径比范围内, IV4项和IV5项对失稳的贡献略微有所不同. 当0.40≤As≤0.59时, 液桥流动失稳的能量主要是由IV4和IV5提供(约占100%), 且二者比例相当; IV2和IV3为负值, 代表其对对流具有稳定作用; 同时扰动引起的热毛细力做功(M项)的贡献占比小于8%(表5). 当0.67≤As≤1.00时, IV5项对液桥流动失稳的贡献在60%左右, 是造成流动失稳的主要原因; IV2 已经为正值, IV3仍然为负值, 热毛细力做功的贡献占比为6.52%~9.86%. 当1.05≤As≤1.18时, 液桥流动失稳主要是由IV4项引起, 其贡献在65%左右, 在这一高径比区间Mz和Mφ的对失稳的贡献有所增加, 热毛细力做功的贡献占比为15.31%~27.35%, 即热毛细机制有所增强. 下面将详细分析上述3个高径比区间的失稳机制.

    当0.40≤As≤0.59时, 由表5可知, 在这一高径比区间中, Mz和Mφ对液桥的流动失稳贡献较小, IV项(IV项表示IV1项到IV5项的总和)是造成流动失稳的主要原因. 由图4(a)扰动能量平衡图可以直观的看出, 失稳主要由IV4和IV5造成, 即基态的轴向速度径向和轴向梯度为引起失稳的扰动动能传递了主要的能量. 由局部扰动动能极值图(图4(b))可知, 液桥热毛细流动最危险的失稳区域出现在距液桥自由液面R/3处并靠近热壁附近. 在液桥自由液面的扰动温度图(图5(b))上分布了5个热极和5个冷极, 其中红色表示热极, 蓝色表示冷极, 自由液面的扰动速度从热极指向冷极, 与扰动热毛细力的方向相同. 从图5(a)中可以看出, 在靠近自由液面约R/3处出现扰动温度强斑, 截面处的扰动速度从热极指向冷极, 并在冷极处形成扰动速度流胞.

    下载: 全尺寸图片 幻灯片
    图  4  高径比As = 0.48, Pr = 0.068的临界状态(Mac = 1936.3)下:(a)扰动能量平衡条状图, (b)局部扰动动能极值图
    Figure  4.  (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 1936.3) at As = 0.48, Pr = 0.068
    下载: 全尺寸图片 幻灯片
    图  5  高径比As = 0.48, Pr = 0.068的临界状态 (Mac = 1936.3) 下:(a) z = 0.51截面扰动温度 (云图) 和扰动速度矢量 (箭头), (b)自由液面扰动温度 (云图) 和扰动速度矢量 (箭头)
    Figure  5.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow): (a) the z = 0.51 cross section, (b) the free surface under the critical condition (Mac = 1936.3) at As = 0.48, Pr = 0.068
    下载: 全尺寸图片 幻灯片

    在0.67≤As≤1区域内, 由图6的扰动能量平衡图可知, 当高径比在这一范围内, 失稳依然是由水动力学惯性机制占主导, 但与0.40≤As≤0.59这一高径比范围内失稳的具体扰动动能项存在一些差异, 该高径比范围内液桥的流动失稳主要是由IV5项主导, 即基态轴向速度的轴向梯度为流动失稳提供了主要能量. 流动最危险的失稳位置出现在靠近热壁处距离自由液面R/2附近(图6(b)). 如图7(a)中扰动温度所示, 该高径比范围内扰动温度的特征与0.40≤As≤0.59时不同, 在这一高径比范围内, 扰动温度的特征是: 在靠近对称轴R/3处和靠近自由液面附近出现强斑. 在自由液面处扰动速度从热极指向冷极, 与扰动热毛细力的方向相同(图7(b)).

    图  6  高径比As = 0.77, Pr = 0.068的临界状态 (Mac = 1083.9) 下: (a)扰动能量平衡条状图, (b)局部扰动动能极值图
    Figure  6.  (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 1083.9) at As = 0.77, Pr = 0.068
    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片
    图  7  高径比As = 0.77, Pr = 0.068的临界状态 (Mac = 1083.9) 下: (a) z = 0.51截面扰动温度 (云图) 和扰动速度矢量 (箭头), (b)自由液面扰动温度 (云图) 和扰动速度矢量 (箭头)
    Figure  7.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow): (a) the z = 0.51 cross section, (b) the free surface under the critical condition (Mac = 1083.9) at As = 0.77, Pr = 0.068
    下载: 全尺寸图片 幻灯片

    当1.05≤As≤1.18时, 由扰动能量平衡图(图8(a))可知, 在这一高径比范围内, 热毛细机制对失稳的贡献增大, 但液桥热毛细流动失稳仍是由水动力学的惯性机制主导. 与上述两种情况不同, 在这一高径比范围内IV4项是造成失稳的主要原因, 但 IV2项对流动失稳有较为明显的抑制作用, 即基态轴向速度的径向梯度为流动失稳提供了主要能量, 而基态径向速度的轴向梯度有抑制流动失稳的作用. 根据图8(b)可知, 最危险的失稳位置出现在热壁附近, 并距离对称轴R/3处. 从自由液面展开图(图9(a))显示: 在靠近冷壁附近扰动速度从热极指向冷极, 而在靠近热壁附近扰动速度从冷极指向热极. 根据z = 0.20,z = 0.51和z = 0.74截面扰度温度和扰动速度矢量图(图9(b)~图9(d))可知, 在z = 0.20截面处扰动速度从冷极指向热极附近, 扰动温度模式的特征为轮辐状强斑, 并且布满整个液桥横截面; 在z = 0.51截面处, 扰动速度在靠近对称轴附近是由热热极指向冷极, 在靠近自由液面附近是由冷极指向热极, 从而扰动速度在热极和冷极之间形成一个扰动流动循环, 扰动温度模式特征为: 在靠近自由液面附近形成椭圆形强斑; 在z = 0.74截面中扰动速度从热极指向冷极, 扰动温度特征为: 在自由液面附近和靠近自由液面R/3处形成强斑.

    图  8  高径比As = 1.11, Pr = 0.068的临界状态 (Mac = 1273.6) 下: (a)扰动能量平衡条状图, (b)局部扰动动能极值图
    Figure  8.  (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 1273.6) at As = 1.11, Pr = 0.068
    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片
    图  9  高径比 As = 1.11，Pr = 0.068的临界状态(Mac = 1273.6)下扰动温度(云图)和扰动速度矢量(箭头)
    Figure  9.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow) under the critical condition (Mac = 1273.6) for the case of As = 1.11 and Pr = 0.068
    下载: 全尺寸图片 幻灯片
    2.3.2   失稳模式II (静态失稳)

    根据表4可知, 在1.20≤As≤2.50区域内, 失稳临界频率等于零, 表示在这一高径比范围内, 热毛细对流从二维轴对称定常对流失稳变为三维定常对流. 随着高径比的增加, 失稳的临界波数kc = 2转变为kc = 1. 由表5和图3可知, 液桥的流动失稳依然是由IV项主导, 这表明失稳模式II仍然是水动力学惯性失稳机制主导, 其中M项在失稳中的占比随高径比的增加先增大, 再减小最后再增大(图3(a)). 在失稳模式II中IV4和IV5是引起液桥流动失稳的主要原因, 但在不同的高径比下 IV4项和IV5项的占比存在差异. 在1.20≤As≤1.54范围内, 液桥的流动失稳是由IV4项主导(约占43%), 但IV5项对失稳的贡献同样较大(约占28%), M项的占比在13.39%~22.18%;在1.67≤As≤2.50, IV4项对失稳的贡献增大, 而IV5项和M项的在失稳中的占比减小, 同时IV2 和IV3变为负值, 所以IV4项是引起流动失稳的主要原因. 下面对失稳模式II进行具体分析.

    在1.20≤As≤1.54范围, 失稳临界波数kc = 2. 根据图10(a)所示, IV4和IV5是造成流动失稳的主原因, 即基态轴向速度的径向梯度和轴向梯度为引起失稳的扰动动能传递了主要的能量, 其中IV4项比IV5项的贡献更大, Mz和Mφ项在流动失稳中的占比也相对较大. 流动失稳最危险的位置出现在液桥中部偏上附近并靠近自由液面处(图10(b)). 从自由液面展开图(图11(a))可知, 扰动速度在靠近冷壁处是由热极指向冷极, 在热壁附近则是由冷极指向热极, 这与图12所示的低普朗特数(Pr = 0.009)液桥在自由液面处扰动速度的指向不同, 在低普朗特数液桥中自由有液面的扰动速度都是从冷极指向热极,且失稳都是水动力学惯性机制作用下的失稳, 热毛细机制对失稳没有贡献. 由图11(b)~图11(d)可知, z = 0.20截面处扰动温度的特征为: 靠近自由液面R/3形成强斑,截面处的扰动速度由冷极流向热极; 在z = 0.51截面处, 扰动温度特征为: 在自由液面R/3形成椭圆形强斑, 并衍生出来的弱斑向自由液面靠近, 截面处的扰动速度在靠近对称轴处是由热级指向冷极, 在靠近自由液面附近是由冷极指向热极; 在z = 0.74截面处, 扰动温度特征为: 扰动温度强斑向自由液面靠近, 截面处的扰动速度是由热极流向冷极.

    图  10  高径比As = 1.43, Pr = 0.068的临界状态 (Mac = 531.9) 下: (a)扰动能量平衡条状图, (b) 局部扰动动能极值图
    Figure  10.  (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 531.9) at As= 1.43, Pr = 0.068
    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片
    图  11  高径比As = 1.43, Pr = 0.068的临界状态(Mac = 531.9)下扰动温度(云图)和扰动速度矢量(箭头)
    Figure  11.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow) under the critical condition (Mac = 531.9) for the case of As = 1.43 and Pr = 0.068
    下载: 全尺寸图片 幻灯片
    图  12  高径比As = 1.43, Pr = 0.009临界状态(Mac = 531.9)下自由液面处扰动温度(云图)和扰动速度矢量(箭头)
    Figure  12.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow) at the free surface under the critical condition (Mac = 531.9) for the case of As = 1.43 and Pr = 0.009
    下载: 全尺寸图片 幻灯片

    在1.67≤As≤2.50区间内, 失稳临界波数kc = 1. 由图13(a)可知, IV4是造成流动失稳的主要原因, 即基态轴向速度的径向梯度为流动失稳提供了最主要能量, 而Mz也对液桥的流动失稳有所影响, 这表明失稳仍然是由水动力学惯性机制主导, 同时热毛细机制对失稳有所影响. 与上述所有失稳机制不同的是, 在该高径比范围内Mφ项起抑制流动失稳的作用. 由图13(b)可知, 最危险的失稳位置出现在液桥中部并靠近自由液面. 在图14(a)中, 扰动速度在靠近液桥冷壁处是由热极指向冷极, 在靠近液桥热壁处扰动速度从冷极指向热极, 在1.05≤As≤1.18和1.20≤As≤1.54这两个高径比范围内都曾出现这一现象, 但在该高径比范围内这一现象更加明显, 扰动速度更多地从冷极指向热极. 由图14(a)和图15对比分析可知, 低普朗特数(Pr = 0.009)液桥自由液面的扰动速度是从冷极指向热极, 这与GaAs熔体液桥自由液面的扰动速度指向不同, 由图14(b)~图14(d)可知, 在z = 0.20, z = 0.51, z = 0.74的横截面的边缘处, 扰动速度都是从冷极指向热极, 同样验证了自由液面处扰动速度更多的是从冷极指向热极; 这三处横截面的扰动温度特征为: 随着z坐标的增大, 扰动温度强斑逐渐向自由液面靠近.

    图  13  高径比 As = 2.00，Pr = 0.068的临界状态 (Mac = 234.6) 下: (a)扰动能量平衡条状图, (b)局部扰动动能极值图
    Figure  13.  (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 234.6) at As= 2.00, Pr = 0.068
    下载: 全尺寸图片 幻灯片
    图  14  高径比As = 2.00, Pr = 0.068的临界状态(Mac = 234.6)下扰动温度(云图)和扰动速度矢量(箭头)
    Figure  14.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow) under the critical condition (Mac = 234.6) for the case of As = 2.00 and Pr = 0.068
    下载: 全尺寸图片 幻灯片
    图  15  As = 2.00 Pr = 0.009时临界状态(Mac = 234.6)下自由液面处扰动温度(云图)和扰动速度矢量(箭头)
    Figure  15.  The disturbance temperature (colored contour) and disturbance velocity vector (arrow) at the free surface under the critical condition (Mac = 234.6) for the case of As = 2.00 and Pr = 0.009
    下载: 全尺寸图片 幻灯片
    2.3.3   基态流函数和温度分布

    由图16和图17可知, 在失稳模式I和失稳模式II中, 随着高径比的增大基态涡胞的位置始终靠近自由液面附近. 从图16和图17可以看出, 当高径比较小时, 温度场在靠近自由液面R/3处向热壁处挤压, 同时温度场也在自由液面处向冷壁处挤压, 从而形成较大温度梯度; 但随着高径比增加, 这种温度场的挤压相对减缓, 从而随着高径比的增加温度场沿轴向分布相对均匀.

    下载: 全尺寸图片 幻灯片
    图  16  失稳模式I在临界状态下的基态流函数(云图)和基态温度场(等值线)
    Figure  16.  The basic flow (colored contour) and temperature (isotherm) of instability mode I under critical condition
    下载: 全尺寸图片 幻灯片
    图  17  失稳模式II在临界状态下的基态流函数(云图)和基态温度场(等值线)
    Figure  17.  The basic flow (colored contour) and temperature (isotherm) of instability mode II under critical condition
    下载: 全尺寸图片 幻灯片
    图  17  失稳模式II在临界状态下的基态流函数(云图)和基态温度场(等值线) (续)
    Figure  17.  The basic flow (colored contour) and temperature (isotherm) of instability mode II under critical condition (continued)
    下载: 全尺寸图片 幻灯片

    3.   结论

    基于已有文献可知,液桥模型中熔体普朗特数对热毛细对流失稳机制和失稳模式具有重要影响: 对典型高普朗特数熔体(例如Pr>1),热毛细对流第一次失稳是从二维轴对称定常对流转变三维周期性振荡对流, 其失稳机制是热毛细失稳机制; 对典型低普朗特数熔体(例如Pr = 0.011), 热毛细对流第一次失稳是从二维轴对称定常对流转变为三维定常对流, 失稳机制是水动力学惯性失稳. 同时文献工作显示液桥高径比对热毛细流对流的失稳临界参数具有重要影响, 但是未发现高径比会影响热毛细对流失稳机制和失稳模式(三维静态失稳和三维振荡失稳).

    本文在基于谱元法线性稳定性分析对GaAs熔体(Pr = 0.068)液桥热毛细对流研究发现, 在同一普朗特数不同高径比的情况下, 随着高径比的改变液桥出现两种失稳模式: 当0.4≤As≤1.18时, 热毛细对流是由二维轴对称定常对流直接转变为三维周期性振荡对流(振荡失稳); 当1.20≤As≤2.5时, 热毛细对流是由二维轴对称定常对流转变为三维定常流动(静态失稳). 基于能量分析我们发现, GaAs熔体液桥的热毛细对流失稳是由水动力学惯性机制占主导, 两种机制(水动力学惯性机制和热毛细机制)共同作用造成. 两种机制对失稳的贡献随着高径比的变化而变化, 其失稳现象更加复杂, 从而随着高径比的改变出现了两种不同的失稳模式. 由IV项各项的占比可知, IV4和IV5是造成流动失稳的主要原因, 即基态轴向速度径向梯度和轴向梯度为扰动场的失稳转递了主要的能量. 在失稳模式I(振荡失稳)中, 最危险的失稳位置出现在热壁附近, 并随着高径比的增加逐渐向对称轴靠近; 在失稳模式II(静态失稳)中, 最危险的失稳位置液桥中部横截面并靠近自由液面附近.

  • 图  1   半浮区液桥模型

    Figure  1.   Floating half-zone model

    下载: 全尺寸图片 幻灯片

    图  2   Mac随高径比(As)的变化(Pr = 0.068)

    Figure  2.   The variation of Mac versus the aspect ratio As for Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  3   扰动能量随高径比As变化

    Figure  3.   The variation of the kinetic energy budgets with respect to the aspect ratio

    下载: 全尺寸图片 幻灯片

    图  3   扰动能量随高径比As变化(续)

    Figure  3.   The variation of the kinetic energy budgets with respect to the aspect ratio (continued)

    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片

    图  4   高径比As = 0.48, Pr = 0.068的临界状态(Mac = 1936.3)下:(a)扰动能量平衡条状图, (b)局部扰动动能极值图

    Figure  4.   (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 1936.3) at As = 0.48, Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  5   高径比As = 0.48, Pr = 0.068的临界状态 (Mac = 1936.3) 下:(a) z = 0.51截面扰动温度 (云图) 和扰动速度矢量 (箭头), (b)自由液面扰动温度 (云图) 和扰动速度矢量 (箭头)

    Figure  5.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow): (a) the z = 0.51 cross section, (b) the free surface under the critical condition (Mac = 1936.3) at As = 0.48, Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  6   高径比As = 0.77, Pr = 0.068的临界状态 (Mac = 1083.9) 下: (a)扰动能量平衡条状图, (b)局部扰动动能极值图

    Figure  6.   (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 1083.9) at As = 0.77, Pr = 0.068

    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片

    图  7   高径比As = 0.77, Pr = 0.068的临界状态 (Mac = 1083.9) 下: (a) z = 0.51截面扰动温度 (云图) 和扰动速度矢量 (箭头), (b)自由液面扰动温度 (云图) 和扰动速度矢量 (箭头)

    Figure  7.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow): (a) the z = 0.51 cross section, (b) the free surface under the critical condition (Mac = 1083.9) at As = 0.77, Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  8   高径比As = 1.11, Pr = 0.068的临界状态 (Mac = 1273.6) 下: (a)扰动能量平衡条状图, (b)局部扰动动能极值图

    Figure  8.   (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 1273.6) at As = 1.11, Pr = 0.068

    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片

    图  9   高径比 As = 1.11，Pr = 0.068的临界状态(Mac = 1273.6)下扰动温度(云图)和扰动速度矢量(箭头)

    Figure  9.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow) under the critical condition (Mac = 1273.6) for the case of As = 1.11 and Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  10   高径比As = 1.43, Pr = 0.068的临界状态 (Mac = 531.9) 下: (a)扰动能量平衡条状图, (b) 局部扰动动能极值图

    Figure  10.   (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 531.9) at As= 1.43, Pr = 0.068

    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片

    图  11   高径比As = 1.43, Pr = 0.068的临界状态(Mac = 531.9)下扰动温度(云图)和扰动速度矢量(箭头)

    Figure  11.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow) under the critical condition (Mac = 531.9) for the case of As = 1.43 and Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  12   高径比As = 1.43, Pr = 0.009临界状态(Mac = 531.9)下自由液面处扰动温度(云图)和扰动速度矢量(箭头)

    Figure  12.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow) at the free surface under the critical condition (Mac = 531.9) for the case of As = 1.43 and Pr = 0.009

    下载: 全尺寸图片 幻灯片

    图  13   高径比 As = 2.00，Pr = 0.068的临界状态 (Mac = 234.6) 下: (a)扰动能量平衡条状图, (b)局部扰动动能极值图

    Figure  13.   (a) The disturbance energy balance and (b) distribution of the local disturbance kinetic energy extreme value under critical condition (Mac = 234.6) at As= 2.00, Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  14   高径比As = 2.00, Pr = 0.068的临界状态(Mac = 234.6)下扰动温度(云图)和扰动速度矢量(箭头)

    Figure  14.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow) under the critical condition (Mac = 234.6) for the case of As = 2.00 and Pr = 0.068

    下载: 全尺寸图片 幻灯片

    图  15   As = 2.00 Pr = 0.009时临界状态(Mac = 234.6)下自由液面处扰动温度(云图)和扰动速度矢量(箭头)

    Figure  15.   The disturbance temperature (colored contour) and disturbance velocity vector (arrow) at the free surface under the critical condition (Mac = 234.6) for the case of As = 2.00 and Pr = 0.009

    下载: 全尺寸图片 幻灯片
    下载: 全尺寸图片 幻灯片

    图  16   失稳模式I在临界状态下的基态流函数(云图)和基态温度场(等值线)

    Figure  16.   The basic flow (colored contour) and temperature (isotherm) of instability mode I under critical condition

    下载: 全尺寸图片 幻灯片

    图  17   失稳模式II在临界状态下的基态流函数(云图)和基态温度场(等值线)

    Figure  17.   The basic flow (colored contour) and temperature (isotherm) of instability mode II under critical condition

    下载: 全尺寸图片 幻灯片

    图  17   失稳模式II在临界状态下的基态流函数(云图)和基态温度场(等值线) (续)

    Figure  17.   The basic flow (colored contour) and temperature (isotherm) of instability mode II under critical condition (continued)

    下载: 全尺寸图片 幻灯片

    表  1   基于不同网格分辨率得到的热毛细对流(Pr = 0.068)失稳临界Marangoni数(Mac)和波数kc

    Table  1   The critical Marangoni number (Mac) and wave number kc for the instability of thermocapillary flow with Pr = 0.068 under different mesh resolution conditions

    Mesh resolution41 × 4649 × 5557 × 64
    Mac / kc1179.7/31182.4/31182.7/3
    下载: 导出CSV

    表  2   计算失稳临界参数(Rec, fc和kc)与文献[14]结果对比

    Table  2   The critical values (Rec, fc, kc) from the present computations against Ref. [14]

    PrRec/ fc/kc (present)Rec/ fc/kc [14]Relative error of Rec/fc
    0.06817379/205.3/317229/206/30.87% / 0.36%
    0.18314021/409.1/314276/415/31.78%/1.42%
    1.02544/64.9/22551/65/20.27%/0.15%
    下载: 导出CSV

    表  3   计算失稳临界参数(Mac和kc)与文献[15]结果对比

    Table  3   The critical values (Mac and kc) from the present computations against Ref. [15]

    AsMac/kc (present)Mac/kc[15] (LSA)Mac/kc[15] (DNS)
    0.632.27/332.21/334.73/3
    0.823.12/223.18/224.71/2
    1.016.76/217.00/219.77/2
    下载: 导出CSV

    表  4   Pr = 0.068时不同高径比下的失稳临界参数

    Table  4   The critical values at different aspect ratios for Pr = 0.068

    AsMacfckcType
    0.402033.1±490.76



    0.431920.0±493.65
    0.481936.3±507.45
    0.531770.8±504.74
    0.591838.4±536.14
    0.671988.8±233.04type I (oscillatory instability)
    0.711278.0±200.34
    0.771083.9±201.74
    0.91995.3±155.23
    1.001182.4±205.33
    1.051246.8±63.62

    1.111273.6±62.72
    1.181393.6±36.72
    1.201348.40.02
    1.251058.60.02type II (stationary instability)
    1.33691.30.02
    1.43531.90.02
    1.54483.60.02
    1.67350.50.01
    1.82262.80.01
    2.00234.60.01
    2.22248.40.01
    2.50326.00.01
    下载: 导出CSV

    表  5   Pr = 0.068时不同高径比能量分析结果

    Table  5   The energy analysis results at different aspect ratios for Pr = 0.068

    AsIV1/%IV2/%IV3/%IV4/%IV5/%(IV4 + IV5)/%Mz/%Mφ/%(Mz + Mφ)/%
    0.402.45−7.90−1.4451.6947.9299.613.543.677.21
    0.432.17−12.31−1.4559.6344.46104.094.043.347.38
    0.484.39−6.12−1.8148.9848.3897.363.112.996.10
    0.534.14−11.26−1.8458.6944.07102.763.552.536.08
    0.596.69−3.11−2.5344.1949.8194.002.592.224.81
    0.670.639.55−1.1118.1963.8782.062.576.178.74
    0.71−0.1810.44−2.6720.4062.9583.352.496.498.98
    0.770.8512.58−4.1619.4163.6583.061.865.757.61
    0.910.339.57−4.2922.3959.0681.453.366.509.86
    1.007.1016.23−5.7514.4361.5375.961.485.046.52
    1.05−1.95−23.35−1.0571.5926.6598.2418.958.4027.35
    1.112.49−16.152.0065.2123.2788.4815.517.0622.57
    1.189.54−6.576.3455.9019.2275.1210.954.3615.31
    1.2010.29−2.227.5350.6020.4171.019.733.6613.39
    1.254.911.527.6644.2827.3371.619.144.9314.07
    1.33−2.441.758.7843.1229.8672.9810.728.1818.90
    1.43−4.711.349.2243.8028.8172.6110.8810.6421.52
    1.54−4.842.519.1843.2727.7170.9810.0412.1422.18
    1.675.33−23.21−13.15117.748.12125.8616.78−11.675.11
    1.825.42−19.76−13.34118.697.33126.0215.32−13.651.67
    2.004.22−14.91−12.64108.4611.71120.1716.04−12.883.16
    2.222.39−9.98−11.0393.1817.04110.2217.92−9.518.41
    2.50−0.11−5.68−8.5777.1319.5796.7020.53−2.8717.66
    下载: 导出CSV
    • 参考文献(32)
  • [1]

    Schwabe D, Scharmann A. Marangoni convection in open boat and crucible. Journal of Crystal Growth, 1981, 52: 435-449 doi: 10.1016/0022-0248(81)90231-1
    [2]

    Azami T, Nakamura S, Eguchi M, et al. The role of surface-tension-driven flow in the formation of a surface pattern on a Czochralski silicon melt. Journal of Crystal Growth, 2001, 233(1): 99-107
    [3]

    Hu KX, Yan CY, Chen QS. Instability of thermocapillary–buoyancy convection in droplet migration. Physics of Fluids, 2019, 31(12): 122101 doi: 10.1063/1.5125846
    [4] 章绍能, 胡开鑫. 黏弹性液滴热毛细迁移的对流不稳定. 力学学报, 2021, 53(5): 1313-1323 (Zhang Shaoneng, Hu Kaixin. Convective instability in thermocapillary migration of a viscoelastic droplet. Chinese Journal of Theoretical and Applied Mechanics, 2021, 53(5): 1313-1323 doi: 10.6052/0459-1879-20-443
    [5]

    Louvis E, Fox P, Sutcliffe CJ. Selective laser melting of aluminium components. Journal of Materials Processing Tech, 2011, 211(2): 275-284 doi: 10.1016/j.jmatprotec.2010.09.019
    [6]

    Dai D, Gu D. Thermal behavior and densification mechanism during selective laser melting of copper matrix composites: Simulation and experiments. Materials and Design, 2014, 55: 482-491 doi: 10.1016/j.matdes.2013.10.006
    [7]

    Smith MK, Davis SH. Instabilities of dynamic thermocapillary liquid layers Part 1. Convective instabilities. Journal of Fluid Mechanics, 1983, 132(1): 145-162
    [8]

    Levenstam M, Amberg G. Hydrodynamical instabilities of thermocapillary flow in a half-zone. Journal of Fluid Mechanics, 1995, 297(1): 357-372
    [9]

    Levenstam M, Amberg G, Carlberg T, et al. Experimental and numerical studies of thermocapillary convection in a floating zone like configuration. Journal of Crystal Growth, 1996, 158(3): 224-230 doi: 10.1016/0022-0248(95)00466-1
    [10]

    Chen QS, Hu WR, Prasad V. Effect of liquid bridge volume on the instability in small-Prandtl-number half zones. Journal of Crystal Growth, 1999, 203(1-2): 261-268 doi: 10.1016/S0022-0248(99)00064-0
    [11]

    Ichiro U, Shiho T, Hiroshi K. Oscillatory and chaotic thermocapillary convection in a half-zone liquid bridge. Physics of Fluids, 2003, 15(2): 408-416 doi: 10.1063/1.1531993
    [12]

    Wanschura M, Shevtsova VM, Kuhlmann HC, et al. Convective instability mechanisms in thermocapillary liquid bridges. Physics of Fluids, 1995, 7(5): 912-925 doi: 10.1063/1.868567
    [13]

    Chen G, Lizée A, Roux B. Bifurcation analysis of the thermocapillary convection in cylindrical liquid bridges. Journal of Crystal Growth, 1997, 180: 638-647
    [14]

    Levenstam M, Amberg G, Winkler C. Instabilities of thermocapillary convection in a half-zone at intermediate Prandtl numbers. Physics of Fluids, 2001, 13(4): 807-816 doi: 10.1063/1.1337063
    [15]

    Li K, Xun B, Imaishi N, et al. Thermocapillary flows in liquid bridges of molten tin with small aspect ratios. International Journal of Heat & Fluid Flow, 2008, 29(4): 1190-1196
    [16]

    Rybicki A, Floryan JM. Thermocapillary Effects in Liquid Bridges I. Thermocapillary Convection. Physics of Fluids, 1998, 30(7): 1956-1972
    [17]

    Velten R, Schwabe D, Scharmann A. The periodic instability of thermocapillary convection in cylindrical liquid bridges. Physics of Fluids A Fluid Dynamics, 1991, 3(2): 267-279 doi: 10.1063/1.858135
    [18]

    Chen QS, Hu WR. Influence of liquid bridge volume on instability of floating half zone convection. International Journal of Heat & Mass Transfer, 1998, 41(6-7): 825-837
    [19]

    Nishino K, Yano T, Kawamura H. Instability of thermocapillary convection in long liquid bridges of high Prandtl number fluids in microgravity. Journal of Crystal Growth, 2015, 420: 57-63 doi: 10.1016/j.jcrysgro.2015.01.039
    [20] 吴勇强, 段俐, 李永强等. 大普朗特数大液桥浮力−热毛细对流地面实验. 力学学报, 2012, 44(6): 981-989 (Wu Yongqiang, Duan Li, Li Yongqiang, et al. Ground experiments of bouyant thermocapillary convection of large scale liquid bridge with large prandtl number. Chinese Journal of Theoretical and Applied Mechanics, 2012, 44(6): 981-989 doi: 10.6052/0459-1879-12-148
    [21] 王佳, 吴笛, 段俐等. 大尺寸液桥热毛细对流失稳性地面实验研究. 力学学报, 2015, 47(4): 580-586 (Wang Jia, Wu Di, Duan Li, et al. Ground experiments on the instability of thermocapillary convection in large scale liquid bridge. Chinese Journal of Theoretical and Applied Mechanics, 2015, 47(4): 580-586 doi: 10.6052/0459-1879-14-309
    [22]

    Kang Q, Wu D, Duan L, et al. The effects of geometry and heating rate on thermocapillary convection in the liquid bridge. Journal of Fluid Mechanics, 2019, 881: 951-982 doi: 10.1017/jfm.2019.757
    [23]

    Wang J, Wu D, Duan L, et al. Ground experiment on the instability of buoyant-thermocapillary convection in large-scale liquid bridge with large Prandtl number. International Journal of Heat & Mass Transfer, 2017, 108: 2107-2119
    [24]

    Liu H, Zeng Z, Yin L, et al. Influence of aspect ratio on the onset of thermocapillary flow instability in annular pool heated from inner wall. International Journal of Heat and Mass Transfer, 2019, 129: 746-752 doi: 10.1016/j.ijheatmasstransfer.2018.10.016
    [25]

    Yin LM, Zhong Z, Qiu Z, et al. Linear stability analysis of thermocapillary flow in a slowly rotating shallow annular pool using spectral element method. International Journal of Heat & Mass Transfer, 2016, 97: 353-363
    [26]

    Lehoucq RB, Sorensen DC. Deflation techniques for an implicitly restarted arnoldi iteration. Siam Journal on Matrix Analysis & Applications, 1996, 17(4): 789-821
    [27]

    Lehoucq RB, Sorensen DC, Yang C. ARPACK users' guide: solution of large-scale eigenvalue problems with implicitly restarted Arnoldi methods. Siam, 1998
    [28]

    Lehoucq RB. Implicitly restarted arnoldi methods and subspace iteration. SIAM Journal on Matrix Analysis and Applications, 23 (2): 551-562, 2001
    [29]

    Boppana V, Gajjar J. Global flow instability in a lid-driven cavity. International Journal for Numerical Methods in Fluids, 2010, 62(8): 827-853
    [30]

    Prange M, Wanschura M, Kuhlmann HC, et al. Linear stability of thermocapillary convection in cylindrical liquid bridges under axial magnetic fields. Journal of Fluid Mechanics, 1999, 394: 281-302 doi: 10.1017/S0022112099005698
    [31]

    Liu H, Zeng Z, Yin L, et al. Effect of the Prandtl number on the instabilities of the thermocapillary flow in an annular pool. Physics of Fluids, 2019, 31(3): 034103 doi: 10.1063/1.5087113
    [32]

    Zeng Z, Mizuseki H, Shimamura K. Marangoni convection in model of floating zone under mircogravity. Journal of Crystal Growth, 2001, 229: 601-604 doi: 10.1016/S0022-0248(01)01236-2
    • 相关文章
    • 施引文献(1)

  • 期刊类型引用(1)

    1. 倪宝玉，熊航，曾令东，孙林华. 移动压力面诱导的固壁冰航道内弯曲重力波传播特性研究. 力学学报. 2024(10): 2838-2853 . 本站查看

    其他类型引用(0)

    • 资源附件(0)
图(24)  /  表(5)
计量
  • 文章访问数:  922
  • HTML全文浏览量:  384
  • PDF下载量:  85
  • 被引次数: 1
出版历程
  • 收稿日期:  2021-05-25
  • 录用日期:  2021-11-10
  • 网络出版日期:  2021-11-11
  • 刊出日期:  2022-02-16

目录

摘要
HTML全文

  • 引言

  • 1.   数学模型和数值方法

  • 1.1   物理模型与基本控制方程

  • 1.2   线性稳定性分析

  • 1.3   扰动能量分析

  • 2.   结果与讨论

  • 2.1   网格依赖性和程序有效性验证

  • 2.2   稳定性曲线与失稳模式

  • 2.3   失稳机制

  • 2.3.1   失稳模式I (振荡失稳)
  • 2.3.2   失稳模式II (静态失稳)
  • 2.3.3   基态流函数和温度分布

  • 3.   结论

参考文献(32)
相关文章
施引文献(1)
资源附件(0)

/

返回文章
返回

下载中心

作者中心

版权所有 © 《力学学报》编辑部备案号：京ICP备05039218号-1

通讯地址：北京市海淀区北四环西路15号邮政编码：100190

联系电话：010-62536271

邮箱：lxxb@cstam.org.cn

微信公众号
RSS

本系统由 北京仁和汇智信息技术有限公司 开发  百度统计