VirtualFlow | 基于粗粒连续接触角模型的微流控装置充注仿真模拟

近年来，微流控技术发展迅速。然而，在建模和仿真工具方面很难有成熟的设计支持。其核心的难点在于需要对传统计算流体动力学(CFD)工具进行许多扩展，例如表面力、三相接触线动力学、非均相化学反应、表面活性剂等小规模的物理效应等。

具有移动接触线的不可压缩两相流在各种应用中很常见，例如涂层和生物过程。模拟这种流动的困难在于两种流体的N-S方程结合无滑移边界条件会导致接触线上存在剪切应力奇点。此外，还需要指定自由表面与固体表面形成的接触角。

在过去几十年里，人们曾多次尝试在连续流体动力学框架内对接触线动力学进行建模。然而，这些模型几乎很少考虑到作用于三线区域的作用力。VirtualFlow中提供的模型是基于接触线的分子动力学建立的，从而不需要给定接触角。在基于杨氏定律的动量方程中加入不平衡杨氏力，而非指定接触角。采用滑移长度对接触线上的应力奇点进行正则化，并对直至滑移长度的对数奇点进行积分，从而解决这一问题。

将模型与Spelt提出的单元模拟进行了比较，并使用其分析具有不同表面特性的液体对原型微器件的填充过程。气泡夹带是微流控器件设计中的一个常见问题。通过精细设计器件的几何形状、控制填充过程和选择材料特性(疏水或亲水性)，可以避免这些气泡的形成。在这一过程中，可以使用可靠的仿真工具对几何形状进行设计修改，以确保设计是有效可行的。

其中动量方程(方程(2))中的RHS项分别为体力、表面张力（方程(4)）及其接触线壁面贡献。不考虑相变时，σij为牛顿应力张量。在level-set方法中，两种流体之间的界面用连续函数φ表示，φ表示到界面的距离，该距离在一侧为正，在另一侧为负。通过这种方式，可以识别两种流体，并且将物理自由界面的位置与0等值面相关联。Level-set控制方程由方程3给出：

其中γ为表面张力，kappa为界面曲率，δ^I为以界面为中心的光滑狄拉克函数。方程在一次迭代之后，level-set函数不再是与界面的有符号距离。为了恢复其在界面附近的正确分布，需要重新赋值距离，将下式积分为稳态:

其中sgn(x)是sgum函数。在VirtualFlow中方程(6)在方程(3)的每个迭代步后，使用非振荡(WENO)三阶格式求解。

 VirtualFlow采用与三线相关物理力的润湿动力学数值处理方法，可以提供对润湿动力学的物理充分描述，从而消除了指定接触角的任何特定边界条件的需要。在本研究中使用的三线力是基于考虑界面自由能的。因此，它只包含两个参数:流体之间的界面张力γ和平衡接触角   θ_eq。

其中，θ_dy是瞬时动态接触角，δ_t是狄拉克函数。

         

松弛到平衡          
         

在这里，我们考虑Spelt（SpeltP. D. M., J. Comp. Phys., 207, 389-404 (2005)）提出的两种情况，即液滴浸入不同的流体中，并附着在矩形域(2 × 1，由2N × N离散)的长边界上。在这些情况下使用的流体性质列于表1中。在Case1中，密度相等(= 1)，而粘度不同(液滴为4.95×10−2，周围的相为4.95 × 10−3)。然而，在case 2中液滴与周围流体的密度比为1:20，粘度恒定为4.95 × 10−2。两个Case的表面张力分别为γ = 0.11和2.21，网格尺寸为128x64。

(a)t=0s

(b)t=3.7s

图1:非润湿液体恢复平衡;θ_eq = 1200

首先，液滴呈现圆帽形(半径为1.66)，图1和图3中的蓝色相的接触角分别为300度和150度，这与平衡接触角θ_eq=120度和600度不同。当前模型与Spelt模型的主要区别在于，后一种情况需要额外规定接触线速度。

(c)t=1.0s

(d)t=3.6s

图1和图3分别显示了case 1和2情况下平衡时的状态。接触角在流动趋于平衡时的变化如图2所示。对于case 2，由于较高的表面张力以及由于较高的外相密度而导致其可以更快地恢复平衡。

为了说明VirtualFlow在实际应用中的实用性，对θ_eq=70°、90°和110°三种不同液体的微型泵原型装置的充液过程进行了仿真。图4显示了最终结果，该结果表明，对于疏水性情况，必须特别小心，才能设计出没有气泡截留的设计。

(b) 接触角=70°

(c) 接触角=90°

(d) 接触角=110°

图4:微器件原型液充

图中蓝色代表液体，灰色代表空气，实线表示相间界面，白色表示固体区域