液桥蒸发动力学过程的格子Boltzmann数值模拟研究

Lattice Boltzmann Simulation of the Evaporation Kinetics of a Liquid Bridge

摘要: 本文采用格子Boltzmann方法研究了相对湿度、基板宽度和液桥高度对液桥蒸发动力学过程的影响。模拟结果表明，提高湿度可以减小液桥的蒸发速度，延长液桥断裂时刻。增大液桥高度可以提高液桥的蒸发速度，减小液桥断裂时刻。基板宽度对液桥蒸发速度可以忽略，因此液桥断裂时刻不变。在蒸发初期液桥颈部与基板边缘钉扎而液桥腰部收缩，使得液桥与基板间接触角减小，液桥以固定接触线的模式蒸发。当液桥接触角减小到足够小后，液桥接触线与基板边缘脱钉，液桥以固定接触角模式蒸发。环境湿度越大，液桥接触线脱钉所需时间越长。液桥高度越高，液桥接触线脱钉所需时间越短。而基板宽度对液桥接触线脱钉所需时间的影响可以忽略。

Abstract: In this paper, the lattice Boltzmann method was used to study the effects of humidity, substrate width and liquid bridge height on the evaporation kinetics of liquid bridges. The simulation results show that increasing the relative humidity can reduce the evaporation velocity of the liquid bridge and prolong the breaking time of the liquid bridge. Increasing the height of the liquid bridge can increase the evaporation rate of the liquid bridge and reduce the breaking time of the liquid bridge. The width of the substrate has a negligible effect on the evaporation rate of the liquid bridge, so the time of fracture of the liquid bridge remains unchanged. In the initial stage of evaporation, the neck of the liquid bridge is pinned to the edge of the substrate and the waist of the liquid bridge is contracted, so that the contact angle between the liquid bridge and the substrate is reduced, and the liquid bridge evaporates in the mode of fixed contact line. When the contact angle of the liquid bridge is reduced to a small enough level, the contact line of the liquid bridge detaches from the edge of the substrate, and the liquid bridge evaporates in a fixed contact angle mode. The higher the ambient humidity, the longer it takes for the contact line of the liquid bridge to detach. The higher the height of the bridge, the shorter the time required for the contact line of the bridge to detach. However, the influence of the width of the substrate on the time required for the detachment of the contact line of the liquid bridge is negligible.

1. 引言

在多相流和微流体动力学的研究中，液桥因其独特的物理性质和广泛的应用而备受关注。液桥是指在两个固体表面之间，由于表面张力的作用而形成的液体柱。在“天宫二号”空间实验室中，科学家们研究了液桥在微重力环境下的行为，这对于理解液体在失重状态下的物理特性和太空探索具有重要意义[1]。在增材制造领域，液桥技术在3D打印技术中发挥着重要作用，Kim等人[2]通过精确控制液桥的蒸发和凝固，成功将电致变色材料打印在AR眼镜中实现了导航的功能。生物技术领域中，Zhang等人[3]通过控制液桥的破裂实现了微液滴的制造，利用该方法制造的微液滴可用于单细胞分析和精准给药，极大地促进了医学研究和治疗实践。液桥技术还可以在固体表面生成图案化的功能涂层，贾若琨等人[4]发现通过控制液桥内溶剂的散逸可以制造表面可控的阵列，进而制造出柔性印章的母模板。基于液桥的曲率，刘逸飞等人[5]通过对液桥边缘进行轮廓拟合，发明了一种液体表面张力系数的测量方法。Ma等人[6]在弯液面表面附加一层薄固体膜来创造固液界面，从而可以通过分析弯液面的曲率来计算固液界面自由能。这些应用中，对液桥蒸发过程的精确控制是一个关键的因素。张晟铭等人[7]通过Surface Evolver模拟研究了在剪切力作用下，两平行板之间液桥的液面形状及其与固体界面之间力的变化。但是模型简化了液面轮廓，理论值与模拟值差异较大。宋志超等人[8]使用VOF-IBM-DEM模型模拟发现小接触角、大质量的液桥有更大的颈部宽度。相对而言，格子Boltzmann方法（LBM）在多相流领域具有显著优势[9]。首先，LBM能够处理复杂的边界条件，如液体与固体板之间的接触角和表面张力效应[10]。其次，LBM基于粒子动力学，能够精确模拟微观尺度下的液体行为[11]。此外，LBM具有出色的并行计算能力，使其能够高效解决大规模的计算问题[12]，尤其是针对复杂几何和动态变化的液桥蒸发过程。本文采用格子Boltzmann方法(LBM)，研究了相对环境湿度$RH$ 和液桥尺寸(宽度$D$ 和高度$L$ )对其蒸发动力学过程的影响。通过研究不同条件下液桥的蒸发速度、蒸发过程中液桥尺寸的变化、接触角的变化和质量的变化，为液桥蒸发行为提供了深入的理解，从而为液桥的控制和应用提供科学依据。

2. 数值方法与问题描述

2.1. 数值方法

LBM通过模拟流体粒子的分布函数在离散格子上的碰撞和迁移来描述流体的流动过程。引入单松弛时间的BGK碰撞算子得到流体粒子密度分布函数的演化方程[13]，可表示为：

$f_i\left(x+e_i\Delta t,t+\Delta t\right)-f_i\left(x,t\right)=-\frac{1}{\tau }\left[f_i\left(x,t\right)-f_i^{eq}\left(x,t\right)\right]+\Delta f_i\left(x,t\right)$ (1)

式中$t$ 和$x$ 分别表示时间步和流体粒子的位置；$i$ 表示流体粒子的运动方向；$f_i\left(x,t\right)$ 表示在时间步$t$ 和位置$x$ 处的流体粒子在$i$ 方向的密度分布函数；$f_i^{eq}\left(x,t\right)$ 表示流体粒子的平衡态密度分布函数；$\tau$ 表示无量纲松弛时间；$\Delta f_i\left(x,t\right)$ 表示体积力项。其中$f_i^{eq}\left(x,t\right)$ 可表示为：

$f_i^{eq}\left(x,t\right)=w_i\rho \left(x,t\right)\left[1+\frac{e_i\cdot u}{c_s{}^2}+\frac{{\left(e_i\cdot u\right)}^2}{2c_s{}^4}-\frac{u\cdot u}{2c_s{}^2}\right]$ (2)

式中$c_s$ 表示格子声速，具体为$c_s=c/\sqrt{3}$ ，其中$c=\Delta t/\Delta x$ 表示格子速度，$\Delta t$ 和$\Delta x$ 分别表示时间步长和格子间距。在LBM中$\Delta t$ 和$\Delta x$ 通常都取1，故$c=1$ ，$c_s=\sqrt{3}/3$ 。$w_i$ 和$e_i$ 分别表示权重系数和离散速度，取值都与方向$i$ 和所使用的格子结构有关。对于本文所采用的D2Q9格子模型，权重系数$w_i$ 和离散速度$e_i$ 分别可表示为：

$w_i=\{\begin{array}{ll}4/9,\hfill & i=0\hfill \\ 1/9,\hfill & i=1,2,3,4\hfill \\ 1/36,\hfill & i=5,6,7,8\hfill \end{array}$ (3)

$e_i=\{\begin{array}{ll}\left(0,0\right),\hfill & i=0\hfill \\ c\left(\cos \left[\frac{\pi \left(i-1\right)}{2}\right],\sin \left[\frac{\pi \left(i-1\right)}{2}\right]\right),\hfill & i=1,2,3,4\hfill \\ \sqrt{2}c\left(\cos \left[\frac{\pi \left(2i-9\right)}{4}\right],\sin \left[\frac{\pi \left(2i-9\right)}{4}\right]\right),\hfill & i=5,6,7,8\hfill \end{array}$ (4)

式(2)中的$\rho \left(x,t\right)$ 和$u$ 分别表示在时间步$t$ 和位置$x$ 处的流体粒子的宏观密度和宏观速度，通过统计流体粒子的密度分布函数$f_i\left(x,t\right)$ 分别可表示为：

$\rho \left(x,t\right)={\displaystyle \sum _{i=0}^8{f}_i\left(x,t\right)}$ (5)

$u=\frac{1}{\rho }{\displaystyle \sum _{i=0}^8{f}_i\left(x,t\right)e_i}$ (6)

在LBM中液相的运动黏度和气相的运动黏度都由无量纲松弛时间$\tau$ 决定，即：

$\upsilon =c_s^2\left(\tau -\frac{1}{2}\right)\Delta t$ (7)

对于气液界面处的运动黏度，可表示为[14]：

$\upsilon ={\upsilon }_{\text{liquid}}\frac{\rho -{\rho }_{\text{vapor}}}{{\rho }_{\text{liquid}}-{\rho }_{\text{vapor}}}+{\upsilon }_{\text{vapor}}\frac{{\rho }_{\text{liquid}}-\rho }{{\rho }_{\text{liquid}}-{\rho }_{\text{vapor}}}$ (8)

许多流体流动问题都会受到力的影响，故必须正确处理这些力，满足正确的流体动力学方程。常用的作用力引入形式有三种：速度修正法[15]、离散作用力法[16] [17]和精确差分法[18]。Gong等人[19]对这三种作用力引入形式进行比较，发现采用精确差分法[18]引入作用力可以保证数值计算结果具有更高的精度和稳定性。通过在平衡态密度分布函数$f_i^{eq}\left(x,t\right)$ 中引入修正平衡态速度$u+\Delta u$ 得到式(1)中的体积力项$\Delta f_i\left(x,t\right)$ ，可表示为：

$\Delta f_i\left(x,t\right)=f_i^{eq}\left(\rho \left(x,t\right),u+\Delta u\right)-f_i^{eq}\left(\rho \left(x,t\right),u\right)$ (9)

式中$\Delta u=F\Delta t/\rho$ 表示在单位时间步长$\Delta t$ 内作用在流体粒子上的作用力$F$ [20]所导致流体粒子速度的变化量。上述宏观速度$u$ 的计算公式是在未加作用力$F$ 的情况下得到的，在引入作用力$F$ 后，流体粒子的宏观速度会发生变化，定义流体粒子碰撞前后速度的平均值为流体粒子的实际宏观速度$U$ ，即

$U=u+\frac{\Delta t}{2\rho }F$ (10)

其中作用力$F$ 主要由流体粒子之间的相互作用力$F_{int}$ 和流体粒子与固体粒子之间的相互作用力$F_s$ 两部分组成，即：

$F=F_{int}+F_s$ (11)

其中流体粒子仅考虑与之相邻格子上流体粒子间的相互作用力$F_{int}$ ，并采用Gong等人[19]提出的粒子间相互作用力形式，可表示为：

$F_{int}\left(x\right)=-\beta \psi \left(x\right){\displaystyle \sum _{x\text{'}}G\left(x,x^{\prime }\right)\psi \left(x^{\prime }\right)\left(x^{\prime }-x\right)}-\frac{1-\beta }{2}{\displaystyle \sum _{x\text{'}}G\left(x,x^{\prime }\right){\psi }^2\left(x^{\prime }\right)\left(x^{\prime }-x\right)}$ (12)

式中$\beta$ 表示权重因子，其大小与所选用的状态方程有关。$\psi \left(x\right)$ 表示Shan-Chen伪势模型中的“有效质量”[15]；$G\left(x,x^{\prime }\right)$ 表示为Green函数。其中$\psi \left(x\right)$ 是粒子宏观密度的函数，与所选用的状态方程有关，可表示为：

$\phi \left(\rho \right)=\sqrt{\frac{2\left(p-\rho c_s^2\right)}{c_{0}g}}$ (13)

式中$c_0$ 的大小由格子结构决定，对于本文所采用的D2Q9格子模型，$c_0=6.0$ 。$g$ 只是一个参数，其唯一要求是确保根号内的整个项为正数，一般取$g=-1.0$ 。$p$ 表示流体压力，由给定的状态方程通过计算获得，本文采用Peng-Robinson (P-R)状态方程[21]，可表示为：

$p=\frac{\rho RT}{1-b\rho }-\frac{a{\rho }^2\epsilon \left(T\right)}{1+2b\rho -b^2{\rho }^2}$ (14)

$a=0.45724\frac{R^2{T}_c^2}{p_c}$ (15)

$c_0=6.0$ $b=0.0778\frac{R{T}_c}{p_c}$ (16)

$\epsilon \left(T\right)={\left[1+\left(0.37464+1.54226\omega -0.26882{\omega }^2\right)\left(1-\sqrt{\frac{T}{T_c}}\right)\right]}^2$ (17)

在上述公式中，$T_c$ 和$p_c$ 分别表示临界温度和临界压力；$\omega$ 表示偏心因子，本文选用水的偏心因子，故$\omega =0.344$ 。其他参数分别设置如下：$a=2/49$ ，$b=1/21$ ，$R=1$ [21]。通过式(15)和式(16)计算可得$T_c=0.0729$ 。式(12)中的$G\left(x,x\text{'}\right)$ 满足$G\left(x,x^{\prime }\right)=G\left(x^{\prime },x\right)$ ，反映了相邻流体粒子之间的相互作用强度，其可表示为[22]：

$G\left(x,x^{\prime }\right)=\{\begin{array}{ll}{g}_1,\hfill & \left|x-x^{\prime }\right|=1\hfill \\ g_2,\hfill & \left|x-x^{\prime }\right|=\sqrt{2}\hfill \\ 0,\hfill & \text{otherwise}\hfill \end{array}$ (18)

其取值由格子结构决定，对于所采用的D2Q9格子模型而言，$g_1=2g$ ，$g_2=g/2$ 。流-固相互作用力$F_s$ 采用Li等人[23]提出的改善虚拟密度格式，可表示为：

$F_s\left(x\right)=-G\psi \left(x\right){\displaystyle \sum _i{w}_i}\psi \left({\rho }_w\left(x\right)\right)s\left(x+e_i\right)e_i$ (19)

${\rho }_{ave}\left(x\right)=\frac{{\displaystyle \sum _i{w}_i\rho \left(x+e_i\right)}\left[1-s\left(x+e_i\right)\right]}{{\displaystyle \sum _i{w}_i\left[1-s\left(x+e_i\right)\right]}}$ (20)

${\rho }_w\left(x\right)=\{\begin{array}{ll}\phi {\rho }_{ave}\left(x\right),\hfill & \phi \ge 1\hfill \\ {\rho }_{ave}\left(x\right)-\Delta \rho ,\hfill & \Delta \rho \ge 0\hfill \end{array}$ (21)

为了获得更好的力离散形式来计算流–固相互作用力$F_s$ ，公式(19)可根据公式(12)进行修改，表示为：

$F_s\left(x\right)=-\beta \psi \left(x\right){\displaystyle \sum _{i}G\left(x,x+e_i\right)\psi \left({\rho }_w\right)e_i}-\frac{1-\beta }{2}{\displaystyle \sum _{i}G\left(x,x+e_i\right){\psi }^2\left({\rho }_w\right)e_i}$ (22)

式(19)中$s$ 表示指示函数，与周围格子上是流体粒子还是固体粒子有关。当位置$x+e_i\Delta t$ 处是流体粒子时，$s\left(x+e_i\Delta t\right)=0$ ；当位置$x+e_i\Delta t$ 处是固体粒子时，$s\left(x+e_i\Delta t\right)=1$ 。式(20)计算的是固体粒子的局部加权平均密度${\rho }_{ave}\left(x\right)$ ，式(21)计算的是固体粒子的虚拟密度${\rho }_w\left(x\right)$ ，通过调整$\phi$ 和$\Delta \rho$ 的值改变固体粒子的虚拟密度${\rho }_w\left(x\right)$ 进而调整固体壁面的接触角$\theta$ 。当固体粒子的虚拟密度${\rho }_w\left(x\right)={\rho }_{ave}\left(x\right)$ 时，相应的接触角$\theta ={90}^{\circ }$ [23]。

2.2. 问题描述

图1为液桥自然蒸发的示意图。将两个宽度为$D$ 的固体基板静止放置在计算域为$NX\times NY=300\times 300$ 的空间中，并在两基板之间初始化高度为$L$ 的液桥。具体模拟参数设置如下，系统温度设置为$T=0.9T_c$ ，根据Maxwell 构造的汽液共存曲线[24]，此系统温度下对应的饱和液相密度和饱和蒸汽密度分别为${\rho }_l=5.91$ 和${\rho }_v=0.58$ 。取松弛时间$\tau =1.0$ ，对应的液相运动粘度和汽相运动粘度分别为${\upsilon }_l=0.17$ 和${\upsilon }_v=1.0$ 。本文采用P-R状态方程，对应的权重因子$\beta =1.16$ 。固体表面的接触角设置为$\theta ={60}^{\circ }$ 。另外，将计算域的左右两侧设置为周期性边界条件，顶部设置为定压边界条件，底部的平坦固体壁面设置为标准半反弹边界。通过调整计算域顶部定压边界的压力大小，可以改变整个系统的相对湿度$RH$ ，进而分别模拟不同相对湿度$RH\left(85\%,90\%,95\%\right)$ 、不同固体基板宽度$D\left(80,100,120\right)$ 和不同液桥高度$L\left(30,50,70\right)$ 的条件下，液桥自然蒸发过程。通过分析对应条件下液桥颈部宽度$d_1$ ，腰部宽度$d_2$ ，液桥与上基板之间的接触角$\theta$ ，液桥的质量$m$ 随时间$t$ 的变化趋势、对应条件下液桥断裂时间$t_b$ 长短，进而得到不同条件下液桥自然蒸发的动力学规律。

Figure 1. Schematic diagram of natural evaporation of a liquid bridge

图1. 液桥自然蒸发示意图

3. 结果与分析

Figure 2. (a) is the evaporation process diagram of the liquid bridge from $t=0$ to $t=480000$ under the conditions of $RH=95\%,D=100,L=50$ , and (b) is the trend of the contact angle $\theta$ with time $t$

图2. (a)为液桥在$RH=95\%,D=100,L=50$ 的条件下从$t=0$ 到$t=480000$ 时刻的蒸发过程图，(b)为接触角$\theta$ 随时间$t$ 的变化趋势图

图2(a)为液桥在$RH=95\%,D=100,L=50$ 的条件下从$t=0$ 到$t=480000$ 时刻的蒸发过程图，图2(b)为在$RH=95\%,D=100,L=50$ 的条件下接触角$\theta$ 随时间$t$ 的变化趋势图。从图2(a)可以看出，在$t=0$ 时刻，液桥近似为一个规则的矩形形状，测得此时$d_1=d_2=101$ ，接触角$\theta ={90}^{\circ }$ 。在$20000\le t\le 50000$ 时刻内液桥腰部位置$d_2$ 逐渐减小为85，但是$d_1$ 几乎不变，接触角$\theta$ 从${79}^{\circ }$ 降低到${66}^{\circ }$ 。这是因为基板边缘对液桥接触线的钉扎作用，使得$d_1$ 在蒸发初期不发生变化。在这个过程中液桥以固定接触线的模式蒸发。液桥腰部的不断蒸发，使得$d_2$ 不断减小并使得接触角$\theta$ 相应减小。当接触角足够小时($t=60000,\theta ={60}^{\circ }$ )，在表面张力作用下液桥接触线在基板边缘脱钉，因此液桥的颈部直径($d_1$ )开始收缩，液桥以固定接触角的模式蒸发。从图2(b)也可以看出$\theta$ 随着$t$ 的增大先减小后稳定在${60}^{\circ }$ 附近。我们把接触角$\theta$ 变为恒定值的时刻(颈部宽度$d_1$ 开始改变的时刻)称为临界时刻$t_c$ 。在$RH=95\%,D=100,L=50$ 的条件下，$t_c=60000$ 。从图2(a)还可以发现，在$t=470000$ 时，液桥的腰部出现断裂，我们把液桥从腰部断裂的时刻称为断裂时刻$t_b$ 。在$RH=95\%,D=100,L=50$ 的条件下，$t_b=470000$ 。

3.1. 相对湿度对液桥蒸发的影响

Figure 3. (a) The variation of $d_1$ with time $t$ during the evaporation process of the liquid bridge under different humidity $RH$ , (b) is the enlarged image in the red box area in (a), and (c) the trend of the contact angle $\theta$ of the liquid bridge with time $t$ under different humidity $RH$

图3. (a)不同湿度$RH$ 下液桥蒸发过程中$d_1$ 随时间$t$ 的变化，(b)为(a)中红框区域内的放大图，(c) 不同湿度$RH$ 下的液桥接触角$\theta$ 随时间$t$ 的变化趋势

不同湿度$RH$ 下液桥颈部的宽度$d_1$ 随时间$t$ 变化趋势如图3(a)、图3(b)所示，其中图3(b)为图3(a)中红框区域的放大图。从图3(a)可以发现，湿度$RH$ 越大，单位时间内液桥颈部宽度$d_1$ 的变化量越小，即蒸发速度越小。这是因为湿度$RH$ 增高导致空气与液体间的水汽压差降低，减缓了水分子的蒸发速率。并且从图3(b)可以发现，不同湿度$RH$ 下颈部宽度$d_1$ 开始减小的临界时刻$t_c$ 不同： $t_c\left(RH=95\%\right)=60000>t_c\left(RH=90\%\right)=14000>t_c\left(RH=85\%\right)=8000$ 。图3(c)为不同湿度$RH$ 下接触角$\theta$ 随时间$t$ 的变化趋势。从图中可以看到，随着时间$t$ 的增大，不同湿度$RH$ 下接触角$\theta$ 随时间$t$ 的变化趋势为：$\theta$ 随着$t$ 先减小，而后稳定在${60}^{\circ }$ 附近。且湿度$RH$ 越大，液桥接触线在基板边缘处脱钉所需时间越长。

Figure 4. (a) is the evaporation velocity field of the liquid bridge when $t=20000,RH=85\%,D=100,L=50$ , (b) is the evaporation velocity field of the liquid bridge when $t=20000,RH=95\%,D=100,L=50$ , (c) is the variation trend of the waist width $d_2$ of the liquid bridge with time $t$ under different humidity $RH$ , and (d) is the trend of the total mass $m$ of the liquid bridge with time $t$ under different humidity $RH$

图4. (a)为$t=20000,RH=85\%,D=100,L=50$ 时液桥的蒸发速度场，(b)为$t=20000,RH=95\%,D=100,L=50$ 时液桥的蒸发速度场，(c)为不同湿度$RH$ 下液桥腰部宽度$d_2$ 随时间$t$ 的变化趋势图，(d)为不同湿度$RH$ 下液桥的总质量$m$ 随时间$t$ 的变化趋势图

图4(a)、图4(b)为$t=20000$ 时，$D=100,L=50$ 条件下的液桥在不同湿度下($RH=85\%,RH=95\%$ )的蒸发速度场矢量图。箭头方向为蒸发的方向，箭头的长短表示蒸发速度的大小。从图4(a)和图4(b)中的红框区域可以发现，在蒸发的过程中与液桥颈部的速度场相比，液桥腰部表示速度大小的箭头更长，这说明液桥不同位置处的蒸发速度不均匀，且液桥腰部的蒸发速度更大。与高湿度下的液桥(图4(b))相比，低湿度下(图4(a))中的箭头更长，意味着低湿度下液桥在腰部蒸发更快。这是因为湿度$RH$ 较低时，空气中水汽分压低，与液体表面的饱和水汽压之间的差距较大，因此蒸发更快。这导致液桥腰部的水分子更快地向空气中扩散，从而产生更强的局部空气流动。且由于液桥腰部界面弯曲效应，导致液桥在腰部蒸发速度比其他位置处更快。图4(c)为不同湿度$RH$ 下液桥的腰部宽度$d_2$ 随时间$t$ 的变化趋势图。从图4(c)可以看出，不同湿度$RH$ 下，液桥断裂时刻($t_b$ )的大小关系为： $t_b\left(RH=95\%\right)=470000>t_b\left(RH=90\%\right)=130000>t_b\left(RH=85\%\right)=80000$ 。这说明湿度$RH$ 越大，液桥从开始蒸发到其断裂所需时间越长。这是因为湿度越高，液桥腰部蒸发越慢。图4(d)为不同湿度$RH$ 下液桥的质量$m$ 随时间$t$ 的变化趋势图。从图4(d)可以看出，湿度$RH$ 越高，液桥质量$m$ 变为0的时间越长，即湿度$RH$ 越大液桥整体的平均蒸发速度也越小。

3.2. 液桥尺寸改变对液桥蒸发的影响

3.2.1. 基板宽度对液桥蒸发的影响

Figure 5. (a) is the variation of $d_1$ with time $t$ during the evaporation process of the liquid bridge of different $D$ , (b) is the enlarged view of the red box area in (a), and (c) is the variation trend of the contact angle $\theta$ with time$t$

图5. (a)为$D$ 不同时液桥蒸发过程中$d_1$ 随时间$t$ 的变化，(b)为(a)图中红框区域的放大图，(c)为$D$ 不同时的接触角$\theta$ 随时间$t$ 的变化趋势图

图5(a)和图5(b)为改变基板宽度$D\left(80,100,120\right)$ 后得到的液桥颈部宽度$d_1$ 随时间$t$ 的变化趋势图，其中图5(b)为图5(a)中红框内的放大图。从图5(a)可以看出，改变基板的宽度$D$ 不会影响单位时间内液桥$d_1$ 的减小量。这说明基板宽度$D$ 不会影响液桥的蒸发速度。但是基板宽度$D$ 越大，同一时刻$t$ 下液桥的$d_1$ 越大。从图5(b)中可以发现，不同基板宽度$D$ 下颈部宽度$d_1$ 开始改变的时间($t_c$ )为：$t_c\left(D=120\right)=14000=t_c\left(D=100\right)=14000\approx t_c\left(D=80\right)=12000$ 。这说明基板宽度$D$ 的改变不影响$t_c$ 的大小。图5(c)为不同基板宽度$D$ 的液桥接触角$\theta$ 随时间$t$ 的变化趋势。从图5(c)中可以发现接触角$\theta$ 随时间$t$ 先减小到${60}^{\circ }$ 附近后稳定在这一值。且基板宽度$D$ 不影响临界时刻$t_c$ ，意味着液桥接触线在基板边缘处脱钉所需时间相同。

Figure 6. (a) is the evaporation velocity field of the liquid bridge with different substrate widths $D$ at $t=10000$ and $t=50000$ , (b) the variation of the waist width $d_2$ of the liquid bridge with different substrate widths $D$ with time $t$ , and (c) the variation trend of the liquid bridge mass $m$ with with time $t$ at different substrate widths $D$

图6. (a)为基板宽度$D$ 不同的液桥在$t=10000$ 和$t=50000$ 时的蒸发速度场，(b) 基板宽度$D$ 不同的液桥的腰部宽度$d_2$ 随时间$t$ 的变化图，(c)为基板宽度$D$ 不同的液桥质量$m$ 随时间$t$ 的变化趋势图

图6(a)为不同基板宽度$D$ 的液桥在$t=10000$ 和$t=50000$ 时的蒸发速度场。从图中红框区域可以看出，基板宽度$D$ 不同，模拟得到的液桥不同位置处的蒸发速度场箭头长短在同一时刻是一致的。这说明基板宽度$D$ 对不同位置处的蒸发速度的影响可以忽略。图6(b)为基板宽度$D$ 不同时，液桥的腰部宽度$d_2$ 随时间$t$ 的变化趋势图。从图6(b)可以看出，腰部宽度$d_2$ 减小到0所需时长不同，液桥断裂时刻$t_b$ 的关系为：$t_b\left(D=120\right)=170000>t_b\left(D=100\right)=130000>t_b\left(D=80\right)=100000$ 。这说明增大基板宽度$D$ ，液桥的断裂时刻$t_b$ 增大。图6(c)为液桥的质量$m$ 随时间$t$ 的变化示意图。如图6(c)所示，增大基板宽度$D$ 提高了液桥的初始质量，但单位时间内质量$m$ 的变化量相同，说明不同基板宽度$D$ 下液桥整体平均蒸发速度不变。

3.2.2. 液桥高度对液桥蒸发的影响

图7(a)、图7(b)为改变液桥高度$L\left(30,50,70\right)$ 后得到液桥颈部宽度$d_1$ 随时间$t$ 的变化趋势图，其中图7(b)为图7(a)中红框内的放大图。从图7(a)可以看出，高度$L$ 越大的液桥颈部宽度$d_1$ 在单位时间内的减小量越大。从图7(b)中可以发现，不同高度$L$ 下颈部宽度$d_1$ 开始改变的时刻$t_c$ 的相对大小关系为：$t_c\left(L=70\right)=12000<t_c\left(L=50\right)=140000<t_c\left(L=30\right)=20000$ 。这意味着增大高度$L$ 会减小临界时刻$t_c$ 。图7(c)为不同高度$L$ 的液桥接触角$\theta$ 随时间$t$ 的变化趋势图。从图7(c)中可以发现接触角$\theta$ 随时间$t$ 先减小后稳定在${60}^{\circ }$ 附近。且液桥高度$L$ 越大，液桥接触线在基板边缘处脱钉所需时间越短。

Figure 7. (a) is the change of $d_1$ with time $t$ during the evaporation process of liquid bridges with different heights $L$ , (b) is a magnified view of the red box area in (a), and (c) is the change of the contact angle $\theta$ of liquid bridges with time $t$ at different heights $L$

图7. (a)为高度$L$ 不同的液桥在蒸发过程中$d_1$ 随时间$t$ 的变化，(b)为(a)图中红框区域的放大图，(c)为高度$L$ 不同的液桥接触角$\theta$ 随时间$t$ 的关系

图8(a)为不同高度$L$ 的液桥在$t=50000$ 时的蒸发速度场。从图8(a)中红框区域可以看出，随着高度$L$ 的增大，液桥腰部蒸发速度场的箭头变长。这说明高度$L$ 的增大使得液桥腰部的蒸发速度增大。图8(b)为不同高度$L$ 的液桥的腰部宽度$d_2$ 随时间$t$ 的变化趋势图，可以看出不同液桥高度$L$ 下液桥的断裂时刻$t_b$ 的大小关系为：$t_b\left(L=30\right)=220000>t_b\left(L=50\right)=130000>t_b\left(L=70\right)=110000$ 。这是由于$L$ 越大的液桥，其腰部蒸发速度越快。图8(c)为不同高度$L$ 的液桥的质量$m$ 随时间$t$ 的变化趋势图。从图8(c)可以发现，液桥高度$L$ 越大液桥的初始质量也越大，且液桥整体的平均蒸发速度也越大。

Figure 8. (a) is the evaporation velocity field of the liquid bridge with different heights $L$ at $t=50000$ , (b) the variation of the waist width $d_2$ of the liquid bridge with time $t$ at different heights $L$ , and (c) the trend of the mass $m$ of the liquid bridge with time $t$ at different heights $L$

图8. (a)为高度$L$ 不同的液桥在$t=50000$ 时的蒸发速度场，(b)高度$L$ 不同的液桥的腰部宽度$d_2$ 随时间$t$ 的变化图，(c)为高度$L$ 不同的液桥质量$m$ 随时间$t$ 的变化趋势图

4. 结论

本文采用格子Boltzmann方法(LBM)研究了液桥蒸发过程中环境湿度$RH$ 、液桥宽度$D$ 和高度$L$ 对液桥蒸发动力学的影响。得到的主要结论如下：

(1) 增大湿度$RH$ 可以降低液桥腰部的蒸发速度，液桥蒸发的整体平均速度也会减小，使得液桥颈部宽度$d_1$ ，腰部宽度$d_2$ 和液桥质量$m$ 在单位时间内的变化率降低，液桥的断裂时刻$t_b$ 更大。也会延长接触角$\theta$ 变为恒定值的时刻$t_c$ ，增大了液桥接触线在基板边缘处脱钉所需的时间。

(2) 增大基板宽度$D$ 不会改变液桥腰部的蒸发速度和液桥蒸发的整体平均速度。但是基板宽度$D$ 越大的液桥初始颈部宽度$d_1$ ，腰部宽度$d_2$ 以及初始质量$m$ 更大，导致液桥的断裂时刻$t_b$ 越大。增大基板宽度$D$ 不会影响接触角$\theta$ 变为恒定值的时刻$t_c$ ，因此也不影响液桥接触线在基板边缘处脱钉所需的时间。

(3) 增大液桥高度$L$ 会使得液桥腰部的蒸发速度增大，液桥蒸发时的整体平均速度也越大，使得液桥颈部宽度$d_1$ ，腰部宽度$d_2$ 和液桥质量$m$ 在单位时间内的变化率增大，使液桥的断裂时刻$t_b$ 越小。增大高度$L$ 可以减小接触角$\theta$ 变为恒定值的时刻$t_c$ ，减小液桥接触线在基板边缘处脱钉所需的时间。

NOTES

^*通讯作者。