胡晟, 陈悦江, 吕晓永, 吴东旭
东北大学秦皇岛分校 控制工程学院, 河北 秦皇岛 066004
收稿日期：2020-03-01
基金项目：国家自然科学基金资助项目(61903069，51705070)；河北省自然科学基金资助项目(E2018501038)。
作者简介：胡晟(1984-)，男，云南景洪人，东北大学秦皇岛分校副教授，博士。

摘要：基于Cahn-Hilliard建立Navier-Stokes两相流体动力学和电场Maxwell应力张量法的多物理场耦合模型, 用于平行极板型、针型和圆环型电极的液滴聚结仿真.研究结果表明: 较强的电场强度诱导分散相液滴聚结耗时较短, 液滴链的结构和均匀度受电极形状和电场空间均匀性的影响较大; 均匀电场诱导液滴成链较为均匀, 并不受液滴数量的影响.在指针型和圆环型电极产生的非均匀电场只能在液滴数量较少的条件下, 实现规则液滴链的生成.通过多物理场模型建模, 仿真结果能够为静电纺丝、液滴合并、液泡回收等复杂微流体电学控制提供理论基础.
关键词：介电泳电聚结液滴Cahn-Hilliard方程COMSOL
Process and Distribution of Droplet Coalescence Based on Dielectrophoresis Effect
HU Sheng, CHEN Yue-jiang, LYU Xiao-yong, WU Dong-xu
School of Control Engineering, Northeastern University at Qinhuangdao, Qinhuangdao 066004, China
Corresponding author: HU Sheng, E-mail: husheng@neuq.edu.cn.

Abstract: Based on Cahn-Hilliard, a multiphysics coupling model of Navier-Stokes two-phase fluid dynamics and electric field Maxwell stress tensor method is established, which is used for the simulation of droplet coalescence of parallel plate, pointer and ring electrodes. The research results show that the stronger electric field intensity induces a shorter time for the coalescence of the dispersed phase droplets, and the structure and uniformity of the droplet chain are greatly affected by the electrode shape and the spatial uniformity of the electric field. The uniform electric field induces the droplets to form a more uniform chain and is not affected by the number of droplets. The non-uniform electric field generated by the pointer-type and ring-type electrodes can only realize the generation of regular droplet chains under the condition of a small number of droplets. The multi-physics model and simulation results involved in this article can provide a theoretical basis for complex microfluidic electrical control such as electrospinning, droplet merging, and bubble recovery.
Key words: dielectrophoresiselectrocoalescencedropletCahn-Hilliard equationCOMSOL
借助微流控芯片进行不同粒径、形状与数量的液滴(10^-9~10^-18L)连续快速生成, 能够满足生化物质信息提取处于高检测精度、快速反应和低样本消耗量的使用标准^[1].目前的液滴生成方法一般分为被动式与主动式^[2].被动式生成法主要通过调节输入通道的连续相与离散相流速比率来实现, 主动式生成法多采用外加物理场, 如: 电场、磁场、光场、热力场和离心旋转等.通过改变两相流中的黏度比、惯性力或毛细力的动态平衡, 使液相界面之间达到不稳定状态而生成预期的液珠形状.电场是目前使用最为广泛、研究较为成熟的技术手段.内置金属电极产生诸如电泳、介电泳、电渗流、电聚结、电穿孔等微观效应.其中, 介电泳效应诱导两相介质水-油聚结, 产生水-油-水双层液滴防止生化试剂的交叉污染^[3].对水-油化合物的文献报道, 国内已有****进行了均匀/非均匀电场下液滴聚结与电场关键参量的研究^[4], 进一步揭示了交流电场作用下单/双液滴的形变过程^[5].
上述研究主要集中于毫米级大尺寸液滴聚结行为的理论仿真和分析, 对微米级小尺寸液滴研究不足, 特别对大量无序分布的液滴聚结过程和拓扑结构研究甚少.微液滴封装单细胞能够较好进行RNA测序或蛋白质组成分分析.为了进行液滴回收, 以及细胞排泄外泌体的二次检测, 都需要准确从油相化合物中提取水相物质.本文采用Cahn-Hilliard方程进行油-水相的动态模拟, 分析水相液滴的形态变化、珠链过程和拓扑结构.针对当前主要的电极结构进行相关的电场、流场、液相场的时域仿真计算.理论模型和计算结果有助于预测水相聚结时间, 为油水分离提供理论指导.
1 理论与模型1.1 原理介绍芯片腔体受到激励电极时, 满足Poisson方程:

(1)

式中: ρ_v为体电荷密度, C/m³; D为电位移量(D=εE, ε=ε_rε₀为介电常数, εr为相对介电常数, ε₀为真空介电常数, F/m), C/m².对电中性介质而言, 其净电荷为零.因此Poisson方程简化为Laplace方程▽·D=0.
根据电场强度E(V/m)和麦克斯韦应力张量(Maxwell stress tensor)法^[6], 可得空间单位体积液滴承受的介电泳力:

(2)

式中，T为应力张量, N/m³.

(3)

式中，I为单位矩阵.
Cahn-Hilliard方程是一种数学物理方程, 描述不同相场之间的分离或聚合.相场变量φ反应二元流体分布函数, 并在每个相中形成纯净的畴.根据化学势的描述, 通过迭代演化来模拟相场的时域变化^[7]:

(4)

式中: u为流体速度, m/s; d_pf为两相介质之间的界面厚度, m, 主要与网格稀疏程度相关; σ为表面张力系数, N/m; χ为两相介质的迁移系数, (m·s)/kg; ψ为中间变量:

(5)

油-水两相的体积分数分别为

(6)

η为混合物质的黏度, Pa·s; ρ为密度, kg/m³; ε为介电常数, F/m, 表达式如下:

(7)

1.2 电极模型由式(1)和式(2)可知, 电场受金属电极位置和排列的影响, 主要对3种电极结构的液滴进行仿真模拟, 如图 1所示.图 1a左/右金属极板长度相同, 类似平行电容极板, 可产生均匀电场, 后两种电极形成的电场为非均匀电场.求解Laplace方程后, 使电场与Navier-Stokes流场耦合求解偏微分方程:

(8)

图 1(Fig. 1)

图 1 均匀与非均匀电场的液滴聚结模型Fig.1 Droplet coalescence model of uniform and non-uniform electric fields (a)—平行极板型; (b)—针型; (c)—环型.

式中: p为压强, Pa; F_st为界面张力, N/m³, 如式(9)所示.

(9)

本文的离散相(水相)和连续相(油相)的仿真参数如表 1所示.
表 1(Table 1)

表 1 油-水相物质的关键系数Table 1 Key parameters of oil-water phases 变量 水相(w) 油相(o) ρ /(kg·m-3) 997 885 η /(Pa·s) 10-3 0.474 εrε0/(F·m-1) 80ε0 2.5ε0 σ/(N·m-1) 0.005

表 1 油-水相物质的关键系数 Table 1 Key parameters of oil-water phases

2 仿真与分析2.1 均匀电场根据图 1a的几何描述, 左边界施加300 V直流电压, 右边界接地.矩形宽度和高度分别为500, 800 μm.液滴半径为25 μm.根据式(1), 式(4)和式(8), 采用COMSOL Multiphysics 5.3a有限元软件进行电场、流场和相场的偏微分求解, 可得图 2中20个液滴随时间的聚结变化.为了贴近实验, 模拟液滴位置满足随机分布.由图 2a的A及B虚线框可知, 3个液滴之间的欧拉间距较其他液滴越小, 它们相互聚结现象越明显.距离更近的2个液滴开始相互铰链融合, 如图 2b中虚线框A′和B′所示.尺寸较大的液滴又与周围较近的液滴完成新的聚结, 如图 2c所示.随着时间的增加, 部分液滴因为空间距离较远, 无法完成进一步聚结, 如图 2d所示.
图 2(Fig. 2)

图 2 300 V时液滴随时间的变化(白色箭头代表电场强度和方向)Fig.2 Change of droplet with time at 300V (a)—0 s; (b)—0.002 s; (c)—0.04 s; (d)—0.1 s.

施加300 V电压时, 剩余14个液滴依然处于较高的分散相状态, 主要原因是诱导液滴产生的电偶极矩强度较弱.将激励电压分别提高到600 V和900 V, 模拟结果如图 3、图 4所示.此时液滴体系相互吸引, 并且聚合的数量更多, 所需时间较短.聚合的液滴开始出现长链结构, 且沿着白色电场线的方向延伸.对比图 3和图 4发现, 电场强度越高, 液滴成链越明显, 同时液滴链更与电场传播方向平行.但长链彼此相互排斥, 不能形成更粗的链型结构.
图 3(Fig. 3)

图 3 600 V时液滴随时间的变化Fig.3 Change of droplet with time at 600 V (a)—0.02 s; (b)—0.1 s; (c)—0.15 s; (d)—0.2 s.

图 4(Fig. 4)

图 4 900 V时液滴随时间的变化Fig.4 Change of droplet with time at 900 V (a)—0.02 s; (b)—0.1 s; (c)—0.15 s; (d)—0.2 s.

2.2 非均匀电场针型电极广泛用于细胞和胶体颗粒的俘获和富集, 本文建立的上、下针型电极模型尺寸相同, 两者间距L_e为500 μm, 宽度W_e为20 μm, 针尖曲率半径R_e为10 μm.由平行极板仿真可知300 V激励过低, 600 V和900 V的激励电压结果如图 5、图 6所示.
图 5(Fig. 5)

图 5 电压为600 V时针型电极液滴形态Fig.5 Pointer-type electrode droplet shape when the voltage is 600 V (a)—0 s; (b)—0.2 s.

图 6(Fig. 6)

图 6 电压为900 V时指针型电极液滴形态Fig.6 Pointer-type electrode droplet shape when the voltage is 900 V (a)—0 s; (b)—0.2 s.

上述计算结果表明: 900 V电场强度相比600 V电压形成的电场强度较大, 液滴受正介电泳力作用而迁移至电极两端.液滴从圆球演变为椭圆形, 该形状与电场方向相同.液滴之间相互聚结逐渐实现上、下电极的桥接.与均匀电场不同, 液滴链仅为单链结构.保证激励电压为900 V时, 降低上、下电极间距L_e为300 μm, 液滴成链的空间均匀性降低, 如图 7所示.液链保持两端粗、中间细的形态分布.不规则短链的形成主要由于电场间距较近, 电极两端的电场强度较高, 液滴更倾向于移动到电极两端聚结成粗状链条, 但中间部分的电场强度较弱, 聚结后的液滴链条相对较细.对环型电极的液滴聚结进行研究, 圆环半径R₁和R₂分别为200, 500 μm.液滴分布与上述两种激励模型相同, 保持互不重合的随机分布, 如图 8a所示.
图 7(Fig. 7)

图 7 900 V时上、下电极间距300 μm的聚结液滴分布Fig.7 Distribution of coalescing droplets with distance of 300 μm between top and bottom electrodes at 900 V (a)—0.2 s; (b)—0.3 s.

图 8(Fig. 8)

图 8 电压300 V时环型电极液滴聚结结果Fig.8 Result of droplet coalescence of ring-type electrode at 300 V (a)—0 s; (b)—0.05 s; (c)—0.1 s; (d)—0.2 s.

通过施加300 V电压发现液滴逐渐靠近R₁电极, 轮辐形状开始成链状向R₂电极伸展, 如图 8b~8d所示.偏置电压提高到600 V时, 相同位置的20个液滴仿真结果如图 9所示.
图 9(Fig. 9)

图 9 电压600 V时环型电极液滴聚结结果Fig.9 Result of droplet coalescence of ring-type electrode at 600 V

研究发现电场强度提高形成的液滴链增多, 与图 4d相同.在图 9a的A区域, 铰链的液滴并非形成更粗的液柱, 而受介电泳排斥力分离为两条新的液柱, 如图 9b的A′区域所示.仿真电压增至900 V时, 液滴聚结速度加快, 在0.11 s实现了轮辐状液柱的生成, 但液柱的稳定性较差.3种电极结构的仿真表明针型电极的成链数量较少, 且电极间距对液柱的粗细影响较大.上述结论与文献[8-10]一致.
在0.2 s时, 提取平行极板和环型电极边界的φ_w值, 如图 10所示.分析发现均匀电场形成的链形粗细基本均匀, 但环型电极在[R₁, R₂]链由粗变细.该现象主要归因于内环电场强度大于外环电场, 液滴倾向移动至内环电极边缘.
图 10(Fig. 10)

图 10 电极边界的液相数值Fig.10 Liquid phase values at electrode boundary (a)—平行极板; (b)—环型极板.

如图 10a所示, 对长链A根据半高全宽编程得到平行极板形成的液柱宽度为30.91 μm.对图 10b中的B长链计算得到最大宽度为52.31 μm, 靠近R₂电极的链宽度为47.56 μm.为了研究不同数量液滴形成的单一链柱, 本文将不同数量的液滴设定在一个狭小空间中.电极尺寸不变, 液滴在宽度[25 μm, 475 μm], 高度[25 μm, 125 μm]矩形框内随机分布.液滴数量分别从5~10个进行求解, 电压为900 V, 0.2 s内各分散相的液滴均聚结形成如图 4d的定向链柱, 计算宽度发现液滴数量增多, 其形成的液链厚度增加, 如图 11所示.
图 11(Fig. 11)

图 11 液滴数对液柱宽度的影响Fig.11 Effect of droplets number on liquid column width

液滴在环型电极径向[225 μm, 475 μm], 方位角[-100°, -80°]以内求解的仿真结果显示非均匀电场相比平行极板电场更难以成链, 如图 12所示.均匀液相分布的液链受液滴数量和分布的限制.与针型电极相似, 密度过高的离散液滴受非均匀电场的控制更强, 而周围液滴相互作用的次级介电泳效应更弱.
图 12(Fig. 12)

图 12 环型电极内不同数量液滴的相场Fig.12 Phase field of different numbers of droplets in a ring electrode (a)—N=5; (b)—N=6; (c)—N=7; (d)—N=8; (e)—N=9; (f)—N=10.

3 结论1) 与电极结构无关, 较大的电场强度易于分散相液滴的聚结.
2) 相比针型和环型电极产生的非均匀电场, 平行极板产生的均匀电场能够较好实现液滴的聚结成链, 液滴数量越多产生的链形越宽.
3) 针型和环型电极在低液滴数下也能形成液链结构, 电场较弱区域的液链宽度小于电场强度较大区域的链宽.
参考文献

[1]	Feng S L, Shirani E, Inglis D W. Droplets for sampling and transport of chemical signals in biosensing: a review[J]. Biosensors-Basel, 2019, 9(2): 1-14.
[2]	Zhu P A, Wang L Q. Passive and active droplet generation with microfluidics: a review[J]. Lab on a Chip, 2017, 17(1): 34-75. DOI:10.1039/C6LC01018K
[3]	Jia Y K, Ren Y K, Liu W Y, et al. Electrocoalescence of paired droplets encapsulated in double-emulsion drops[J]. Lab on a Chip, 2016, 16(22): 4313-4318. DOI:10.1039/C6LC01052K
[4]	危卫, 张云伟, 顾兆林. 电场作用下电流变液滴的变形及力学行为[J]. 科学通报, 2013, 58(3): 197-205. (Wei Wei, Zhang Yun-wei, Gu Zhao-lin. Deformation and mechanical behavior of electrohydrodynamic droplet under external electric field[J]. Chinese Science Bulletin, 2013, 58(3): 197-205.)
[5]	陈庆国, 郑天宇, 梁雯, 等. 交流电场下液滴形变及聚结影响因素分析[J]. 电机与控制学报, 2015, 19(11): 51-58. (Chen Qing-guo, Zheng Tian-yu, Liang Wen, et al. Analysis of influencing factors for droplets deformation and coalescence in AC electric field[J]. Electric Machines and Control, 2015, 19(11): 51-58.)
[6]	Kumar S, Hesketh P J. Interpretation of AC dielectrophoretic behavior of tin oxide nanobelts using Maxwell stress tensor approach modeling[J]. Sensors and Actuators B—Chemical, 2012, 161(1): 1198-1208. DOI:10.1016/j.snb.2011.12.010
[7]	韩磊, 张艳玲, 许伟伟. 匀强电场下分散相液滴聚结行为分析[J]. 安全、健康和环境, 2018, 18(12): 62-67. (Han Lei, Zhang Yan-ling, Xu Wei-wei. Analysis of droplet coalescence behavior in dispersed phase under uniform electric field[J]. Safety Environment Health, 2018, 18(12): 62-67. DOI:10.3969/j.issn.1672-7932.2018.12.015)
[8]	Luo S R, Schiffbauer J, Luo T F. Effect of electric field non-uniformity on droplets coalescence[J]. Physical Chemistry, 2016, 18(43): 29786-29796.
[9]	Hadidi H, Kamali R, Manshadi M K D. Numerical simulation of a novel non-uniform electric field design to enhance the electrocoalescence of droplets[J]. European Journal of Mechanics B—Fluids, 2020, 80: 206-215. DOI:10.1016/j.euromechflu.2019.10.010
[10]	Zhou Y T, Li B, Zhang M Y, et al. Effect of dielectrophoresis on the coalescence of binary droplets under a non-uniform electric field[J]. Chemical Engineering Science, 2020, 224: 1-13.