任晓飞1,魏长智2,刘飞飞1,魏守水1
(1. 山东大学 控制科学与工程学院,济南 250061; 2.济南大学 信息科学与工程学院,济南 250002)
摘要:
为增进对液体微透镜变焦的动力学特征的理解,将晶格玻尔兹曼方法与电动力模型相结合,提出了一种晶格玻尔兹曼-电动力(LB-ED)方法研究介电润湿(EWOD)驱动的变焦液体微透镜原理.采用晶格玻尔兹曼方法求解Navier-Stokes方程以研究透镜的变焦过程,引入新的分布函数求解电场分布以计算驱动透镜变焦的电场力.首先数值分析了EWOD效应,并与理论分析及实验结果进行对比,验证数值方法的准确性；然后研究了电压对EWOD驱动的变焦液体透镜焦距的影响；分析了透镜变焦的动态过程；最后讨论了绝缘液体黏度对透镜响应时间及系统稳定性的影响.研究表明:不仅低电压下接触角变化与Lippmann-Young方程吻合良好,且高电压时出现接触角饱和现象,与实验结果一致,数值方法正确；根据数值计算和理论推导,建立了电压与焦距的关系；施加电压的初始时刻,电场力引起接触角突变,透镜需要延迟时间响应接触角的变化；发现绝缘液体黏度过小,系统处于振荡状态,黏度过大,系统处于过阻尼状态.合适的液体黏度可以使系统性能达到最佳.
关键词: 电流体动力学 晶格玻尔兹曼方法 介电润湿 变焦液体透镜 多组分多相流
DOI：10.11918/j.issn.0367-6234.201807068
分类号:O361.4
文献标识码:A
基金项目:国家自然科学基金(3,4);山东省自然科学基金(ZR2014EEM033)
Electrodynamic study of variable-focus liquid lens driven by electrowetting on dielectric
REN Xiaofei1,WEI Changzhi2,LIU Feifei1,WEI Shoushui1
(1.School of Control Science and Engineering, Shandong University, Jinan 250061, China; 2.School of Information Science and Engineering, University of Jinan, Jinan 250002, China)
Abstract:
With the development of micro-optical electromechanical systems, micro-lenses with adjustable focal length have become a research hotspot in the industry. In order to improve our understanding of the dynamic feature of variable-focus liquid micro-lens, a lattice Boltzmann-electrodynamic (LB-ED) method combining the lattice Boltzmann method and electrodynamic model was proposed to study the transient process of lens driven by electrowetting on dielectric (EWOD). The multi-component lattice Boltzmann method (LBM) was used to study the motion of the lens during a focusing action which is governed by Navier-Stokes equation, and then a new distribution function was introduced to calculate the electrical force. Firstly, the EWOD effect was numerically analyzed and compared with theoretical values and experimental results. Then the influence of voltage on the focal length was studied, and the dynamic process of zoom of the lens was analyzed. Finally, effects of the viscosity of the insulating liquid on the response time and stability of the lens were investigated. Results show that the change of contact angle agreed well with Loppmann-Young equation under low voltage and saturation occurred under high voltage, which indicated that the change of contact angle was consistent with the theoretical values and experimental results and verified the correctness of the numerical method. The relationship between applied voltage and focal length was established. When the voltage was switched, the contact angle changed abruptly but it cost the system some time to respond to the change in the force at the contact point. It was found that the system was in an oscillating state when the viscosity of the insulating liquid was small, and in an over-damped state when the viscosity was large. A suitable liquid viscosity could optimize the system performance.
Key words: electrohydrodynamic lattice Boltzmann method electrowetting on dielectric variable-focus liquid lens multi-component multiphase fluids
任晓飞, 魏长智, 刘飞飞, 魏守水. 介电润湿驱动的变焦液体透镜电动力学分析[J]. 哈尔滨工业大学学报, 2019, 51(10): 61-67. DOI: 10.11918/j.issn.0367-6234.201807068.
REN Xiaofei, WEI Changzhi, LIU Feifei, WEI Shoushui. Electrodynamic study of variable-focus liquid lens driven by electrowetting on dielectric[J]. Journal of Harbin Institute of Technology, 2019, 51(10): 61-67. DOI: 10.11918/j.issn.0367-6234.201807068.
基金项目 国家自然科学基金(51075243, 11672164);山东省自然科学基金(ZR2014EEM033) 作者简介 任晓飞(1989—)，男，博士研究生;
魏守水(1965—)，男，教授，博士生导师 通信作者 魏守水，sswei@sdu.edu.cn 文章历史 收稿日期: 2018-07-08
Contents -->Abstract Full text Figures/Tables PDF
介电润湿驱动的变焦液体透镜电动力学分析
任晓飞1, 魏长智2, 刘飞飞1, 魏守水1
1. 山东大学 控制科学与工程学院, 济南 250061;
2. 济南大学 信息科学与工程学院, 济南 250002
收稿日期: 2018-07-08
基金项目: 国家自然科学基金(51075243, 11672164);山东省自然科学基金(ZR2014EEM033)
作者简介: 任晓飞(1989—)，男，博士研究生;
魏守水(1965—)，男，教授，博士生导师
通信作者: 魏守水，sswei@sdu.edu.cn
摘要: 为增进对液体微透镜变焦的动力学特征的理解，将晶格玻尔兹曼方法与电动力模型相结合，提出了一种晶格玻尔兹曼-电动力(LB-ED)方法研究介电润湿(EWOD)驱动的变焦液体微透镜原理.采用晶格玻尔兹曼方法求解Navier-Stokes方程以研究透镜的变焦过程，引入新的分布函数求解电场分布以计算驱动透镜变焦的电场力.首先数值分析了EWOD效应，并与理论分析及实验结果进行对比，验证数值方法的准确性；然后研究了电压对EWOD驱动的变焦液体透镜焦距的影响；分析了透镜变焦的动态过程；最后讨论了绝缘液体黏度对透镜响应时间及系统稳定性的影响.研究表明：不仅低电压下接触角变化与Lippmann-Young方程吻合良好，且高电压时出现接触角饱和现象，与实验结果一致，数值方法正确；根据数值计算和理论推导，建立了电压与焦距的关系；施加电压的初始时刻，电场力引起接触角突变，透镜需要延迟时间响应接触角的变化；发现绝缘液体黏度过小，系统处于振荡状态，黏度过大，系统处于过阻尼状态.合适的液体黏度可以使系统性能达到最佳.
关键词: 电流体动力学 晶格玻尔兹曼方法 介电润湿 变焦液体透镜 多组分多相流
Electrodynamic study of variable-focus liquid lens driven by electrowetting on dielectric
REN Xiaofei1, WEI Changzhi2, LIU Feifei1, WEI Shoushui1
1. School of Control Science and Engineering, Shandong University, Jinan 250061, China;
2. School of Information Science and Engineering, University of Jinan, Jinan 250002, China
Abstract: With the development of micro-optical electromechanical systems, micro-lenses with adjustable focal length have become a research hotspot in the industry. In order to improve our understanding of the dynamic feature of variable-focus liquid micro-lens, a lattice Boltzmann-electrodynamic (LB-ED) method combining the lattice Boltzmann method and electrodynamic model was proposed to study the transient process of lens driven by electrowetting on dielectric (EWOD). The multi-component lattice Boltzmann method (LBM) was used to study the motion of the lens during a focusing action which is governed by Navier-Stokes equation, and then a new distribution function was introduced to calculate the electrical force. Firstly, the EWOD effect was numerically analyzed and compared with theoretical values and experimental results. Then the influence of voltage on the focal length was studied, and the dynamic process of zoom of the lens was analyzed. Finally, effects of the viscosity of the insulating liquid on the response time and stability of the lens were investigated. Results show that the change of contact angle agreed well with Loppmann-Young equation under low voltage and saturation occurred under high voltage, which indicated that the change of contact angle was consistent with the theoretical values and experimental results and verified the correctness of the numerical method. The relationship between applied voltage and focal length was established. When the voltage was switched, the contact angle changed abruptly but it cost the system some time to respond to the change in the force at the contact point. It was found that the system was in an oscillating state when the viscosity of the insulating liquid was small, and in an over-damped state when the viscosity was large. A suitable liquid viscosity could optimize the system performance.
Keywords: electrohydrodynamic lattice Boltzmann method electrowetting on dielectric variable-focus liquid lens multi-component multiphase fluids
随着微光机电系统(micro-optical electro mechanical system, MOEMS)的进展，变焦液体微透镜成为研究的热点，其改变了传统机械变焦原理，由透镜位置变为透镜曲面的变化实现变焦，克服了机械可动部件、结构复杂、易磨损、难微型化的缺点.1995年，Gorman等[1]利用介电润湿(electrowetting on dielectric, EWOD)原理，通过外加电压改变液滴表面形状；2000年，法国科学家Berge等[2]最早提出将电润湿技术应用于变焦液体透镜；Varioptic公司[3-4]先后推出了多款液体镜头聚焦照相模块；2004年，Philips公司[5]发布了实用化的变焦液体透镜，其响应时间约为几十毫秒.目前对这些微型变焦液体透镜的研究主要是采用实验手段，然而，由于受空间结构和尺度的限制，这类透镜的微制作工艺复杂、实验周期长、难度大.另一方面，EWOD的机制尚不完全清楚.对于EWOD的研究，最具代表的是电化学模型[6-7]，即微液滴与电介质层之间电荷积累产生的电容效应引起的界面能变化导致接触角变化，虽然Lippmann-Young方程预测了接触角随电压(低压)的变化情况，但通过实验，接触线上表面张力在水平方向上并不满足力学平衡[8].Peykov等[9]认为微液滴双电层中电荷密度极低，不足以驱动如此大的接触角变化，没有考虑三相边界上的电荷分布及静电力效应，造成理论与实验结果的差距较大.
晶格玻尔兹曼方法(LBM)[10]是介观尺度上多组分/多相流分析的有力工具，以其较高的计算效率和灵活的边界处理较好地解决了非连续介质界面的演变问题.多个研究多组分/多相流的模型：如颜色模型[11]、伪势模型[12]、自由能模型[13]及He-Chen-Zhang模型[14]等成功地用于微尺度流的动力学行为研究[15-20].其中，伪势模型由Shan等[12]在1993年提出，该模型的最大特点是直接对微观相互作用力进行描述，能够反映多组分/多相流体动力学的物理本质，因而得到了广泛的应用.本文结合LBM伪势模型和电动力(electrodynamic, ED)模型，提出晶格玻尔兹曼-电动力(LB-ED)方法对变焦液体微透镜进行研究.采用LBM伪势模型求解Navier-Stokes方程以研究透镜的变焦过程，引入新的分布函数求解电场分布以计算驱动透镜变焦的电场力.分析了EWOD效应，发现不仅低电压下与Lippmann-Young方程吻合良好，且高电压时与实验结果一致[21-22]，揭示了接触角饱和现象；进一步研究了电压与透镜焦距的关系、分析了透镜变焦的动态变化过程、讨论了绝缘液体黏度系数对透镜响应时间的影响.方法可以用于不同表面张力、壁面润湿性及液体黏度等复杂情况的液体透镜系统研究.
1 理论方法与数值模型 1.1 EWOD原理如图 1所示，在金属电极上涂覆电介质层以避免电解的发生.双流体界面与固体表面的接触角由接触点处的力平衡决定.当所加电压为零时，平衡接触角由Young’s方程[23]给出:
${\gamma _{12}}\cos {\theta _0} = {\gamma _{{\rm{s2}}}} - {\gamma _{{\rm{s1}}}}.$ (1)
Fig. 1
图 1 EWOD原理图 Fig. 1 Schematic diagram of EWOD 式中：γs1、γs2、γ12分别为固体表面和导电液滴、固体表面和绝缘液体及导电液滴与绝缘液体间的表面张力；θ0为导电液滴与绝缘电介质层间的初始平衡角.
当所加电压不为零时，接触角会在电场的作用下发生变化.由于液滴的体积不变，所以液滴界面曲率会在外加电压的作用下发生变化.这时Young’s方程修正如下[24]:
${\gamma _{12}}\cos \theta = {\gamma _{{\rm{s2}}}} - {\gamma _{{\rm{s1}}}} + \frac{{{\varepsilon _r}{U^2}}}{{2d}}.$ (2)
式中：θ为施加电压后的平衡接触角; U为所加电压; εr、d分别为电介质层的介电常数和厚度.联立方程(1)、(2)可得
$\cos \theta = \cos {\theta _0} + \frac{{{\varepsilon _r}{U^2}}}{{2{\gamma _{12}}d}}.$ (3)
1.2 变焦液体透镜原理图 2为变焦液体透镜的结构示意图，在透镜腔的内壁依次附着一层电极、一层电介质层以及疏水层.在镜头腔中充入两种组分的液态物质(不相溶且一种为导电体，一种绝缘体).两种液体具有不同的折射率，但具有相同的密度，以致液体界面在任何方位都能保持球面形状.通过施加电压可以有效地改变导电液体和内壁之间的接触角，从而影响透镜的形状，起到变焦的效果[2].
Fig. 2
图 2 变焦液体透镜的结构示意 Fig. 2 Schematic cross section of a variable-focus liquid lens 1.3 LB-ED数值方法 1.3.1 LBM伪势模型LBM伪势模型中，若存在S个不同组分，则伪势模型包含S个分布函数的演化方程，第k类组分的分布函数的演化方程为
$\begin{array}{*{20}{c}}{f_i^{(k)}\left( {\mathit{\boldsymbol{x}} + {\mathit{\boldsymbol{e}}_i},t + 1} \right) - f_i^{(k)}\left( {\mathit{\boldsymbol{x}},t} \right) = }\\{ - \frac{1}{{{\mathit{\boldsymbol{\tau }}^{(k)}}}}\left[ {f_i^{(k)}\left( {\mathit{\boldsymbol{x}},t} \right) - f_i^{k({\rm{eq}})}\left( {{n^{(k)}},\mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)}} \right)} \right].}\end{array}$ (4)
式中：上标k为第k类组分；fi(k)(x, t)为具有速度ei的流体粒子在t时刻，x位置处的分布函数；τ(k)为第k类组分的量纲一的松弛时间，组分k的运动黏度ν(k)与松弛时间τ(k)之间满足：ν(k)=(cs(k))2(τ(k)-0.5)Δt；fik(eq)为相应的平衡态分布函数[25]为
$f_i^{k\left( {{\rm{eq}}} \right)} = \left\{ {\begin{array}{*{20}{l}}{{\alpha _k}{n_k} - \frac{2}{3}{n_k}{{\left( {\mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)}} \right)}^2},i = 0;}\\{\frac{{\left( {1 - {\alpha _k}} \right){n^{(k)}}}}{5} + \frac{1}{3}{n^{(k)}}{\mathit{\boldsymbol{e}}_i} \cdot \mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)} + }\\{\frac{1}{2}{n_k}{{\left( {{\mathit{\boldsymbol{e}}_i} \cdot \mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)}} \right)}^2} - \frac{1}{6}{n^{(k)}}{{\left( {\mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)}} \right)}^2},i = 1,2,3,4;}\\{\frac{{\left( {1 - {\alpha _k}} \right){n^{(k)}}}}{{20}} + \frac{1}{{12}}{n^{(k)}}{\mathit{\boldsymbol{e}}_i} \cdot \mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)} + }\\{\frac{1}{8}{n_k}{{\left( {{\mathit{\boldsymbol{e}}_i} \cdot \mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)}} \right)}^2} - \frac{1}{{24}}{n^{(k)}}{{\left( {\mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)}} \right)}^2},i = 5,6,7,8,}\end{array}} \right.$ (5)
式中：αk与组分k的格子声速cs(k)有关，定义为
本文采用常用的D2Q9格子模型[27]，离散格子速度为
${\mathit{\boldsymbol{e}}_i} = \left\{ {\begin{array}{*{20}{l}}{{\bf{0}},i = 0;}\\{\left( {\cos \frac{{i - 1}}{2}{\rm{ \mathsf{ π} }},\sin \frac{{i - 1}}{2}{\rm{ \mathsf{ π} }}} \right),i = 1,2,3,4;}\\{\sqrt 2 \left( {\cos \frac{{2i - 1}}{4}{\rm{ \mathsf{ π} }},\sin \frac{{2i - 1}}{4}{\rm{ \mathsf{ π} }}} \right),i = 5,6,7,8.}\end{array}} \right.$
式中：
伪势模型通过平衡态速度ueq(k)来体现作用力的影响，即
${\rho ^{(k)}}\mathit{\boldsymbol{u}}_{{\rm{eq}}}^{(k)} = {\rho ^{(k)}}\mathit{\boldsymbol{u'}} + {\mathit{\boldsymbol{\tau}} ^{(k)}}{\mathit{\boldsymbol{F}}^{(k)}},$ (6)
其中:u′要保证在没有作用力时碰撞过程动量守恒[28]，即
$\sum\limits_k {\sum\limits_i {\frac{{{c_i}}}{{{\mathit{\boldsymbol{\tau}} ^{(k)}}}}} } \left[ {{f_{ki}} - f_{ki}^{({\rm{eq}})}\left( {{\rho _k},\mathit{\boldsymbol{u'}}} \right)} \right] = 0,$
其解为
$\mathit{\boldsymbol{u'}} = \left( {\sum\limits_k {\frac{{{\rho ^{(k)}}{\mathit{\boldsymbol{u}}^{(k)}}}}{{{\mathit{\boldsymbol{\tau}} ^{(k)}}}}} } \right)/\left( {\sum\limits_k {\frac{{{\rho ^{(k)}}}}{{{\mathit{\boldsymbol{\tau}} ^{(k)}}}}} } \right).$
F(k)是第k类流体粒子受到的总的作用力，包括：
${\mathit{\boldsymbol{F}}^{(k)}} = \mathit{\boldsymbol{F}}_{\rm{f}}^{(k)} + \mathit{\boldsymbol{F}}_{\rm{s}}^{(k)} + \mathit{\boldsymbol{F}}_{\rm{e}}^{(k)}.$
式中：Ff(k)为流体粒子间的相互作用力；Fs(k)为流体粒子与固体壁面间的相互作用力；Fe(k)为流体粒子受到的电场力.
混合流体的宏观速度定义为碰撞前后速度的平均值[29]为
$\rho \mathit{\boldsymbol{u}} = \sum\limits_k {{\rho ^{(k)}}{\mathit{\boldsymbol{u}}^{(k)}}} + \frac{1}{2}\sum\limits_k \mathit{\boldsymbol{F}} (k).$
式中：混合流体的质量密度为
压力可以通过下式求得：
$p = \frac{1}{3}\sum\limits_k {{n^{(k)}}} + \frac{3}{2}\sum\limits_{k\bar k} {{G_{k\bar k}}{n^{(k)}}{n^{(\bar k)}}} .$
在伪势模型中，流体粒子间的相互作用力为
$\mathit{\boldsymbol{F}}_f^{(k)} = - {\psi ^{(k)}}\left( \mathit{\boldsymbol{x}} \right)\sum\limits_i {\sum\limits_k {{G_{k\bar k}}\left( {\mathit{\boldsymbol{x}} + {\mathit{\boldsymbol{e}}_i}} \right){\psi ^{(\bar k)}}\left( {\mathit{\boldsymbol{x}} + {\mathit{\boldsymbol{e}}_i}} \right){e_i}} } ,$ (7)
式中:ψ(k)为组分k的有效密度，为了计算方便，本文中取
${\psi ^{(k)}} = {n^{(k)}},{G_{k\bar k}}\left( {\left| {\mathit{\boldsymbol{x}} - \mathit{\boldsymbol{x'}}} \right|} \right) = \left\{ \begin{array}{l}{g_{k\bar k}},\left| {\mathit{\boldsymbol{x}} - \mathit{\boldsymbol{x'}}} \right| = 1;\\{g_{k\bar k}}/4,\left| {\mathit{\boldsymbol{x}} - \mathit{\boldsymbol{x'}}} \right| = \sqrt 2 ;\\0,{\rm{otherwise}},\end{array} \right.$
gkk决定了组分k和k粒子之间相互作用的强度.
类似地，将固态壁面看作一个密度为常数的相，从而引入流体粒子与固体壁面间的相互作用力为
$\mathit{\boldsymbol{F}}_{\rm{s}}^{(k)} = - {\psi ^{(k)}}\left( \mathit{\boldsymbol{x}} \right)\sum\limits_i {{G_{ks}}s\left( {\mathit{\boldsymbol{x}} + {\mathit{\boldsymbol{e}}_i}} \right){\mathit{\boldsymbol{e}}_i}} .$ (8)
式中：s(x)为指示函数，当x是固体点时，s(x)=1，当x是液体点时，s(x)=0；Gks为组分k与固态壁面之间的相互作用，Gks>0表示壁面为非润湿性，Gks < 0表示壁面为润湿性.通过改变Gks，可以获得不同的壁面润湿性.
1.3.2 电动力模型流体粒子所受的电场力可以表示为Maxwell应力张量或体积力[30-31]为
${\mathit{\boldsymbol{F}}_{\rm{e}}} = - \frac{1}{2}\mathit{\boldsymbol{E}} \cdot \mathit{\boldsymbol{E}}\nabla \varepsilon + q\mathit{\boldsymbol{E}} + \frac{1}{2}\nabla \left[ {\mathit{\boldsymbol{E}} \cdot \mathit{\boldsymbol{E}}\frac{{\partial \varepsilon }}{{\partial \rho }}\rho } \right].$ (9)
式中：E为外加电场强度; ε为流体介电常数; q为自由电荷密度.式(9)右侧第1项源于极化应力，并沿着界面法线方向作用；第2项为沿着电场方向作用的静电力；第3项是由液体密度的改变引起的，本文中，由于液体的不可压缩性，该项可以忽略.
因此电场力可以简化为
${\mathit{\boldsymbol{F}}_{\rm{e}}} = - \frac{1}{2}\mathit{\boldsymbol{E}} \cdot \mathit{\boldsymbol{E}}\nabla \varepsilon + q\mathit{\boldsymbol{E}},$ (10)
电场强度为
$\mathit{\boldsymbol{E}} = - \nabla U,$ (11)
自由电荷密度为
$q = \nabla \cdot (\varepsilon \mathit{\boldsymbol{E}}).$ (12)
为了计算液滴所受电场力，首先通过以下控制方程计算EWOD系统中电势分布：
$\nabla \cdot \left[ {\sigma \nabla U} \right] = 0,\nabla \cdot \left( {\varepsilon \nabla U} \right) = 0,$
其中σ为流体电导率.
在不同组分界面处，Laplace方程满足连续性边界条件：
$\vec n \cdot \left( {{\mathit{\boldsymbol{D}}_1} - {\mathit{\boldsymbol{D}}_2}} \right) = {q_{\rm{s}}},\vec n \times \left( {{\mathit{\boldsymbol{E}}_1} - {\mathit{\boldsymbol{E}}_2}} \right) = 0.$
式中：qs为表面电荷密度; D1、D2分别为界面两侧的电位移；E1、E2分别为界面两侧的场强.
微液滴和电介质层表面可认为没有自由电荷，因此两种介质界面处的边界条件为：
${\varepsilon _1}\frac{{\partial {U_1}}}{{\partial \mathit{\boldsymbol{n}}}} = {\varepsilon _2}\frac{{\partial {U_2}}}{{\partial \mathit{\boldsymbol{n}}}},{U_1} = {U_2}.$
在此采用LBM方法求解上述控制方程.为此，引入一个新的分布函数hi，其演化方程为[32]
$\begin{array}{*{20}{l}}{{h_i}\left( {\mathit{\boldsymbol{x}} + {\mathit{\boldsymbol{e}}_i},t + 1} \right) - {h_i}\left( {\mathit{\boldsymbol{x}},t} \right) = }\\{ - \frac{1}{{{\mathit{\boldsymbol{\tau}} ^h}}}\left[ {{h_i}\left( {\mathit{\boldsymbol{x}},t} \right) - h_i^{{\rm{eq}}}\left( {\mathit{\boldsymbol{x}},t} \right)} \right],}\end{array}$ (13)
其中，平衡态分布函数hieq为
$h_i^{{\rm{eq}}} = {\omega _i}U.$ (14)
权重系数ω0=4/9，ω1-4=1/9，ω5-8=1/36.松弛时间τh=3σ+0.5，其中
电势可由下式获得：
$U = \sum\limits_i {{h_i}} .$ (15)
每种组分的流体粒子所受的电场力为
$\mathit{\boldsymbol{F}}_{\rm{e}}^{(k)} = \frac{{{\mathit{\boldsymbol{\rho }}^{(k)}}}}{\rho }{\mathit{\boldsymbol{F}}_{\rm{e}}}.$ (16)
1.4 LB-ED算法流程LB-ED算法流程如图 3所示，其中：模型初始化为计算区域网格划分、流场初始化、电场初始化、设置边界条件；输入为电极电压U、最大迭代次数Tmax、时间步t=0.
Fig. 3
图 3 LB-ED算法流程 Fig. 3 Flow chart of LB-ED algorithm 2 结果与讨论 2.1 算法验证为了EWOD系统有效工作，约束条件为Bo=ρgL2/γ < 1[33]，其中:ρ为液体密度；g为重力加速度；γ为表面张力；L为特征尺寸.下文各参数单位均为格子单位(lattice unit, lu).
根据EWOD原理，在壁面润湿性一定时，通过改变所加电压，可以改变液滴与壁面的接触角.为了验证LB-ED方法的正确性，首先探究电压对接触角的影响.计算区域网格设定为201×101，四周边界为固态壁面，采用无滑移反弹边界条件；形成液滴的流体记为组分1，液滴外的流体记为组分2.液滴介电常数取ε1=5，电导率σ1=0.1，电介质层与绝缘液体取相同的介电常数ε2=εr=2，电介质层厚度取d=1，其他参数取n(1)=n(2)=1，m(1)=m(2)=1，τ(1)=τ(2)=1，G1s=0，G2s=0，gkk=0.2.初始时刻，半圆形液滴置于底部壁面上，初始平衡态接触角为90°.
对于施加电压后达到稳定状态的液滴，其半径r可以通过下式求得[17]:
$r = \frac{{{a_0}}}{2} + \frac{{b_0^2}}{{8{a_0}}}.$
液滴接触角满足以下关系:
$\tan (\theta ) = \frac{{{b_0}}}{{2\left( {r - {a_0}} \right)}}.$
式中:a0、b0分别为液滴顶部距离壁面的高度和液滴与壁面接触线的长度(如图 1所示).
图 4给出了电压分别为0、0.2、0.4、0.5、0.6、0.7时液滴达到稳态时的形状.随着电压的增大，液滴与壁面的接触角减小，当电压大于0.6 V时，接触角不再随着电压的增大而减小.
Fig. 4
图 4 电压U=0、0.2、0.4、0.5、0.6、0.7时，液滴的稳态形状 Fig. 4 Steady-state shape of droplets under the conditions of U=0, 0.2, 0.4, 0.5, 0.6, 0.7 图 5、6分别给出了本文数值方法得到的不同电介质层厚度d及不同介电常数εr条件下，电压U与Δcos θ=cos θ-cos θ0的关系，并与式(3)得到的理论值进行比较.离散点为不同条件下的数值计算结果，实线为相应条件下的理论值.由图 5中d=1、d=2及d=3条件下的离散点可以看出，当电压分别达到0.55、0.80、0.95时，Δcos θ值趋于稳定.由图 6发现，εr=4、εr=3和εr=2条件下，当电压分别达到0.40、0.45、0.55时，Δcos θ值趋于稳定.接触角随着电压的增大而减小，并且，电介质层越厚，要达到同一接触角所需施加的电压越大；电介质层的介电常数越大，获得相同接触角所需的电压越小，且接触角的变化也越快.结果符合Lippmann-Young方程在低压时的表现；同时超过一定电压时，接触角不再随电压的增大而改变的饱和现象，解决Lippmann-Young方程无法预测高电压的情况.Vallet等[21-22]对介电润湿的实验中也发现了接触角饱和现象，验证了本文方法的合理性.
Fig. 5
图 5 不同电介质层厚度条件下，电压对接触角的影响 Fig. 5 Influence of voltage on the contact angle for different thicknesses of the dielectric layer Fig. 6
图 6 不同介电常数条件下，电压对接触角的影响 Fig. 6 Influence of voltage on the contact angle for different permittivity of dielectric layer 2.2 变焦液体透镜的分析将计算区域网格设定为Lx×Ly=201×101，四周边界为固态壁面，采用无滑移反弹边界条件；导电液体记为组分1，绝缘液体记为组分2.粒子数密度取值为n(1)=n(2)=1，分子质量为m(1)=m(2)=1，松弛时间为τ(1)=τ(2)=1，对应的运动黏度系数为ν(1)=ν(2)=0.167，G1s=0，G2s=-0.15，液-液相互作用强度gkk=0.2，由Laplace公式可以得到界面张力γ12=0.409，导电液体介电常数取ε1=5，电导率σ1=0.1，电介质层与绝缘液体取相同的介电常数ε2=εr=2，电介质层厚度取d=1，透镜腔半径R=(Lx-2d)/2，导电液体和绝缘液体的折射率分别取值为n1=1.394 0，n2=1.560 2.
2.2.1 电压对透镜焦距的影响透镜焦距与曲率半径的关系为
$\frac{1}{f} = \left( {{n_2} - {n_1}} \right)\left[ {\frac{1}{{{r_1}}} - \frac{1}{{{r_2}}} + \frac{{\left( {{n_2} - {n_1}} \right){d_l}}}{{{n_1}{n_2}{r_1}{r_2}}}} \right].$ (17)
式中：r1、r2分别为透镜两侧的曲率半径；n1、n2分别为导电液体和绝缘液体的折射率；dl为透镜厚度.
在此，透镜一侧为平面，所以r1为无穷大，又由几何关系知道，r2=R/sin(
$f = \frac{R}{{\left( {{n_2} - {n_1}} \right)\cos \theta }}.$
图 7为初始角为135°时，电压与焦距的关系.从图 7可以发现，当电压较小时，透镜焦距变化缓慢，随着电压的增大，透镜的焦距迅速减小，该过程对应的透镜为凹透镜.曲线开始变化缓慢的原因是，由于疏水层对导电液体的约束，低电压时，电场力对接触角的影响较小. Berge等[2]于2000年利用介电润湿原理设计的液体透镜，将Na2SO4溶液和α-氯萘油滴封装成3 mm的变焦透镜.通过实验测量发现，当电压低于一定值时，透镜焦距变化很小，在电压达到一定值(原文80 V)时出现拐点，透镜焦距迅速变化，本文的数值结果与其实验结果吻合良好.工作电压约为0.53时，液体界面形状为平面.透镜焦距从负无穷跃变为正无穷，透镜从凹透镜转变为凸透镜.随着电压的继续增大，焦距从无穷大迅速减小，直至达到一个稳定值.从图 7可以发现，当电压达到0.80时，焦距几乎不再随电压的增大而变化，原因是在高电压下，导电液体与透镜腔内壁的接触角达到饱和.
Fig. 7
图 7 初始角为135°时，电压与焦距的关系 Fig. 7 Relationship between the voltage U and the focal length f under the condition that the initial angle is 135° 因此，根据数值计算结果，结合电压对接触角的影响，可以得到液体变焦透镜中电压与焦距的关系：
$f = \frac{R}{{\left( {{n_2} - {n_1}} \right)\cos {\theta _0} + \frac{{\left( {{n_2} - {n_1}} \right)}}{{{\gamma _{12}}}} \cdot \frac{{{\varepsilon _{\rm{r}}}}}{{2d}}{U^2}}}.$
2.2.2 透镜性能分析本文研究了透镜在变焦过程中的运动情况，图 8给出了透镜切换过程中不同时步(time step, ts)的速度矢量图.图中箭头表示速度方向，颜色表示速度大小.图 8(a)~8(d)为透镜从凸面变为平面的过程，从图 8(a)可以看出，施加电压的初始时刻，由于电场力的作用，透镜与壁面的接触点获得向上运动的速度，接触角突然改变，此时透镜其他部位速度为零.为了响应接触角的变化，靠近接触点的部位开始运动，随后引起透镜中心部位向下运动.透镜中心速度首先迅速增大，当t=2 400 ts时，速度开始减小.当t=3 600 ts时，透镜变为平面. 图 8(d)~8(f)为透镜从平面变为凹面的过程.该过程中，透镜形状缓慢变化直至达到稳定状态.由图 8可以看出，在外加电压的作用下，由于液-液界面的易变性，两种液体的界面可以从凸面变成凹面，透镜不仅可以改变焦距，而且也能改变透镜的正负类型，焦距可调范围广.
Fig. 8
图 8 t=10、1 200、2 400、3 600、4 800、6 000 s时刻的透镜形状及速度分布 Fig. 8 Shapes of the lens and the velocity distribution at t=10, 1 200, 2 400, 3 600, 4 800, 6 000 s 图 9给出了不同的绝缘液体黏度条件下，透镜中心处的位移随时间的变化关系.从图 9可以看出，曲线开始变化缓慢，随后迅速减小，最终经过振荡或平滑达到稳定状态.原因是施加电压的初始时刻，电场力与张力导致接触角突变，但是透镜中心需要延迟时间来响应接触角的变化.分析发现，当ν(2)>4ν(1)时，系统处于过阻尼状态，透镜进入稳定状态需要较长时间；当ν(2)≤3ν(1)时，系统处于欠阻尼状态，透镜要经过多次振荡才能达到稳定状态；当ν(2)=4ν(1)时，系统接近临界阻尼状态.从图 9可以看出，当系统处于临近阻尼状态时，透镜达到稳态所需时间最短，系统性能最佳.可以发现，黏度系数对界面位移的改变没有影响，但对系统的响应时间及稳定性能有较大影响.因此合理选择液体黏度，可以使系统即具有令人满意的响应快速性，又具有比较好的响应平稳性.
Fig. 9
图 9 不同的绝缘液体黏度条件下，透镜中心处位移随时间的变化 Fig. 9 Location of the center of the lens as a function of time for different values of the viscosity of the insulating fluid 3 结论1) 数值分析了介电润湿效应，发现不仅低电压下与Lippmann-Young方程吻合，且在高电压时与实验结果一致，揭示了接触角饱和现象，验证了数值方法的正确性.
2) 研究了电压对透镜焦距的影响，建立了电压与焦距的关系.同时研究了透镜从凹透镜到凸透镜的转变，分析了透镜变焦的动态演变过程，发现施加电压的初始时刻，电场力与张力导致接触角突变，透镜需要延迟时间来响应接触角的变化；
3) 讨论了绝缘液体黏度对透镜响应时间的影响，当绝缘液体黏度系数为导电液体黏度系数4倍时，系统处于临界状态，透镜响应时间最短且平稳性最佳；当绝缘液体黏度系数大于导电液体黏度系数4倍时，系统处于过阻尼状态，透镜达到稳定状态需要较长时间，无振荡；当绝缘液体黏度系数小于导电液体黏度系数4倍时，系统处于欠阻尼状态，透镜在振荡中达到稳态.
4) 该方法可以用于不同表面张力、壁面润湿性及液体黏度等复杂情况的液体透镜系统研究，增进了对透镜变焦的动力学特征的理解.
参考文献
[1] GORMAN C B, BIEBUYCK H A, WHITESIDES G M. Control of the shape of liquid lenses on a modified gold surface using an applied electrical potential across a self-assembled monolayer[J]. Langmuir, 1995, 11(6): 2242. DOI:10.1021/la00006a063
[2] BERGE B, PESEUX J. Variable focal lens controlled by an external voltage: An application of electrowetting[J]. The European Physical Journal E, 2000, 3(2): 159. DOI:10.1007/s101890070029
[3] GABAY C, BERGE B, DOVILLAIRE G, et al. Dynamic study of a Varioptic variable focal lens[C]//Proceedings of the SPIE 4767, Current Developments in Lens Design and Optical Engineering Ⅲ. Seattle, WA: SPIE Press, 2002: 159. DOI: 10.1117/12.468224
[4] BERGE B. Liquid lens technology: Principle of electrowetting based lenses and applications to imaging[C]//Proceedings of the 18th IEEE International Conference on Micro Electro Mechanical Systems. Miami Beach, FL: IEEE Press, 2005: 227. DOI: 10.1109/MEMSYS.2005.1453908
[5] KUIPER S, HENDRIKS B H W. Variable-focus liquid lens for miniature cameras[J]. Applied Physics Letters, 2004, 85(7): 1128. DOI:10.1063/1.1779954
[6] KUMARI N, BAHADUR V, GARIMELLA S V. Electrical actuation of dielectric droplets[J]. Journal of Micromechanics and Microengineering, 2008, 18(8): 085018. DOI:10.1088/0960-1317/18/8/085018
[7] WALKER S W, SHAPIRO B. Modeling the fluid dynamics of electrowetting on dielectric (EWOD)[J]. Journal of Microelectromechanical Systems, 2006, 15(4): 986. DOI:10.1109/JMEMS.2006.878876
[8] DIGILOV R. Charge-induced modification of contact angle: The secondary electrocapillary effect[J]. Langmuir, 2000, 16(16): 6719. DOI:10.1021/la991308a
[9] PEYKOV V, QUINN A, RALSTON J. Electrowetting:A model for contact-angle saturation[J]. Colloid and Polymer Science, 2000, 278(8): 789. DOI:10.1007/s003960000333
[10] CHEN Shiyi, DOOLEN G D. Lattice Boltzmann method for fluid flows[J]. Annual Review of Fluid Mechanics, 1998, 30(1): 329. DOI:10.1146/annurev.fluid.30.1.329
[11] ROTHMAN D H, KELLER J M. Immiscible cellular-automaton fluids[J]. Journal of Statistical Physics, 1988, 52(3/4): 1119. DOI:10.1007/BF01019743
[12] SHAN Xiaowen, CHEN Hudong. Lattice Boltzmann model for simulating flows with multiple phases and components[J]. Physical Review E, 1993, 47(3): 1815. DOI:10.1103/PhysRevE.47.1815
[13] SWIFT M R, OSBORN W R, YEOMANS J M. Lattice Boltzmann simulation of nonideal fluids[J]. Physical Review Letters, 1995, 75(5): 830. DOI:10.1103/PhysRevLett.75.830
[14] HE Xiaoyi, CHEN Shiyi, ZHANG Raoyang. A lattice Boltzmann scheme for incompressible multiphase flow and its application in simulation of Rayleigh-Taylor instability[J]. Journal of Computational Physics, 1999, 152(2): 642. DOI:10.1006/jcph.1999.6257
[15] FAROKHIRAD S, MORRIS J F, LEE T. Coalescence-induced jumping of droplet: Inertia and viscosity effects[J]. Physics of Fluids, 2015, 27(10): 102102. DOI:10.1063/1.4932085
[16] ASHOKE RAMAN K, JAIMAN R K, LEE T S, et al. A numerical study on electrowetting-induced jumping and transport of droplet[J]. International Journal of Heat and Mass Transfer, 2016, 99: 805. DOI:10.1016/j.ijheatmasstransfer.2016.04.038
[17] XU Zhiyuan, LIU Haihu, VALOCCHI A J. Lattice Boltzmann simulation of immiscible two-phase flow with capillary valve effect in porous media[J]. Water Resources Research, 2017, 53(5): 3770. DOI:10.1002/2017WR020373
[18] DONG B, YAN Y Y, LI W Z. LBM simulation of viscous fingering phenomenon in immiscible displacement of two fluids in porous media[J]. Transport in Porous Media, 2011, 88(2): 293. DOI:10.1007/s11242-011-9740-y
[19] LIU Haihu, JU Yaping, WANG Ningning, et al. Lattice Boltzmann modeling of contact angle and its hysteresis in two-phase flow with large viscosity difference[J]. Physical Review E, 2015, 92(3): 033306. DOI:10.1103/PhysRevE.92.033306
[20] 周超英, 葛家, ISLAMShams Ul. 平板间方柱绕流的格子Boltzmann方法模拟[J]. 哈尔滨工业大学学报, 2010, 42(11): 1753.
ZHOU Chaoying, GE Jia, ISLAM Shams Ul. Numerical simulation of confined flow around a square cylinder in a channel with lattice Boltzmann method[J]. Journal of Harbin Institute of Technology, 2010, 42(11): 1753. DOI:10.11918/j.issn.0367-6234.2010.11.016
[21] CHO S K, MOON H, KIM C J. Creating, transporting, cutting, and merging liquid droplets by electrowetting-based actuation for digital microfluidic circuits[J]. Journal of Microelectromechanical Systems, 2003, 12(1): 70. DOI:10.1109/JMEMS.2002.807467
[22] VALLET M, VALLADE M, BERGE B. Limiting phenomena for the spreading of water on polymer films by electrowetting[J]. The European Physical Journal B-Condensed Matter and Complex Systems, 1999, 11(4): 583. DOI:10.1007/s100510051186
[23] YOUNG T. Ⅲ. An essay on the cohesion of fluids[J]. Philosophical Transactions of the Royal Society of London, 1805, 95: 65. DOI:10.1098/rstl.1805.0005
[24] LIPPMANN M G. Relation entre les phénomènes électriques et capillaries[J]. Annales de Chimie et de Physique, 1875, 5: 494.
[25] CHEN Hudong, CHEN Shiyi, MATTHAEUS W H. Recovery of the Navier-Stokes equations using a lattice Boltzmann method[J]. Physical Review A, 1992, 45(8): R5339. DOI:10.1103/PhysRevA.45.R5339
[26] KANG Qinjun, ZHANG Dongxiao, CHEN Shiyi. Displacement of a two-dimensional immiscible droplet in a channel[J]. Physics of Fluids, 2002, 14(9): 3203. DOI:10.1063/1.1499125
[27] QIAN Y H, D'HUMIèRES D, LALLEMAND P. Lattice BGK models for Navier-Stokes equation[J]. Europhysics Letters, 1992, 17(6): 479. DOI:10.1209/0295-5075/17/6/001
[28] SHAN Xiaowen, CHEN Hudong. Simulation of nonideal gases and liquid-gas phase transitions by the lattice Boltzmann equation[J]. Physical Review E, 1994, 49(4): 2941. DOI:10.1103/PhysRevE.49.2941
[29] SHAN Xiaowen, DOOLEN G. Multicomponent lattice-Boltzmann model with interparticle interaction[J]. Journal of Statistical Physics, 1995, 81(1/2): 379. DOI:10.1007/BF02179985
[30] MELCHER J R, TAYLOR G I. Electrohydrodynamics: A review of the role of interfacial shear stresses[J]. Annual Review of Fluid Mechanics, 1969, 1(1): 111. DOI:10.1146/annurev.fl.01.010169.000551
[31] SAVILLE D A. Electrohydrodynamics: The Taylor-Melcher leaky dielectric model[J]. Annual Review of Fluid Mechanics, 1997, 29(1): 27. DOI:10.1146/annurev.fluid.29.1.27
[32] ZHANG Junfeng, KWOK D Y. A 2D lattice Boltzmann study on electrohydrodynamic drop deformation with the leaky dielectric theory[J]. Journal of Computational Physics, 2005, 206(1): 150. DOI:10.1016/j.jcp.2004.11.032
[33] MUGELE F, BARET J C. Topical review: Electrowetting from basics to applications[J]. Journal of Physics Condensed Matter, 2005, 17(28): R705. DOI:10.1088/0953-8984/17/28/R01
