4.4.6 建模设计案例分析
例4.6 管道流动的三维模拟。
微流控芯片中，一般情况下都会涉及通道流动操作，因此通道流(也称泊肃叶流)是建立微流控芯片模型的基础。这里以建立一个通道流模型为例。
1.通道流的概念与基本原理
所谓通道流，是指无限长直圆管中的不可压缩黏性层流，此流动可由Navier-Stokes方程组描述。19世纪G.H.L.哈根和J.L.M.泊肃叶从实验角度归纳出圆管内的流动规律，后来证实与理论计算符合[5]，即圆管截面上的速度分布为绕中心线的旋转抛物面，相应表达式为
式中，ΔP —— 在长度为L的圆管通道上的压力差(或压力损失)；
D —— 圆管通道的直径；
μ—— 流体的动力黏度；
r —— 距离圆管通道中轴线的距离。
流量Q和压力差ΔP的公式如下：
式中，g —— 重力加速度。
由式(4.116)可以看出，在其他条件不变的情况下，压力差与流量线性相关。
2.建立通道流模型
拟建立通道流模型的主要信息：圆管的直径为100 μm，长度为600 μm，环境温度为293 K，流体为水。采用压力驱动，压力差为10 Pa。
1)进入开发模式
按照如下导航流程选择应用模式：【新建】→【模型向导】→【三维】→【流体流动】→【单相流】→【层流】。双击【层流】，即可将层流模式添加进“添加的物理场接口”。此时，可以在右侧窗口中看到“检查物理场接口”下显示的三个速度分量“u，v，w”和一个压力分量“p”，这是本案例中要求解的变量，也是Navier-Stokes方程组的自变量。由于是不可压缩流动，流体的密度是常量，用稳态模型描述流动即可。在【研究】中，双击【稳态】即可进入开发模式界面。
2)具体建模流程
将工程文件名设为“Flow”，默认的COMSOL Multiphysics工程文件扩展名为.mph。因此可从＜模型开发器＞中看到“Flow.mph(root)”的字样。接下来，按照前面介绍的建模流程进行操作。
第1步：设置参数与变量
设置环境温度。操作步骤：＜模型开发器＞→【Flow.mph(root)】→【全局定义】→【参数】，在＜参数设置区＞设置【名称】为“T”，设置【表达式】为“293”，设置【描述】为“温度”。T作为温度通用参数，在模型中可以自动识别。
第2步：构建几何模型
操作步骤：＜模型开发器＞→【Flow.mph(root)】，右击【几何】，在弹出的快捷菜单中选择【圆柱体】，进入＜参数设置区＞。将【长度单位】选择“μm”。在模型树中，单击【圆柱体1】(如果不喜欢这个命名，可以右击【圆柱体1】后重命名，如命名为“管道”)。此时，在＜参数设置区＞中可见圆柱体的尺寸、位置、旋转角度等空间位置定义。几个关键的参数设置如下：
【大小和形状】半径：50；高度：600。
【位置】x：0；y：0；z：0。位置坐标是圆柱底面圆心的位置。
【轴】轴指圆柱的轴心线，默认为“z轴”。
【工作平面】选择“xy平面”为工作平面。
设置完毕后，单击＜参数设置区＞顶部的【构建所有对象】，即可在＜模型可视化区＞得到如图4.17所示的图形。
图4.17 三维圆形管道几何模型(书后附彩插)
图4.18 弹出的快捷菜单
第3步：设置物理场属性
设置流动物质属性。＜模型开发器＞→【Flow.mph(root)】→【组件1】→【材料】，右击【材料】，弹出的快捷菜单如图4.18所示。其中，通过【从库中添加材料】选项可以直接导入需要添加材料的主要属性。对于流体来说，密度和黏度是最关键的两个参数。本例流动介质选择“水”。操作步骤：单击【从库中添加材料】，转至＜自定义菜单区＞，选择【液体气体】→【Liquids】→【Water】，双击【Water】即可将其添加至“模型树”中的【材料】项。此时单击【材料】→【Water】，即可在＜参数设置区＞看到与水相关的属性。可以看到，激活了【动力黏度】和【密度】两项，且这两项均是环境温度“T”的函数。当环境温度“T”在计算过程中发生变化时，水的动力黏度和密度也会发生变化。
如果不使用材料库中提供的材料属性，也可以选择图4.18所示菜单中的【空材料】，手动定义动力黏度和密度的数值。例如，定义动力黏度为0.001 Pa·s，密度为1 000 kg/m3，这样设置的参数就与环境温度不相关。
本例中，选择材料库中的【Water】。
第4步：设置边界条件
对于本例，管道的边界条件分为三种——入口、出口和壁面。设置方法如下：＜模型开发器＞→【Flow.mph(root)】→【组件1】→【层流】，右击【层流】，在弹出的快捷菜单中分别选择【入口】、【出口】，以添加边界条件种类。在流动模拟中，壁面边界是默认边界，在设置其他边界类型时，相应壁面边界会被新设置的边界类型覆盖。
单击【入口1】，在＜参数设置区＞可对入口边界进行设置。在＜模型可视化区＞单击z=0处的圆柱端面，即可在【边界选择】窗口中增加一个数字“3”，这里的数值是选择边界的序号。在【边界条件】的下拉菜单中选择【压力】，并在下面的【压力条件】中设置数值为“10”，表示入口压力为10 Pa。类似地，单击【出口1】，然后在＜模型可视化区＞单击z=600处的圆柱端面，该端面的序号为4。出口边界默认为压力边界条件，压力默认值为0，因此不用进行人为设置。圆柱的其他边界面已经默认为壁面，所以也不用再进行设置。
第5步：网格划分
单击“模型树”中的【网格】，在＜参数设置区＞进行单元格大小设置，其选项包括“极细化”“超细化”“较细化”“常规”“较粗化”“超粗化”“极粗化”。此处选择【常规】，单击＜参数设置区＞顶部的【全部构建】，即可在＜模型可视化区＞查看构建结果。需要说明的是，在圆柱通道内部的网格的最大网格尺寸不大于管道直径的15%为宜，若过于粗化则解的精度不够，若过于细化则计算量会非常大。理论上，不论网格如何划分，解的结果都应该是相同的，此为数值解的“网格无关性”。得到的边界网格如图4.19所示。
图4.19 圆柱通道端面的有限元网格
可见，最大网格的边长为管道直径的10%左右，且壁面边界处的网格自动加密，这样的网格划分结果是可以接受的。
第6步：配置求解器
配置方法：在＜模型开发器＞→【Flow.mph(root)】→【研究】中进行配置。这里的【研究】和【全局】、【组件1】的层次是并列的，属于“模型树”的主干。本例不涉及多场耦合，是最简单也是最基础的模型，求解器的选择及相关参数采用默认值即可。
第7步：求解
求解方法：＜模型开发器＞→【Flow.mph(root)】→【研究】，右击其下选项【步骤1：稳态】，在弹出的快捷菜单中选择【计算选定步骤】，即可开始计算。
第8步：后处理
一旦求解完成，在＜模型开发器＞→【Flow.mph(root)】→【结果】中，可以进行多种类型的后处理操作。在这一步骤，可以观察求解结果的分布情况，如可以查看在任意点、线或截面的解的分布。
1)查看某截面的速度分布
例如，要查看在z轴向长度为200和400处的速度大小分布，可采用如下操作：【结果】→【速度】→【切面1】，进入＜参数设置区＞，在【平面数据】中设置【平面类型】为【常规】，在【平面定义方法】中选择【三点】(即三点 确 定 一个平 面)，设 置 下 面 的“x，y，z”坐 标。首 先，设 置z=200 μm的切面。在设置三个点的位置数据时，确保z坐标均设置为200，x、y坐标可任意设置，但三个点不能共线。设置完毕后，单击＜参数设置区＞上端的【绘制】，即可得到z=200 μm的速度切面。然后，建立另一个速度切面，找到【结果】→【速度】，右击【速度】，在弹出的快捷菜单中选择【切面】，即可再建立一个切面图。与切面1比较，需要将z方向的坐标均设为400。最终结果如图4.20所示。
2)查看z=200 μm处直径上的速度分布
本操作需要两个步骤：
第1步，建立三维截线。找到【结果】→【数据集】，右击【数据集】，在弹出的快捷菜单中选择【三维截线】，三维截线的定义方法采用【两点】。选择在“yz平面”的直径为截线，则定义截线的两点的(x，y，z)坐标分别为(0，-50，200)和(0，50，200)。在＜参数设置区＞顶部单击【绘制】，即可完成设置，获得的对象是【三维截线1】。
图4.20 圆管通道速度场截面图(书后附彩插)
第2步，显示沿着截线的速度分布。右击【结果】，在弹出的快捷菜单中选择【一维绘图组】，即可发现在【结果】下面增加了【一维绘图组】选项，右键单击，在弹出的快捷菜单中选择【线图】，即可建立一个对象【线图1】，进入【线图1】设置，将【数据集】选择为【三维截线1】，其他默认设置，单击【绘制】，即可得到【三维截线1】上的速度分布曲线，如图4.21(a)所示。
采用式(4.115)获得的理论值曲线如图4.21(b)所示，其中实线部分为理论值，圆圈点为从图4.21(a)中提取的计算值，可见两者十分接近，说明理论计算结果与数值模拟结果吻合。提取数据的方法：【结果】→【导出】，右键单击，在弹出的快捷菜单中选择【数据】，即可生成对象【数据1】；在＜参数设置区＞，将【数据集】选择为【三维截线1】。在【输出】项设置数据文件的导出地址，然后在＜参数设置区＞顶端单击【导出】，即可完成数据提取。
图4.21 速度剖面曲线(书后附彩插)
(a)速度剖面曲线；(b)速度理论值与计算值比较
理论值计算的MATLAB脚本如下：
例4.7 复杂管道流动的三维模拟。
本任务是建立一个具有复杂形状的三维流道模型，通道的形状和尺寸如图4.22所示，本例中的环境温度为293 K，流体为水。流动采用出入口压力差驱动，相应压力差为10 Pa。
图4.22 通道模型
其几何模型的建立思路：首先构建通道截面的二维模型，然后通过“拉伸”操作，形成一个三维模型。其中，通道截面的二维模型可以采用多个矩形进行布尔运算得到。本模型依然需要采用流体流动模块，采用导航流程【新建】→【模型向导】→【三维】→【流体流动】→【单相流】→【层流】，选用稳态求解器，单击【完成】即可进入模型开发界面，开始建模。保存此工程文件为ComplexFlow.mph。
具体建模流程如下。
第1步：设置参数与变量
所需的参数可以直接在“模型树”的【全局定义】→【参数】中进行定义，如表4.1所示。
表4.1 计算参数表
其中，“名称”“表达式”和“描述”需要手动输入。在“表达式”一栏，当输入参数值时，后面用中括号[]来标识物理单位，物理单位之间通过数学运算符连接。“值”根据“表达式”的输入自动显示。关于参数设置，还可以借助参数设置栏下面的功能图标，采用导入、导出的方式进行批量操作。
第2步：构建几何模型
采用7个矩形进行布尔运算，获得通道模型的纵截面，矩形的坐标位置及编号如图4.23所示。各矩形的形状和位置参数如表4.2所示。
图4.23 几何模型的二维平面规划
表4.2 几何模型中的矩形设置参数 μm
续表
由于首先要建立一个可供拉伸的二维模型，所以要先建立一个工作平面。建立流程：＜模型开发器＞→【ComplexFlow.mph(root)】→【组件1】→【几何1】，在＜参数设置区＞选择【长度单位】为“μm”。右击【几何1】，在弹出的快捷菜单中选择【工作平面】，即可建立一个工作平面对象【工作平面1】。在＜参数设置区＞选择【平面】为【xy平面】，【z坐标】默认为“0 μm”。单击【工作平面1】下的【平面几何】，＜模型可视化区＞进入几何设计模式。
建立7个矩形：右击【平面几何】，在弹出的快捷菜单中选择【矩形】，这样操作7次；然后按照表4.3中的数据，在＜参数设置区＞设置矩形的尺寸和位置，其中“xw”表示矩形左下角x方向坐标，“yw”表示矩形左下角y方向坐标。设置完成后，单击＜参数设置区＞顶端的【全部构建】即可。
为了建立一个整体的二维通道模型，还需要对这7个矩形进行布尔运算，分为两部分进行：建立差集、建立并集。
1)建立差集
单击“模型树”中【几何1】→【工作平面1】下的【平面几何】。单击＜主菜单选项卡＞中的【几何】，此时在相应的＜主菜单功能选择区＞中单击【布尔操作和分割】，在其下拉菜单中选择【差集】，即可在【平面几何】下建立对象【差集1】。
对【差集1】进行设置：在＜参数设置区＞的【差集】设置区域，单击【要添加的对象】旁边的【激活】，使之切换为“on”状态，之后选择“r2”(在＜模型可视化区＞中单击【矩形r2】)即可；然后，单击【要减去的对象】旁边的【激活】，之后选择“r4”“r5”“r6”“r7”，再单击＜参数设置区＞顶端的【构建选定对象】即可。
2)建立并集
单击“模型树”中【几何1】→【工作平面1】下的【平面几何】。单击＜主菜单选项卡＞中的【几何】，此时在相应的＜主菜单功能选择区＞中单击【布尔操作和分割】，在其下拉菜单中选择【并集】，即可在【平面几何】下建立对象【并集1】。
对【并集1】进行设置：在＜参数设置区＞的【并集】设置中，单击【输入对象】旁边的【激活】，使之切换为“on”状态，选择“r1”“dif1”和“r3”，之后取消选中【保留内部边界】复选框，再单击＜参数设置区＞顶端的【构建选定对象】即可。
右击【工作平面1】，在弹出的快捷菜单中选择“拉伸”，在＜参数设置区＞单击【距离】→【指定：到面的距离】，将【距离】设置为“150”，再单击＜参数设置区＞顶端的【构建选定对象】，即可获得图4.24所示的几何模型。
图4.24 三维几何对象(书后附彩插)
至此，本例的三维几何对象建立完毕。
第3步：设置物理场属性
本例中，流质的设置不采用材料库中的【Water】，而是直接采用“空材料”定义动力黏度和密度。定义方法：＜模型开发器＞→【ComplexFlow.mph(root)】→【组件1】→【材料】，右击【材料】，在弹出的快捷菜单中选择【空材料】；在＜参数设置区＞的【几何实体选择】中选择所有域，设置【材料属性明细】中的【密度】(变量rho)为“Rhos”，设置【黏度】(变量mu)为“Mus”。
第4步：设置边界条件
＜模型开发器＞→【ComplexFlow.mph(root)】→【组件1】→【层流】，右击【层流】，分别单击【入口】、【出口】，添加边界条件，建立的对象为【入口1】和【出口1】。
单击【入口1】，在＜参数设置区＞可对入口边界进行设置。在＜模型可视化区＞的【边界选择】窗口中，单击x=0处端面(编号为1)，添加入口边界。在【边界条件】的下拉菜单中选择【压力】，并在下面的【压力条件】中输入“Pin”。类似地，单击【出口1】，然后在＜模型可视化区＞单击“x=1200”处的边界端面，该端面的序号为28。出口边界默认为压力边界条件，数值处输入“Pout”。其他边界面默认为壁面，不用进行设置。
第5步：网格划分
单击“模型树”中的【网格1】，在＜参数设置区＞进行网格大小设置。然而，如果按照上例的方式，采用自动网格设置，会出现一些问题。例如，采用“常规”网络大小划分的网格，在边界和拐角处的网格会比较密，而在通道内部的网格又会比较粗，显然，默认的网格划分模式不能较好地满足需求。在这种情况下，需要在＜参数设置区＞的【网格设置】→【序列类型】中选择【用户控制网格】，选定后发现在“模型树”的【网格1】下面出现了一些新的对象，主要对其中的【大小】进行设置。
单击【大小】，进入＜参数设置区＞，在【单元大小】中选择【定制】，将【最大单元大小】设置为“15”(通道宽度的10%)，将【最小单元大小】设置为“5”，其他设置不变。最后，单击＜参数设置区＞顶端的【全部构建】，即可完成网格划分。
第6步：配置求解器
此例不涉及多场耦合，属于求解层流流动的基本情况，求解器的选择及相关参数采用默认值即可。
第7步：求解
求解方法：＜模型开发器＞→【ComplexFlow.mph(root)】→【研究】，右击其下选项【步骤1：稳态】，在弹出的快捷菜单中选择【计算选定步骤】，即可开始计算。
第8步：后处理
本例的后处理任务有两个：任务1，显示出通道纵截面的速度分布；任务2，找到通道中速度最大值出现的位置和数值大小。
任务1：单击【结果】→【速度】→【切面1】，进入＜参数设置区＞。在【平面数据】中设置【平面类型】为“快速”，设置【平面】为“xy平面”，设置【定义方法】为“平面数”，设置【平面数】为“1”，表示只选择一个显示平面。设置完毕后，单击＜参数设置区＞上端的【绘制】按钮，即可得到在z=75处(通道中央)的速度场。
任务2：单击【结果】→【速度spf】，在＜主菜单功能选择区＞找到【更多绘图】，在其下拉菜单中选择【体图】→【体最大值/最小值】，即可在【结果】→【速度spf】下建立【体最大值/最小值1】的对象。单击该对象，进入＜参数设置区＞，在【高级】中设置如下：将【单元细化】设置为“4”(进一步细化内部网格，提高最大值的定位精度)；将【显示精度】设置为“3”(保留3位有效数字)；将【显示】选择“最大值”；选中【包含单位】复选框；在【着色和样式】中，将【颜色】选择“白色”；将【背景色】选择“黑色”。最后，单击＜参数设置区＞上端的【绘制】，即可在原速度场上增加关于最大值的显示。
速度截面和最大值显示如图4.25所示。
图4.25 速度截面和最大值(书后附彩插)
例4.8 基于电泳原理的十字通道样品分离模型二维模拟。
进样，是微流控芯片中典型的样品操作步骤。设置本例的主要目的是让大家不仅了解电泳进样的基本原理，而且理解多场耦合的概念。因此，本例采用二维模型，以减少在三维建模和计算过程中的开销。本例涉及的三维模型在COMSOL Multiphysics中几乎所有版本的案例库中都包含，有兴趣的读者可以在COMSOL官方网站(http：//cn.comsol.com)下载可供学习的建模文档，模型名称为Transport in an Electrokinetic Valve，或者在COMSOL Multiphysics软件的案例库中找到该文档。本例除了需要用到流体流动模块以外，还涉及电场及溶质的对流扩散，因此是一个典型的多场耦合模型。由于溶质的电泳是一个时变过程，因此需要用到瞬态求解器。此外，本例还涉及多个时序操作步骤，对各步骤的初始值设置需要有清晰的思路，这也是初学者容易出错之处，是本模型中的一个主要难点。
1)几何模型
本模型需要建立的通道几何形状如图4.26所示。所需的主要参数：缓冲液为水；在室温20℃下操作；样品浓度为1 mol/m3；一个样品离子带2个正电荷；样品的扩散系数为10×10-12 m2/s。
2)样品分离步骤
样品分离主要通过两个步骤来实现。为了对不同步骤进行的仿真有所区别，在此将每个步骤的设置、计算和分析称为一个“仿真模式”。
第1步，夹流聚焦(仿真模式1)。如图4.27所示，样品从a口随着流动进入通道，同时b、d口往内部输入缓冲液，c口是唯一的出口。通过控制a口或b、d口的流动速度，就可以调节“聚焦”样品的分布形状。图中指示的聚焦区用不同的颜色表示，是一个位于通道交叉口的一个锥形区域。
图4.26 进样通道结构与尺寸
图4.27 夹流聚焦示意图(书后附彩插)
图4.28 液体流动停止，依靠电泳分离聚焦区域样本(书后附彩插)
第2步，电泳进样(仿真模式2)。在这一步，一方面需要撤销a、b、c、d四个边界的速度边界，流动被禁止；另一方面要布置电场，溶质在合适的静电场作用下，在流体中发生定向迁移，即完全依靠电泳的作用实现样品的迁移。图4.28展示了样品在电场作用下一个时刻的迁移结果，为了较好地将聚焦区的样品分离出来，电场的设置需要一定技巧。在本例中，设置图4.28所示的电势值，所形成的电场在向下迁移聚焦区样本的同时，左侧通道的样品以一定速度向左迁移，右侧通道的样品以一定的速度向右迁移。这样分离出的样品才会独立分布于一段缓冲液中。此外，在迁移过程中，由于溶质的扩散作用，下侧通道中被分离出来的样品浓度将随着时间的推进而下降，相应地其扩散的区域会逐渐增大。
为了定量显示分离效果，可以在分离通道的某个纵截面观察浓度随时间的分布。这一步操作将在后面的演示中完成。
需要说明的是，尽管两个样品的操作步骤对应两个不同的仿真模式，但其中有大量建模步骤都是相同的，只有少量参数设置、求解步骤和分析方法有区别，因此在后续的建模流程中，除了特别说明某个设置对应特定的“仿真模式”以外，其他设置均为通用设置。
3)建模流程
【新建】→【模型向导】→【二维】，进入物理场选择界面，需要选择3个物理场。选择【流体流动】→【单相流】→【层流spf】，单击【添加】，即可将【层流spf】导入【添加的物理场接口】。类似地，添加【AC/DC】→【电流ec】；添加【化学物质传递】→【稀物质传递tds】。
接下来，要添加【研究】，单击【研究】→【稳态】，然后单击【完成】，进入开发模式。需要说明的是，由于求解的问题是一个时变问题，因此在后续操作中将加入【瞬态】求解器。
保存此工程文件为SampleSeparation.mph。下面介绍具体的建模过程。
1.“仿真模式1”的建模过程
第1步：设置变量与参数
所需设置的一些物理场参数如表4.3所示。在＜模型开发器＞→【SampleSeparation.mph(root)】→【全局定义】→【参数】处进行设置。
表4.3 模拟参数设置
续表
第2步：构建几何模型
操作步骤：＜模型开发器＞→【Samp leSeparation.mph(root)】→【组件1】→【几何1】，右击【几何1】，在弹出的快捷菜单中选择【矩形】，进入＜参数设置区＞，将【长度单位】设置为“μm”。
在＜主菜单选项卡＞单击【几何】，进入＜主菜单功能选择区＞，单击【矩形】，在＜模型可视化区＞按住左键拖动，即可建立一个矩形。此操作重复一次，即可建立第二个矩形。
设置矩形1的参数如下：
【大小和形状】宽度：1 000；高度：50。
【位置】→【居中】x：100；y：0。
设置矩形2的参数如下：
【大小和形状】宽度：50；高度：1 000。
【位置】→【居中】x：0；y：-100。
建立两个矩形后，需要合并两个矩形为一个整体的十字通道。操作方法：在＜主菜单功能选择区＞选择【布尔操作和分割】→【并集】，即可在【几何1】中建立并集对象【并集1】。单击＜模型可视化区＞中的两个矩形对象，即可在＜参数设置区＞设置【输入对象】r1和r2。取消选中【保留内部边界】复选框，然后单击＜参数设置区＞顶端的【构建选定对象】，即可完成十字通道的建立。
通过组合键“Ctrl+S”保存本工程文件。
第3步：设置物理场属性
设置【材 料】。操 作 方 法：＜模 型 开 发 器＞→【SampleSeparation.mph(root)】→【组件1】→【材料】，右击【材料】，在弹出的快捷菜单中选择【从库中添加材料】，转至＜自定义菜单区＞，选择【液体气体】→【Liquids】→【Water】，双击【Water】即可将其添加至“模型树”中的【材料】项。此时单击【材料】中的【Water】，即可在＜参数设置区＞看到与水相关的属性。需要说明的是，由于本例中施加了电场，因此水的材料属性中自动增加了【相对介电常数】和【电导率】这两个选项，这里的【相对介电常数】需要手动设置，在20℃(293 K)温度下，水的相对介电常数为80，直接输入“80”即可。
在添加【材料】之后，需要设置各物理场的参数。
(1)设置溶质的物理特性。＜模型开发器＞→【SampleSeparation.mph(root)】→【组件1】→【稀物质传递(tds)】→【传递属性1】。在＜参数设置区＞，将【温度】设置为在【参数】中设置的“T”。对流扩散现象描述了溶质在流动作用下的空间迁移效应，因此在溶质自然扩散的同时，要考虑流体的运载作用。在【对流】→【速度场】→【u】处，选择【速度场(spf)】，这一步操作很重要，主要是将浓度场和流动耦合，读者可以通过分析偏微分方程的形式来加深对这一步操作的理解。单击【稀物质传递(tds)】→【传递属性1】，查看＜参数设置区＞中的【方程】可见：
式中，J=-DΔc为菲克定律[1]；D为扩散系数(m2/s)；ci为物质组分的体积浓度(mol/m3或kg/m3)； Δ为梯度；“-”号表示扩散方向为浓度梯度的反方向，即扩散组元由高浓度区向低浓度区扩散。扩散通量J的单位是mol/m2·s或kg/m2·s。 Δ·J+u·Δc=R称为对流扩散方程，这里的u是流动的速度场，R为单位时间、单位体积空间内因化学反应生成物质组分的量[6]，本例不考虑化学反应。
通过分析式(4.117)可知，应明确浓度场和流场之间的关系，以及扩散系数对描述浓度变化迁移的影响。单击【扩散】→【材料】→【Water(mat1)】；在【扩散系数】→【Dc】处，将【用户定义】设置为“D”即可。
(2)设置电场属性。单击【电流ec】→【电流守恒1】，查看在＜参数设置区＞的参数设置，除了将默认的【温度】设置为“T”之外，其他参数值不改动。
(3)设置流动属性，单击【层流sp f】→【流体属性1】，查看在＜参数设置区＞的参数设置，除了将默认的【温度】设置为“T”之外，其他参数值不改动。
在各物理场设置下，【初始值1】均为默认设置，不作改动。
第4步：设置边界条件
“仿真模式1”只涉及【稀物质传递(tds)】和【层流(spf)】两个物理场的耦合作用，【电流(ec)】在“仿真模式2”才激活，因此可以暂时不改动。
首先，设置【稀物质传递(tds)】的边界条件。右击【稀物质传递(tds)】，重复添加3个【浓度】边界，再添加1个【通量】边界。例4.6和例4.7已介绍过如何设置边界条件，因此这里只陈述设置要求，具体操作步骤省略。将边界【浓度1】设置为a口处流入溶液的溶质浓度，浓度值为“Cin”，将边界【浓度2】和边界【流入3】分别设置为b、d口的缓冲液浓度，其值均设置为“0”。注意要选中【浓度】→【物质“c”】前面的复选框，方能设置浓度值。将边界【通量1】设置为c口的流出边界。
其次，设置【层流(sp f)】的边界条件。右击【层流(sp f)】，重复添加3个【入口】边界，再添加1个【出口】边界。将边界【入口1】设置为a口的速度边界，将边界【入口2】设置为b口的速度边界，将边界【入口3】设置为d口的速度边界。以上3个边界的速度已经在【参数】中设置，因此只需输入相应变量即可。将边界【出口1】设置为c口的零压力边界。
第5步：网格划分
由于求解域中涉及细长通道及拐角，因此默认的网格划分疏密程度不能满足正常的计算需求。在此需要用户自定义划分。
单击“模型树”中的【网格1】，在＜参数设置区＞进行网格划分。在【网格设置】→【序列类型】中，选择【用户控制网格】。这时在“模型树”中的【网格1】下，增加了一些可编辑对象。单击【大小】，在＜参数设置区＞激活网格设置。【单元大小】→【定制】：将【最大单元格大小】设为“5”，这是因为通道的宽度为50 μm，最大单元格边长为通道宽度的10%；将【最小单元格大小】设置为“1.5”；将【最大单元增长率】设置为“1.2”；将【曲率因子】设为“0.3”；将【狭窄区域分辨率】设置为“1”。然后，单击＜参数设置区＞顶部的【全部构建】。可以在＜模型可视化区＞查看网格划分结果。
第6步：配置求解器
在＜模型开发器＞→【Samp leSeparation.mph(root)】→【研究1】中进行配置。单击【研究1】→【步骤1：稳态】，在＜参数设置区＞的【研究设置】中，将物理场【电流(ec)】的复选框取消选中，意为在求解中不考虑电场作用。
第7步：求解
求解方法：＜模型开发器＞→【Samp leSeparation.mph(root)】→【研究1】，右击其下选项【步骤1：稳态】，在弹出的快捷菜单中选择【计算选定步骤】，即可开始计算。
第8步：后处理
单击【结果】→【速度(spf)】→【表面】，即可显示速度场，如图4.29(a)所示；单击【结果】→【浓度(tds)】→【表面】，即可显示浓度场，如图4.29(b)所示。
图4.29 通道内速度场与浓度场(书后附彩插)
(a)速度场；(b)浓度场
其中，在浓度场中的通道交叉区域得到了锥形的样本聚焦区，这样的样品分布即样本分离步骤(仿真模式2)的初始值。
最后，按组合键“Ctrl+S”保存本工程文件。
2.“仿真模式2”的建模过程
“仿真模式2”是在“仿真模式1”的基础上，实现样品的电泳迁移，其中绝大部分的设置沿用“仿真模式1”中的设置。
第1步：设置参数与变量
无变动设置。
第2步：构建几何模型
无变动设置。
第3步：设置物理场参数
由于样品的迁移只需用到电泳作用，因此需要消除流动产生的溶质迁移效应。具体方法：＜模型开发器＞→【SampleSeparation.mph(root)】→【组件1】→【稀物质传递(tds)】→【传递属性1】，在这里主要设置溶质的物理特性。在＜参数设置区＞，修改【对流】→【速度场】→【u】的设置，选择【用户定义】，将x和y方向的速度均设置为“0”，表示无载流作用。
此外，由于要施加电场作用，因此需要激活【稀物质传递(tds)】中电泳迁移功能。操作方法：单击＜模型开发器＞→【SampleSeparation.mph(root)】→【组件1】→【稀物质传递(tds)】，在＜参数设置区＞，选中【传递机理】→【附加传递机理】下的【电场迁移】。此时，单击【稀物质传递(tds)】下的【传递属性1】，可在＜参数设置区＞发现多了【电场迁移】的设置选项。在【电场迁移】→【电势】“V”处选择【电势(ec)】，将【电荷数】“zc”设置为“Nzc”。
注意：在激活【电场迁移】后，单击【稀物质传递(tds)】→【传递属性1】，查看＜参数设置区＞中的【方程】可见，原来的菲克定律描述Ji=-Di
Δci改变为Ji=-Di Δci-zi um，i Fci
ΔV，公式中多了一项：-zi um，i Fci ΔV。其中，zi为离子电荷数；um，i为离子迁移率；F为外力项；ci为浓度， ΔV为电势梯度。这样的改变，就将电泳导致的离子迁移作用加入了。
第4步：设置边界条件
＜模型开发器＞→【SampleSeparation.mph(root)】→【组件1】→【电流(ec)】，分别添加4个电 势边界，在a、b、c、d四个边界分 别 设 置“Va”“Vb”“Vc”“Vd”。
找到＜模型开发器＞→【SampleSeparation.mph(root)】→【组件1】→【稀物质传递(tds)】，首先禁用原来a、b、c、d四个边界的设置。操作方法：右击相应边界对象，在弹出的快捷菜单中选择【禁用】即可。然后，建立4个通量边界，分别设为a、b、c、d的四个边界。在＜参数设置区＞，选中【物质“c”】复选框，然后在“J0，c”输入“-tds.nm flux_c”，这是一个内置参数，主要用于描述跨边界的自由通量。
第5步：网格划分
无变动设置。
第6步：配置求解器
在这一步，需要用到瞬态求解器，以在一定的时间段内完成溶质的电泳迁移。操作方法：单击＜主菜单选项卡＞中的【研究】，在＜主菜单功能选择区＞单击【添加研究】，此时在＜自定义菜单区＞显示要添加研究的种类，这里选择【一般研究】→【瞬态】，双击【瞬态】即可完成瞬态求解器的添加。在“模型树”中，可以看到新的求解器对象【研究2】→【步骤1：瞬态】。
继续配置瞬态求解器。【研究2】→【步骤1：瞬态】，＜参数设置区＞→【研究设置】：将【时间步】设置“range(0，0.1，10)”，意为计算时间为0～10 s，时间步长为0.1 s；选中物理场【电流(ec)】复选框，意为考虑电场作用。
接下来，设置溶质浓度的初始值，即【研究1】中的浓度分布结果。找到【因变量值】→【求解变量的初始值】，这里的求解变量包括溶质浓度和电场。依次选择【设置】→【用户控制】、【方法】→【解】、【研究】→【研究1：稳态】、【选择】→【自动】，这样就将溶质浓度的初始值设置为“仿真模式1”中的计算结果了。
第7步：求解
求解方法：＜模型开发器＞→【Samp leSeparation.mph(root)】→【研究2】，右击其下的选项【步骤1：瞬态】，在弹出的快捷菜单中选择【计算选定步骤】，即可开始计算。
第8步：后处理
注意：对于“仿真模式1”下的求解结果，在【结果】中也是保留的，在后处理中要进行区分。例如，【浓度(tds)】存储的是“仿真模式1”下的求解结果，而【浓度(tds)1】存储的是“仿真模式2”下的求解结果。
单击【结果】→【浓度(tds)1】→【表面】，即可显示浓度场。在＜参数设置区＞找到【数据】，将【时间(s)】选择为“10”，然后单击＜参数设置区＞顶端的【绘制】，即可得到如图4.30所示在浓度场演化时间为10 s的浓度分布情况。
图4.30 样品分离效果图(书后附彩插)
此外，为了定量观测某个时刻在分离通道中心处的浓度分布变化，需要用到“截线显示”方法。具体操作：单击【结果】→【浓度(tds)1】→【表面】，在＜主菜单功能区＞单击【选择】→【截线起点】，这时就会发现＜模型可视化区＞的浓度图上多了一条红线，之后在“十字”通道中央单击一次，即设置了截线的起点；然后，单击【选择】→【截线终点】，在“十字”通道的下端分离通道的中线处单击一次，这样就会发现红线穿越了分离的样本区域，并且在【结果】中添加了一个新的对象【一维绘图组8】，单击【结果】→【一维绘图组8】，在＜参数设置区＞→【数据】→【时间选择】→【手动】区域，在【时间索引(1-101)】处输入“range(1，10，101)”，意为每隔10个时间步(1 s)显示一次浓度分布。最后，单击＜参数设置区＞顶端的【绘制】，结果如图4.31所示。图中的10条曲线从左到右依次为1 s、2 s、…、10 s时的浓度分布曲线。
图4.31 分离通道中线处浓度动态分布曲线
可以看出，随着时间的推移，溶质样品逐渐被分离，并且浓度的峰值逐渐降低，展宽增加。这说明在电泳迁移时，溶质也在发生扩散现象。
例4.9 微混合器优化设计案例。
将两种溶液中的溶质相互混合，是微流控芯片中经常要进行的操作。流体在微尺度下通常以层流方式流动，在这种情况下，溶质通常仅能通过扩散来实现混合。通过这样的方式来实现充分混合存在两个问题：其一，要耗费足够长的时间等待溶质扩散，这极大地降低了微流控芯片的工作效率；其二，两种溶液之间需要较大的接触面积，这不利于微流控芯片的微型化。
为了快速地在有限空间内实现溶质的充分混合，就需要改进其混合过程，设计专门的混合单元——微混合器。按照混合机理的差异，微混合器分为主动式和非主动式两种。主动式微混合器主要通过增加外力作用(如气、磁、声、电等作用)来引导溶质的迁移或对层流进行干扰。增加外力作用来引导溶质的迁移，本质上是采用外界作用来提高溶质的扩散效率；若对层流进行干扰，本质上则是在单位时间、单位空间内增加了两种流体之间的混合面积，从而促进扩散。被动式微混合器主要通过改变通道形状来增加两种溶质在单位空间内的接触面积，常见的有T形、L形和人字形微混合器等。主动式微混合器的混合效率一般比较高，但需要增加外围设备以施加外力作用，并且容易出现局部发热、产生气泡等负效应；被动式微混合器的混合效率相对低一些，但既不需要外围设备，也不用担心外力作用带来的负效应，因此成本要低一些，结构和工况也要简单一些。
本例主要介绍一种被动式微混合器的设计，通过改变通道中的结构设计来对比研究提出方案的混合效率。本例共提出三种通道结构设计，通道的结构和尺寸如图4.32所示。
图4.32 通道结构设计
各通道中央处的黑线是区域分界线，而不是无滑移壁面，其不影响流动，用于初始化时设置不同的溶液类型。通过初始化不同类型溶液的分布，可以减少从初始值演化到稳态混合过程的计算步骤。通道的左侧是流动入口，右侧为出口。在入口中点设置分隔点，用于区分入口处流入的溶液类型。入口上半侧设置的是缓冲液(浓度为0)，下半侧设置的是浓度为“Cin”的溶液。为提高有限空间内的混合效率，在结构2和结构3中，设置一些列规律分布的“挡块”结构。其中，结构2和结构3中挡块的倾斜角度不同。为了比较不同通道结构的混合效率，可直接观测在达到稳态混合的情况下出口处的浓度分布，易知通道上下两侧的浓度差异越小就代表混合越充分。作为一个完全理想的混合结果，出口位置从上到下的浓度分布均无差异，但事实上不可能得到这个结果。下面来完成这个实例的建模设计。
采用导航流程【新建】→【模型向导】→【二维】，进入“物理场”选择界面，选择两个物理场：选择【流体流动】→【单相流】→【层流spf】，单击【添加】即可将【层流sp f】导入【添加的物理场接口】；类似地，添加【化学物质传递】→【稀物质传递tds】。
接下来，要添加【研究】。单击【研究】→【瞬态】，然后单击【完成】，进入开发模式。
保存此工程文件为MicroM ixer.mph。下面进入具体的建模过程。
第1步：设置参数与变量
涉及的参数设置如表4.4所示。
表4.4 模拟参数设置
续表
第2步：构建几何模型
操作步骤：＜模型开发器＞→【MicroMixer.mph(root)】→【组件1】→【几何1】，将＜参数显示区＞中的【长度单位】设置为“μm”。
采用3个【矩形】对象构建主通道。其尺寸和布局参照表4.4中的要求设置。在控制通道之间的距离时，【位置】→【基】选择【角】，将“x”均设为“0”，使三个通道左边对齐。对于通道2，将“y”设为“Ch_W+Dis_C”；对于通道3，将“y”设为“2*Ch_W+2*Dis_C”。这样，即可确定三个通道的位置和尺寸。
建立7个【矩形】对象，设置通道2内的挡块。采用参数中的设置来控制挡块的大小和位置。按照挡块的序号顺序，关键参数如表4.5所示。
表4.5 通道2挡块位置设置
建立完成后，采用布尔运算方法，用通道2减去上述7个挡块的区域，即可得到所需的通道结构。具体步骤已在前面的例子演示过，这里不再赘述。
然后，建立7个【矩形】对象，设置通道3内的挡块。按照挡块的序号顺序，关键参数如表4.6所示。
表4.6 通道3挡块位置设置
建立完成后，采用布尔运算方法，用通道3减去上述7个挡块的区域，即可得到所需的通道结构。
接下来，绘制通道内部的“分界线”。操作方法：单击＜主菜单选项卡＞→【几何】→【绘图】→【线】，单击相应图标。以作通道1的分界线为例，将鼠标指针移至＜模型可视化区＞，当指针接近通道1入口边界的中点时，指针的十字聚焦中心将会“磁力对准”该中点，单击，将指针移至通道1出口中点，再次单击，即可确定需要的分界线线段。最后，右击，完成操作。采用同样的方法，可以完成通道2、通道3中的分界线设置。
通过按组合键“Ctrl+S”保存本工程文件。
第3步：设置物理场属性
设置【材料】。操作方法：＜模型开发器＞→【MicroMixer.mph(root)】→【组件1】→【材料】，右击【材料】，在弹出的快捷菜单中选择【从库中添加材料】，转至＜自定义菜单区＞，选择【液体气体】→【Liquids】→【Water】，双击【Water】即可将其添加至“模型树”中的【材料】项。在＜参数设置区＞，要保证所设置的材料在【所有域】均适用。
在添加【材料】之后，需进行各物理场的参数设置。【层流(spf)】的【流体属性1】采用默认设置即可。在【稀物质传递(tds)】的【传递属性1】中，需要将【对流】→【速度场】→【u】→【速度场(spf)】；将【扩散】→【扩散系数】→【Dc】设置为“D”。
右击【稀物质传递(tds)】，在弹出的快捷菜单中选择【初始值】，增加对象【初始值2】，单击【初始值2】，在＜参数设置区＞选择所有通道的下半区域，在【初始值】→【浓度】→【c】中输入“Cin”。这样就在各通道设置了具有特定分布的浓度初始值，每个通道的下半侧初始分布为待稀释溶液，而上半侧则是浓度为0的缓冲液。
第4步：设置边界条件
【稀物质传递(tds)】的边界条件设置。右击【稀物质传递(tds)】，先添加2个【浓度】边界，再添加1个【通量】边界。其中，浓度边界【浓度1】为缓冲液，施加在三个通道左侧入口的上半侧边界。在＜参数设置区＞选中【物质“c”】复选框，在【c0，c】中输入“0”。【浓度2】为待稀释溶液的浓度，施加于三个通道左侧入口的下半侧边界。在＜参数设置区＞选中【物质“c”】复选框，在【c0，c】中输入“Cin”。对于出口设置，选择三个通道的右侧边界即可，设置边界类型为“通量”。所设置的边界条件可被三个通道公用，只要在＜参数设置区＞的【边界选择】中添加相应的通道边界即可。
【层流(spf)】的边界条件设置。右击【层流(spf)】，先添加2个【入口】边界，再添加1个【出口】边界。在【入口1】设置“Vin_1”，在＜参数设置区＞的【边界选择】中选择所有缓冲液入口边界(通道入口的上半侧)，在【入口2】设置“Vin_2”，在＜参数设置区＞的【边界选择】中选择所有待稀释溶液入口边界(通道入口的下半侧)。对【出口1】，在＜参数设置区＞的【边界选择】中选择所有通道的右侧出口边界，默认的“0”压力出口无须改动。
第5步：网格划分
由于求解域中涉及细长通道及拐角，因此默认的网格划分疏密程度不能满足正常的计算需求。在此需用户自定义划分。
单击“模型树”中的【网格1】，在＜参数设置区＞进行网格划分。在【网格设置】→【序列类型】中，选择【用户控制网格】。这时，在“模型树”中的【网格1】下增加了一些可编辑对象。单击【大小】，在＜参数设置区＞激活网格设置。【单元大小】→【定制】：将【最大单元格大小】设为“12”，最大单元格边长为通道宽度的6%；将【最小单元格大小】设为“10”；将【最大单元增长率】设为“1.2”；将【曲率因子】设为“0.4”；将【狭窄区域分辨率】设为“1”。然后，单击＜参数设置区＞顶部的【全部构建】，即可在＜模型可视化区＞查看网格划分结果。
以上划分主要针对求解域内部区域，在划分结构中发现，边界层附近的网格很密，网格疏密差异太大，对此可以通过【网格1】下的对象【大小1】来专门对边界层附近的网格进行自定义修改。操作方法：单击【大小1】，进入＜参数设置区＞，可以看到【几何实体选择】→【几何实体层】中默认为【边界】，所有壁面边界都已经被作为被选择的对象。在【单元大小】下选择【定制】，勾选下面的所有定制项，设置【最大单元格大小】为“12”；将【最小单元格大小】设为“8”；将【最大单元增长率】设为“1.2”；将【曲率因子】设为“0.4”；将【狭窄区域分辨率】设为“1”。然后，单击＜参数设置区＞顶部的【全部构建】，即可重划网格。观察网格疏密情况，可发现边界层附近的网格与求解域内部的网格大小基本匹配。
第6步：配置求解器
在＜模型开发器＞→【M icroMixer.mph(root)】→【研究1】中进行配置。单击【研究1】→【步骤1：瞬态】，在＜参数设置区＞的【研究设置】中，将【时间步】设置为“range(0，0.1，5)”，意为流场演化时间为0～5 s，时间步长为0.1 s。其他设置无须改动。
第7步：求解
求解方法：＜模型开发器＞→【MicroMixer.mph(root)】→【研究1】，右击其下选项【步骤1：瞬态】，在弹出的快捷菜单中选择【计算选定步骤】，即可开始计算。
第8步：后处理
在计算完成后，如果在【结果】中没有关于“速度”“浓度”和“压力”的绘图组，那么可以进行手动添加。例如，计算完成后，发现没有“压力”绘图组，那么可以通过如下步骤添加：＜模型开发器＞→【MicroMixer.mph(root)】，右击【结果】，在弹出的快捷菜单中选择【二维绘图组】，此时发现在【结果】中增加了一个【二维绘图组×】对象(这里的×代表一个数值，根据已有的绘图组数量自动编号)。然后右击【二维绘图组×】，将其重命名为“压力”，再右击【压力】，在弹出的快捷菜单中选择【表面】，即可生成一个二维的表面对象【表面1】。然后，进入＜参数设置区＞，将【数据】→【数据集】选择为【研究1/解1(sol1)】，找到【表达式】，在【表达式】同一栏的右侧，有两个下拉菜单按钮，任意单击一个按钮可展开下拉菜单，找到【模型】→【组件1】→【层流】→【速度和压力】→【P-压力-Pa】，双击【P-压力-Pa】即可自动设置表达式。在＜参数设置区＞顶端单击【绘制】，即可显示压力云图。
(1)观察速度分布，如图4.33所示。
从速度分布上看，通道1中的速度分布相对平均一些，通道2中的局部速度差异有一定的加大，通道3中挡块顶端速度有大幅提高。由此可见，挡块的存在让局部流动的速度分布变得具有很大的差异。
图4.33 流场演化时间为5 s时的速度场分布(书后附彩插)
(2)观察压力分布。如图4.34所示。
图4.34 流场演化时间为5 s时的压力场分布(书后附彩插)
从压力分布来看，随着流动的进行，压力逐渐降低。由于三个通道的入口和出口的速度分布都是相同的，所以各通道的流量相同。此外，由于出口边界均为0压力设置，因此可知产生通道1所需的压力差是最小的，而通道3需要最大的压力差来保证产生相同的流量。
(3)观察浓度分布。
由于是层流扩散，因此从初始值到稳态分布有一个过渡过程，而我们需要观察最后的稳态分布，所以需要判断出口处浓度的时变情况。操作方法：＜模型开发器＞→【M icroMixer.mph(root)】，右击【结果】，在弹出的快捷菜单中选择【一维绘图组】，建立一个【一维绘图组×】对象；右击【一维绘图组×】，可在弹出的快捷菜单中选择【上移】或【下移】来调整一维绘图组所在的上下位置。右击【一维绘图组×】，在弹出的快捷菜单中，选择【点图】，在【一维绘图组×】下面建立一个点图对象【点图1】，单击【点图1】，在＜参数设置区＞选择要监测的点，然后单击【选择】→【手动】，选择3个通道出口的中点。在【y轴数据】同一栏的右侧，有两个下拉菜单按钮，任意单击一个按钮可展开下拉菜单，找到【模型】→【组件1】→【稀物质传递】→【物质c】→【通量】→【c-浓度-mol/m3】，双击【c-浓度-mol/m3】即可自动设置表达式。在＜参数设置区＞顶端单击【绘制】，即可显示选择的3个观察点浓度随时间的变化，如图4.35所示。
图4.35 3个观察点的浓度变化曲线(书后附彩插)
可见，在2 s之后，三个通道的混合就已经达到稳定了，这说明所设置的计算时间范围可满足观察混合达到稳态后浓度分布的需求。
浓度分布情况如图4.36所示。浓度在出口位置上的分布可以反映通道结构对混合效果的影响，这是本设计研究的核心内容。通过色度图不易定量比较差异，因此可以用“截线显示”的方法来观察3个通道入口处和出口处的浓度分布。需要建立两条截线，并将其绘制在一幅图上对比显示。“截线显示”的方法在例4.8中已经介绍过，这里不再赘述，只给出在流场演化时间为5 s时的浓度分布。
图4.36 流场演化时间为5 s时的浓度分布(书后附彩插)
图4.37左侧曲线为通道1的入口处和出口处的浓度分布，中间曲线为通道2的入口处和出口处的浓度分布，右侧曲线为通道3的入口处和出口处的浓度分布。其中，点划线描述的是入口处的浓度分布，3个通道入口处的浓度分布相同；实线代表出口附近的浓度分布。由其分布规律可见，通道1基本没有实现混合功能，通道2的出口处浓度的差异有减小的趋势，通道3的出口处浓度的差异最小，说明混合效果是最好的。
图4.37 三个通道出口处的稳态浓度分布(书后附彩插)
这个例子表明，在通道流量相同的情况下，较大的通道压力差会导致较大的局部速度差异，也可增强溶质的混合效率。
例4.10 基于弹性瓣膜的流体定向输送装置设计。
在微流控芯片中，流体的定向输送是流动操控的重要内容。在低雷诺数下(Re＜40)，流动产生的原因主要是压力差。因此，通过控制流动演化过程中的压力差就可以控制流动的方向。本部分设计通过设置弹性瓣膜，利用其流固耦合机理，在驱动端施加往复运动的流动模式来实现主通道中的流体输送。本例借鉴了COMSOL案例库中的“微泵机理”，在这个案例基础上进行进一步的机理分析。由于在主通道中没有涉及复杂的流动驱动结构，因此能够最大限度地降低堵塞风险，并保护细胞、组织等生物材料不受到极端压力和高剪切力的损害。流体输送装置的结构和主要尺寸如图4.38所示。
图4.38 微泵结构与机理示意图(书后附彩插)
整个通道的尺度在微米级，驱动端和主通道结合形成一个倒“T”形的结构，两个初始倾斜安置的阀瓣分布固定于主通道下侧，两个阀瓣中间的区域与驱动端通道直接连通，形成一个整体区域。
当驱动端流体流入时，两个阀瓣处于不同的工作状态：一方面，左侧阀瓣的右下侧压强增大，会迫使其往左运动，并且有减小通道流动截面的趋势，这样会形成某种程度的“液封”，产生的阻力会减少流体向左侧流动；另一方面，右侧阀瓣的左上侧压强增大，会迫使其往右运动，并且有增大通道流动截面的趋势，这样会减小流动阻力，增加流体向右侧流动的趋势。
当驱动端流体流出时，左侧阀瓣的右下侧压强减小，会迫使其往右运动，并且有增大通道流动截面的趋势，产生的阻力会减小，有利于流体的通过。右侧阀瓣的左上侧压强减小，会迫使其往左运动，有减小通道流动截面的趋势，这样也会形成“液封”，减少液体回流。
以上两个过程交替进行，虽然在驱动端的流动有进有出，但在主通道形成的净流量是向右的，从而形成了流体的定向输送。该装置不但可以用于溶液的空间转移，而且可用于构建具有连续回路的循环系统来为微电子器件散热。
这个案例涉及流体和弹性体之间的力学耦合问题，其设置、求解都有一定的代表性，并且与微流控芯片设计密切相关，因此被作为一个典型的教学案例供读者研究和学习。
下面来完成这个实例的建模设计。
采用导航流程【新建】→【模型向导】→【二维】，进入“物理场”选择界面后，选择【流体流动】→【流-固耦合】→【流-固耦合】，单击【添加】，即可将【流-固耦合】导入【添加的物理场接口】。
然后要添加【研究】。由于本例主要研究对象的时变特征，因此选择【研究】→【瞬态】，然后单击【完成】，进入开发模式。
保存此工程文件为MicroPump.mph。下面进入具体的建模过程。
第1步：设置参数与变量
涉及的参数如表4.7所示。
表4.7 模型参数设置
第2步：构建几何模型
(1)进入几何设置模式：＜模型开发器＞→【MicroPump.mph(root)】→【组件1】→【几何1】，在＜参数显示区＞将【长度单位】设置为“μm”。
(2)建立1个【矩形】对象设置主通道。矩形的参数设置：
【大小和形状】宽度：L；高度：H。
【位置】→【基】→【角】x：0；y：0。
(3)建立1个【矩形】对象设置左阀瓣。矩形的参数设置：
【大小和形状】宽度：2*rp；高度：2*hp。
【位置】→【基】→【居中】x：x0-hp*sin(beta)/2；y：0。
在【旋转角度】→【旋转】输入“-beta”。
(4)建立主通道的矩形副本。操作方法：在＜主菜单选项卡＞单击【几何】，在＜主菜单功能选择区＞单击【变换】→【复制】，即可建立1个复制对象，默认名为“复制1(copy1)”。在＜参数设置区＞选择“r1”为复制对象，然后在＜参数设置区＞顶端单击【构建选定对象】。在此，将r1的复制对象作为主通道几何。
(5)采用布尔运算方法，获得阀瓣的几何区域。操作方法：在＜主菜单功能选择区＞单击【布尔操作和分割】，在其下拉菜单中选择【交集】，即可在【平面几何】下建立对象【交集1(int1)】。对【交集1】进行设置，在＜参数设置区＞的【交集】设置中，单击【要添加的对象】旁边的【激活】，使之切换为“on”状态，之后选择“r1”“r2”。取消勾选【清除内部边界】，再单击＜参数设置区＞顶端的【构建选定对象】，即可构建。这里由于r1的副本存在，所以在【输入对象】不能直接在＜模型可视化区＞选择“r1”，此时可在【输入对象】列表显示栏右击【粘贴选择】，在弹出的文本框中输入“r1”即可。
(6)在获得阀瓣的几何区域【交集1(int1)】后，需对阀瓣的顶端进行圆角处理。在＜主菜单选项卡＞单击【几何】，在＜主菜单功能选择区＞中单击【操作】→【圆角】。在＜参数设置区＞选择要制造圆角的点，也就是Int1对象的第3、4点，在【半径】处输入“rp”，在＜参数设置区＞顶端选择【构建选定对象】即可。
(7)复制1个左侧的阀瓣，构建右侧阀瓣。在＜主菜单选项卡＞中单击【几何】，在＜主菜单功能选择区＞展开【变换】，单击【复制】，建立一个新对象，名为“复制2(copy2)”。在＜参数设置区＞选择【圆角1(fil1)】为复制对象。然后在【位移】→【x】输入“2*x0”。在＜参数设置区＞顶端选择【构建选定对象】，即可完成右侧阀瓣的几何构建。
(8)建立1个【矩形】对象，设置驱动端通道。矩形的参数设置：
【大小和形状】宽度：0.8*H；高度：H。
【位置】→【基】→【角】x：0.43*L；y：H。
在＜参数设置区＞顶端选择【构建选定对象】即可。
到此，几何模型建立完毕，按组合键“Ctrl+S”保存本工程文件。
第3步：设置物理场属性
1)设置材料
＜模型开发器＞→【MicroPump.mph(root)】→【组件1】，右击【材料】，在弹出的快捷菜单中选择【空材料】。在【材料】设置窗口中，在【标签】文本框中输入“流体”。在＜参数设置区＞，【密度】→【rho】，设置为“dens”；【动力黏度】→【mu】，设置为“visc”；其他无关项输入“0”。
采用上述同样操作，建立一个【空材料】对象，命名为“固体”。在＜参数设置区＞，将【几何实体选择】设置为“2，4”；【密度】→【rho】，设置为“rho_fb”；【杨氏模量】→【E】，设置为“E_fb”；【泊松比】→【nu】，设置为“Psr_fb”；其他无关项输入“0”。
2)设置其他属性
添加两个积分组件耦合算子，将其中一个添加到形成左侧出口的边界，将另一个添加到形成右侧出口的边界。这两个算子将用于计算流出每个出口的流率。
(1)设置入口积分算子。＜模型开发器＞→【MicroPump.mph(root)】→【组件1】→【定义】，右击【定义】，在弹出的快捷菜单中选择【组件耦合】→【积分】可建立【积分1(intop1)】，将其命名为“intopL”。在＜参数设置区＞，【源选择】→【几何实体层】→【边界】，将边界设置为“1”。
重复以上步骤，设置出口(边界17)积分算子，将其命名为“intopR”。
(2)创建几个变量来计算通道中的净流率。其中，变量UoutL、UoutR分别为主通道左侧出口和右侧出口的流率(单位为m3/s)。计算方法分别为UoutL=-intopL(u_fluid)*W和UoutR=intopR(u_fluid)*W，其中u_fluid为通道流动横截面处与界面垂直的速度分量。请注意符号规则：为从左向右的流动指定正值，为从右向左的流动指定负值。变量UoutNet用于计算每个出口的平均流率差，UoutNet=UoutR-UoutL，正值对应于从左向右的净流量。
设置方法：＜模型开发器＞→【MicroPump.mph(root)】→【组件1】，右击【定义】，在弹出的快捷菜单中单击【变量】可建立【变量1】对象，定位到＜参数设置区＞，分别输入“UoutL”“UoutR”“UoutNet”三个变量及相应表达式。
为了对一段时间内的平均净流率求积分来计算泵送的流体体积，需要在模型中添加一个“全局常微分和微分代数方程”物理场接口。添加方法：＜模型开发器＞→【MicroPump.mph(root)】→【组件1】，右击，在弹出的快捷菜单中选择【添加物理场】，定位到＜自定义菜单区＞，在“模型树”中选择【数学】→【常微分和微分代数方程接口】→【全局常微分和微分代数方程(ge)】，双击【添加】即可。
单击【组件1】→【全局常微分和微分代数方程(ge)】→【全局方程1】，在＜参数设置区＞设置【名称】为“Vpump”，设置【f(u，ut，utt，t)(1)】为“Vpumpt-UoutNet”。这里实际上是设置的一个微分方程，主要用于描述从左向右泵送的净体积V pump，其中V pump为要求解的自变量，V pumpt为V pump的一阶导数(二阶导数则为V pumptt)。之后，定位到【单位】栏，定义因变量V pump的单位。【单位】→【因变量物理量】→【定制单位】，输入“m^3”；【单位】→【源项物理量】→【定制单位】，输入“m^3/s”。
3)配置流-固耦合物理场接口
最后，是定义变形域。＜模型开发器＞→【MicroPump.mph(root)】→【组件1】→【定义】→【动网格】，单击【变形域1】，在＜参数设置区＞的【域选择】中选择“1”，在【平滑处理】→【硬化因子】下的“C2”文本框中输入 “100”。
第4步：设置边界条件
(1)设置变形域边界。【组件1】→【定义】，右击【动网格】，在弹出的快捷菜单中选择【固定边界】，建立对象【固定边界1】。在＜参数设置区＞的【边界选择】中选择“1”“2”“7”“9”“16”和“17”，这些边界在变形域发生变形时，是固定不动的。【组件1】→【定义】，右击【动网格】，在弹出的快捷菜单中选择【指定法向网格位移】，建立对象【指定法向网格位移1】，在＜参数设置区＞的【边界选择】中选择“3”和“12”，表示该边界上的节点只做法向运动。在【指定法向网格位移】的“dn”处输入“0”，其作用相当于固定，只是指定方式不一样。
(2)设置流动边界。＜模型开发器＞→【MicroPump.mph(root)】→【组件1】→【层流(sp f)】，单击【层流(spf)】，进入＜参数设置区＞，在【域选择】选择“1”和“3”(如果有“2”和“4”，则在列表中将其选中，单击列表右侧的“-”即可将其删除)。右击【层流(sp f)】，添加1个【入口】边界。在＜参数设置区＞的【边界选择】中，选择“10”，即驱动端入口。在【速度】→【法向流入速度】的“U0”处输入“U*6*s*(1-s)*sin(2*pi*t/(1[s]))”，注意：这里的“6*s*(1-s)”本质上是一个抛物线分布的速度，s是从0到1变化的一个一维空间算子；“sin(2*pi*t/(1[s]))”是一个正弦函数，描述了流动连续、周期性地改变方向。与上述建立入口边界的方法类似，建立一个出口边界，选择主通道两端的边界“1”和“17”，默认的0压力出口无须改动，注意到这里的出口边界可能存在“倒灌”情况，因此要取消【抑制回流】复选框中的勾选。
(3)设置阀瓣的固体力学边界。＜模型开发器＞→【MicroPump.mph(root)】→【组件1】→【固体力学(solid)】，单击【固体力学(solid)】，进入＜参数设置区＞，在【域选择】中选择“1”“2”和“4”(如果存在“3”，则将其删除)。这里域1是变形域，它虽然在列表中，但不起作用。在【固体力学(solid)】节点下，单击【线弹性材料1】，在【线弹性材料】的设置窗口中，定位到【线弹性材料】栏，从【使用混合公式】列表中选择【压力公式】，这表示采用流体压力作为阀瓣的外力。最后，设置阀瓣的固定端。右击【固体力学(solid)】，在弹出的快捷菜单中选择【固定约束】，单击新建的【固定约束】对象，在＜参数设置区＞的【边界选择】中，选择“4”和“13”两个边界。注意：这两个边界很小，不易用鼠标点选，可以单击“粘贴选择”图标，在弹出的对话框中直接输入。
到此，边界设置完成，按组合键“Ctrl+S”保存本工程文件。
第5步：网格划分
为了建立合适的网格分布，这里需要自定义网格大小和分布。单击“模型树”中的【网格1】，在＜参数设置区＞进行网格划分。在【网格设置】→【序列类型】中，选择【用户控制网格】。这时在“模型树”中的【网格1】下，增加了一些可编辑对象。右击【网格1】，可以设置若干【大小】对象，在此建立3个“大小”对象，其中【大小】指建立域内网格，【大小1】指除了阀瓣网格之外的边界网格，【大小2】指阀瓣边界(不含固定边界)的网格。具体【定制】参数设置如表4.8所示。
表4.8 网格参数设置
设置完毕后，单击＜参数设置区＞顶部的【全部构建】即可。
第6步：配置求解器
配置方法：在＜模型开发器＞→【MicroPump.mph(root)】→【研究1】中进行配置。单击【研究1】→【步骤1：瞬态】，在＜参数设置区＞的【研究设置】中，将【时间步】设置为“range(0，0.02，2)”；在【相对容差】中输入“0.001”。
第7步：求解
求解方法：＜模型开发器＞→【MicroPump.mph(root)】→【研究1】，右击其下选项【步骤1：瞬态】，在弹出的快捷菜单中选择【计算选定步骤】，即可开始计算。
第8步：后处理
(1)显示动画。单击＜主菜单选项卡＞中的【结果】，定位到动画，选择【文件】，可以在【结果】→【导出】添加【动画1】对象。在＜参数设置区＞将【目标】选择为【播放器】，在【帧】中选择【全部】。之后在【结果】→【导出】节点下右击【动画1】，在弹出的快捷菜单中选择【播放】，即可观看动画。动画演示了阀瓣对流体做出反应时的被动运动是如何导致流体从左向右流动的。在驱动端流体“下行”(流入)过程中，右侧阀瓣朝通道底部弯曲，左侧阀瓣朝远离底部的方向弯曲。相对于右侧出口来说，朝左侧出口的流动受到限制。在驱动端流体“上行”(流出)过程中，两个阀瓣向相反方向弯曲，从右侧口向通道中心的向内流动受到限制。这将导致从左向右的净流动，即流体从左侧口吸入，从右侧口被推出。
(2)显示净流量的演化曲线。在【结果】节点下，单击【一维绘图组4】，在一维绘图组的设置窗口中，在【标签】文本框中输入“净输送体积”。定位到【图例】栏，从【位置】列表中选择【左上角】，即可得到如图4.39所示的曲线。该曲线展示了在主通道中输送净流体积随时间变化的情况。
图4.39 从左到右的净流体积
例4.11 液滴操控案例。
微液滴技术是在微尺度通道内，利用流动剪切力与表面张力之间的相互作用，将连续流体分割成离散的纳升级(及以下)体积的液滴的一种微纳技术。它是近年来发展起来的一种全新的操纵微小液体体积的技术[7]。微液滴类型主要有气-液相液滴和液-液相液滴两种。气-液相液滴由于容易在微通道中挥发和造成交叉污染而限制了其应用。液-液相液滴根据连续相和分散相的不同，分为水包油(O/W)、油包水(W/O)等。液-液相微液滴有体积小、液滴样品间无扩散、可避免样品间的交叉污染、反应条件稳定、适当操控下可实现迅速混合等优点，是一种十分理想的微反应器，已经被用于研究微尺度条件下的众多反应及其过程，并在化学和生命科学等领域拓展出了重要应用[8]。目前，液滴技术已广泛应用于DNA、蛋白质、酶等生物大分子的分析检测，以及药物传递等生物医学领域。
液滴微流控芯片设计的关键在于通道结构和通道尺寸的设计。根据液滴生成方式的不同，可以将微通道的结构分为T形微通道、流动聚焦型微通道、毛细管共轴型微通道、阶梯型微通道等。在这些通道中，T形微通道具有结构简单、生成的液滴大小易控等特点，是液滴操控中最基本的应用。本部分主要完成T形微通道生成液滴示例，演示多相流建模方法(水平集方法，Level Set)的应用。通道结构与各相流速设置如图4.40所示。
图4.40 通道结构与各相流速设置(书后附彩插)
美国数学家Stanley Osher和James Sethian于20世纪80年代开发出了水平集(Level Set，LS)方法。这一方法最初主要应用于智能控制、图像处理、材料微细加工等领域，经进一步发展后用于两相流动的数值模拟。从几何角度上看，气-液(或液-液)相界面对应的是一个高阶曲面与一个低阶曲面的交界线(或交界面)；从数学角度上看，气-液(或液-液)相界面对应于一个高阶方程的零等值面。使用方法追踪气液相界面位置时，需要在整个计算区域内定义一个光滑的连续相函数Φ；两种流体的交界面可用相函数的零值点表示，相界面表示如图4.41所示，通过在计算区域上求解一个一阶偏微分方程来更新相界面的位置。
图4.41 相界面(书后附彩插)蓝色-油相；红色-水相
LS方法中，随着界面运动，封闭界面内的流体质量会损失，即造成质量不守恒。为解决这个问题，需使用特殊方法对计算结果进行修正。此外，水平集函数通常不保持为演化界面的有向距离函数，这是该方法的不足。
尽管LS方法的应用存在以上难点和不足，但其优点也是显而易见的。该方法属于相界面捕捉方法的范畴，通过相函数的等值点可以方便地追踪相界面位置，从而避免方法中相界面构造的复杂运算；由于相函数空间分布连续，且具有距离函数的特性，因此可以精确、方便地计算出相界面的曲率以及与曲率相关的物理量，从而解决该方法难以精确计算相界面曲率及与曲率相关的物理量的问题[9]。
在COMSOL Multiphysics中，除了LS方法，还提供了相场(Phase Field，PF)法，这两种方法均可以用于描述同一种多相流现象。下面采用LS方法来完成本实例的建模设计。
首先，采用导航流程【新建】→【模型向导】→【二维】，进入“物理场”选择界面，选择【流体流动】→【多相流】→【两相流，水平集】→【层流两相流，水平集】，单击【添加】，即可将其导入【添加的物理场接口】。
然后，添加【研究】。由于本案例主要研究对象的时变特征，因此要选择【研究】→【包括相初始化的瞬态】，然后单击【完成】，进入开发模式。
保存此工程文件为Droplet-LS.mph。下面进入具体的建模过程。
第1步：设置参数与变量
本次建模涉及的参数设置如表4.9所示。将这些参数输入＜模型开发器＞→【MicroPump.mph(root)】→【全局定义】→【参数】。
表4.9 模拟参数设置
续表
注意：
(1)相界面宽度epsilon。该参数描述了相界面两侧相函数的过渡区宽度。epsilon越大，界面就越宽，相界面之间的色度分隔就越模糊；epsilon越小，界面之间的色度分隔就越明显，但是太小容易引起计算不稳定。根据经验，epsilon一般取通道宽度中最窄宽度的2%～5%为宜。
(2)重新初始化系数r。重新初始化是LS方法中非常重要的步骤。相界面函数Φ在经历几个时间步的运动后，将不再维持距离函数的性质，这会导致相函数Φ分隔的区域出现质量不守恒，重新初始化就是要解决这个问题，以确保水平集函数始终近似为距离函数。根据经验，r一般取流场中的最大速度值。
(3)接触角θ。其定义可参考图4.42。
(4)滑移长度β。这是描述润湿壁面的参数，太大和太小都会改变相界面和壁面之间的作用规律。根据经验，这个参数在数值上可取求解域中网格边长的平均尺度。滑移长度β的具体定义可参考图4.43。
图4.42 接触角
图4.43 滑移长度
第2步：构建几何模型
(1)进入几何设置模式：＜模型开发器＞→【Droplet-LS.mph(root)】→【组件1】→【几何1】，在＜参数显示区＞的【长度单位】选择“μm”。
(2)建立1个【矩形】对象构建主通道。矩形的参数：
【大小和形状】宽度：L1；高度：W。
【位置】→【基】→【角】x：0；y：0。
(3)建立1个【矩形】对象设置侧方通道。矩形的参数：
【大小和形状】宽度：W；高度：L2。
【位置】→【基】→【角】x：L2；y：-L2。
(4)采用布尔运算方法，将主通道和侧方通道合并。在＜主菜单功能选择区＞中单击【布尔操作和分割】，在其下拉菜单选择【并集】，即可在【平面几何】下建立对象【并集1(int1)】。对【并集1】进行设置，在＜参数设置区＞的【并集】设置中，单击【要添加的对象】旁边的【激活】，使之切换为“on”状态，之后选择“r1”“r2”。选中【保留内部边界】复选框，再单击＜参数设置区＞顶端的【构建选定对象】，即可进行构建。
第3步：设置物理场属性
(1)设置多向流计算参数。＜模型开发器＞→【Droplet-LS.mph(root)】→【组件1】→【水平集(ls)】，单击下面的节点【水平集模型1】，在＜参数设置区＞的【域选择】中，选择“1”和“2”。在【水平集模型】→【重新初始化参数】的“γ”下输入“r”，在【界面厚度控制参数】的“εls”后输入“epsilon”。
(2)设置多向流物理参数。定位到＜模型开发器＞→【Drop let-LS.mph(root)】→【组件1】→【多物理场】，单击其下的【两相流，水平集1(tpf1)】，在＜参数显示区＞→【流体1的属性】→【流体1】→【流体1的密度】的“rho1”下选择【用户定义】，输入“rho_w”，【流体1的动力黏度】的“μ1”下选择【用户定义】，输入“niu_w”；【流体2的属性】→【流体2】→【流体2的密度】的“rho2”下选择【用户定义】，输入“rho_o”，【流体2的动力黏度】的“μ2”下选择【用户定义】，输入“niu_o”。在【表面张力】栏选中【包含表面张力】复选框，在激活的【表面张力系数】中选择【用户定义】。在【用户定义】下面的【表面张力系数】后输入“delta”。
这样，就完成了耦合流场及物理参数的设置。
第4步：设置边界条件
(1)设置层流边界。＜模型开发器＞→【Droplet-LS.mph(root)】→【组件1】→【层流(spf)】，右击【层流(spf)】，建立2个【入口】边界和1个【出口】边界。单击【入口1】边界对象，在＜参数设置区＞的【边界选择】中选择“1”(即主通道入口)，【速度】→【法向流入速度】，在“U0”处输入“Vo”；单击【入口2】边界对象，在＜参数设置区＞的【边界选择】中，选择“5”(即主通道入口)，【速度】→【法向流入速度】，在“U0”处输入“Vw”。对【出口1】，在＜参数设置区＞的【边界选择】中，选择主通道的右侧出口边界“9”，默认的0压力出口无须改动。
(2)设置水平集边界。＜模型开发器＞→【Droplet-LS.mph(root)】→【组件1】→【水平集(ls)】，右击【水平集(ls)】，建立2个【入口】边界、1个【出口】边界和1个【初始值】对象。单击【入口1】边界对象，在＜参数设置区＞的【边界选择】中，选择“1”，在【水平集条件】选择【流体2(Φ=1)】；单击【入口2】边界对象，在＜参数设置区＞的【边界选择】中，选择“5”，在【水平集条件】选择【流体1(Φ=0)】。对【出口1】，在＜参数设置区＞的【边界选择】中，选择主通道的右侧出口边界“9”即可。
单击【初始值1】对象，在＜参数设置区＞的【域选择】中，选择“1”，在【初始值】→【水平集变量】中选择【指定相】，然后选择【流体2(Φ=1)】。这表示在主通道中，初始时充盈流体2。
单击【初始值2】对象，在＜参数设置区＞的【域选择】中，选择“2”，在【初始值】→【水平集变量】中选择【指定相】，然后选择【流体1(Φ=0)】。这表示在侧通道中，初始时就充盈了流体1。
(3)设置多物理场边界。＜模型开发器＞→【Droplet-LS.mph(root)】→【组件1】→【多物理场】，右击【多物理场】，在弹出的快捷菜单中选择【润湿壁】。单击建立的对象【润湿壁1(ww1)】，在＜参数设置区＞的【边界选择】中选择“2”“3”“4”“7”“8”。在【润湿壁】→【接触角】的“θw”处输入“pi-ceta”，在【滑移长度】的“β”处输入“beta”。
第5步：网格划分
单击“模型树”中的【网格】，在＜参数设置区＞进行单元格大小设置。由于本例中通道的结构简单，因此选择【网格设置】→【序列类型】→【物理场控制网格】→【单元大小】→【细化】，单击＜参数设置区＞顶部的【全部构建】，即可实现网络划分。
第6步：配置求解器
＜模型开发器＞→【Droplet-LS.mph(root)】→【研究1】→【步骤2：瞬态】，在＜参数设置区＞的【研究设置】中，在【时间步】设置“range(0，0.1，40)”。在【求解时显示结果】下选中【绘制】复选框，在绘图组中选择【流体1的体积分数(ls)】。
第7步：求解
(1)＜模型开发器＞→【Droplet-LS.mph(root)】→【研究1】，右击下面的节点【步骤1：相场初始化】，在弹出的快捷菜单中选择【计算选定步骤】。
(2)＜模型开发器＞→【Droplet-LS.mph(root)】→【研究1】，右击下面的节点【步骤2：瞬态】，在弹出的快捷菜单中选择【计算选定步骤】，开始计算。
第8步：后处理
(1)调整色度显示。“模型树”→【结果】→【流体1的体积分数(ls1D)】，单击其下节点【流体1的体积分数】，在＜参数设置区＞→【着色和样式】→【着色】处选择【颜色表】，在【演示表】处选择【Thermal】，分别选中【颜色表反序】和【对称颜色范围】复选框，即可获得图4.44所示的液滴分离效果。
图4.44 液滴分离效果图(书后附彩插)
(2)显示动画。在＜主菜单选项卡＞中单击【结果】，定位到动画，选择【文件】，可以在【结果】→【导出】添加【动画1】对象。在＜参数设置区＞将【目标】选择为【播放器】，在【帧】中选择【全部】。之后在【结果】→【导出】节点下右击【动画1】，在弹出的快捷菜单中选择【播放】，即可观看液滴生成的动画。