本发明涉及微流控领域，更具体的说，它涉及终端封闭微通道中粒子移动耦合模型控制方法。
背景技术：
在20世纪，人们借助于电子在半导体或金属中流动得到的“信息”，成就了具有战略意义的信息科学技术；而在21世纪，通过带有可溶性生物分子或悬浮细胞的水溶液在移动控制芯片通道或平面上流动以研究生命、理解生命，以至于部分地改造生命，将有可能同样成就一种新的具有战略意义的科学技术：微流控学。微通道中粒子移动是微流控学研究的热点之一。而在终端封闭的微通道中，由于其一端封闭的几何形状是不允许有净流体流量产生，如何将粒子有效地传输到终端封闭的几何形状中对许多生命领域是非常重要的。在药物传输中，药物通过主血管通道传输，然而连接着的旁支(封闭支路)，都是一些需要药物治疗的重要部位，比如坏死的细胞区域。在大脑细胞外间隙中，由于不规则的细胞排列产生“死区”，分子在死区内运动不仅与细胞外间隙的扩散紧密相关，而且影响死区附近，细胞对分子的吸收。
众所周知，布朗运动很自然地作为一种粒子移动的机制，然而单靠布朗运动是非常低效的。压力驱动微流体运动，流体携带着粒子移动的方案并不可取，由于终端封闭的微通道中没有净流量，施加的流体压力不仅没能驱动粒子，反而有可能破坏终端封闭结构。电动驱动作为微流控学近十几年来主要研究热点之一，其在粒子移动方面有较成熟的研究，但在终端封闭这些较偏的“支路”结构中施加电场，显得力不从心。而通过化学的方法，设计粒子移动的化学条件，对于现有的科技提出很大的挑战。
技术实现要素：
本发明克服了现有技术的不足，提供了快速参数辨识、最优精准控制边界溶液进行移动扩散的终端封闭微通道中粒子移动耦合模型控制方法。
本发明的技术方案如下：
终端封闭微通道中荧光粒子移动耦合模型控制方法，包括模型通道和控制装置；所述模型通道包括主微通道和支路通道，所述主微通道与支路通道的一端连接固定，支路通道的另一端封闭，该主微通道的入口设置微泵控制溶液进入，支路通道为真空状态；所述溶液包括带荧光粒子的nacl溶液、nacl溶液和油滴；所述控制装置包括微泵、移动控制芯片和高速相机，微泵、高速相机与移动控制芯片连接，具体移动芯片处理步骤如下：
101)获得支路通道的灰度图步骤：通过高速相机获取支路通道的灰度图，荧光粒子的实际流动速度，该实际流动速度包括荧光粒子扩散电泳速度和流体运动速度两者之和；
102)获得荧光粒子浓度步骤：灰度值正比于荧光粒子浓度值，而获得i(x,z,t)∝n(x,z,t)，其中i为像素点的灰度值，n为高速相机拍摄时刻的荧光粒子浓度值，该灰度值、荧光粒子浓度值与x方向即x轴线方向和z方向即高度方向有关，t为时间变量；
103)获取荧光粒子溶液浓度变化步骤：在支路通道中，采用一维的扩散模型来获得荧光粒子溶液在n个基函数的变化叠加组成，其中n为10，具体公式如下：
其中c(x,t)为支路通道中的nacl溶液浓度，an＝2(1-cosλn)/λn和λn＝(2n+1)π/2，其将溶液浓度变化在时间上和空间上分离，dc为nacl双极性溶液的扩散系数，初始条件c(x,0)＝ci，其中ci为一常数；边界条件c(0,t)＝c0,
104)获得荧光粒子扩散电泳速度步骤：令c(x,z,t)＝c(x,t)，c(x,z,t)表示为荧光粒子的移动变化，此时c(x,z,t)与z方向无关，即在特定的x方向上，所有z方向的溶液浓度都相等；从而得到荧光粒子在溶液浓度梯度驱动下移动，产生的荧光粒子扩散电泳速度公式如下：
其中，γp为荧光粒子扩散电泳的系数，vp(x,z,t)为荧光粒子扩散电泳速度，此速度是在溶液浓度梯度驱动下产生的；
105)获得流体速度步骤：流体的速度为vf(x,z,t)＝[vxf(x,z,t),vzf(x,z,t)]·，具体流体在x方向z方向的两个维度上的速度如下：
其中γw为壁面滑移系数，h表示z方向上的高度，是常量；
106)获得荧光粒子移动速度步骤：由流体速度、荧光粒子扩散电泳速度得到荧光粒子移动的实际速度v(x,z,t)＝vp(x,z,t)+vf(x,z,t)，并令
v(x,z,t)＝[vx(x,z,t),vz(x,z,t)]·；从而建立封闭终端微通道中荧光粒子在x方向z方向上的移动公式如下：
其中dp是荧光粒子的扩散系数，其边界条件为如下公式：
n为荧光粒子在x＝0处的浓度值，为一常数；后三项为在x＝l,z＝h,z＝0处满足的自然边界条件；
107)最优策略步骤：设置荧光粒子在一定的时间内要达到的分布量，结合荧光粒子的相关数据，而采用近似梯度优化方法和形式梯度优化法来制定最佳溶液边界控制方法，所述荧光粒子的相关数据包括荧光粒子溶液浓度、荧光粒子浓度变化、荧光粒子扩散电泳速度、流体速度和荧光粒子移动速度。
进一步的，所述步骤107)具体的最优策略步骤如下：
701)确定控制形式步骤：移动控制芯片采用边界溶液浓度作为实际控制溶液移动扩散的控制量；
702)近似梯度优化步骤：由荧光粒子移动模型推演得到多个不同方向的优值，并比较所有的优值，获得荧光粒子移动相对最优值，而防止了优化方法陷入局部最优解；
703)形式梯度优化步骤：对相对最优值进行处理，具体处理为通过变分方法得到伴随模型，并结合有限元方法，获得荧光粒子移动模型，最终得到边界溶液最佳控制方法。
进一步的，确定控制形式的控制目标为：
nd(x,z)为荧光粒子浓度分布达到期望的分布。
进一步的，基于fenicsproject的近似梯度优化方法步骤如下：
1)将需要优化的量c0在猜测区间中分成n个，那么就有n组优化猜测值；
2)加一个扰动项系数β，取[0,0.05]之间的数；
3)计算每个目标函数之差
4)在第k+1步，采用平移变换算子
5)计算
6)最后，将n个最优值中选取最小值，并记录此时的优化量c0。
进一步的，形式梯度优化步骤的具体如下：采用变分方法，并基于fenicsproject，结合有限元方法，得到荧光粒子移动模型，荧光粒子移动模型为：
其中m＝1,2,3,...,m，m为[0,t]之间的离散数，δt＝t/m为时间间隔量，w(x,z)为fenicsproject设定的测试函数，满足第一边界条件时，其形成fenics形式的伴随模型如下：
本发明相比现有技术优点在于：
1、本发明精简巧妙的支路通道荧光粒子移动控制芯片拓扑结构设计，通过主微通道中油滴的注入，防止不同浓度的nacl溶液在主微通道中混合，从而在支路通道中产生溶液浓度梯度。对支路通道和主微通道采用真空处理，从而使初始的nacl溶液能自发地注满整个荧光粒子移动控制芯片。
2、本发明的终端微通道中荧光粒子移动耦合模型包括荧光粒子移动时空演化模型和溶液浓度时空演化模型，其中后者耦合与前者。分析物理过程，将溶液浓度模型弱化为一维扩散模型。后者采用二维荧光粒子移动时空演化模型，精准描述荧光粒子移动的过程。控制器设计在溶液浓度时空演化模型的边界上，而控制目标只跟荧光粒子移动时空演化模型有关。在控制器设计时，通过弱化溶液浓度时空演化模型为一维，从而只需考虑荧光粒子移动时空演化模型的控制器设计，无需考虑溶液浓度时空演化模型的。
3、本发明对于二维荧光粒子移动模型中的荧光粒子移动速度计算，考虑荧光粒子扩散电泳速度和流体运动速度，更精准地完成参数辨识工作。
4、本发明将二维荧光粒子移动模型的参数辨识问题，转化为两个子问题，一是荧光粒子移动速度优化问题，二是参数辨识问题。在荧光粒子移动速度优化问题处理中，通过设定荧光粒子移动性能指标和加入正则项，转化为一组代数优化问题，容易求解。
5、本发明在参数辨识问题求解中，考虑到跟荧光粒子扩散电泳系数，壁面滑移系数和nacl双极性溶液的扩散系数都有关。而只与壁面滑移系数和nacl双极性溶液的扩散系数相关。先用智能优化方法求解得到壁面滑移系数和nacl双极性溶液的扩散系数，然后用最小二乘方法直接得到得到荧光粒子扩散电泳系数。分别设定方向和方向的目标函数，加快参数辨识的求解。
6、本发明给定边界溶液浓度控制形式，不仅将控制轨迹优化问题转化为参数优化问题，而且只需考虑荧光粒子移动模型的变分形式，大大减少了伴随模型的复杂度。引入fenicsproject平台和有限元方法，给出二维荧光粒子移动模型和二维伴随模型的fenics形式，提升计算的稳定性和准确度。
7、本发明采用近似梯度优化方法与形式梯度优化方法结合的方式设计边界溶液浓度的控制器，尽可能地避免陷入局部最优解，得到全局最优解。在近似梯度优化方法中，结合随机和布点思想，尽可能找到次优解。由于近似梯度优化不够精确。将近似优化的结果当成形式梯度优化的初始值，精确的形式梯度求解达到最优值。
附图说明
图1为本发明支路通道荧光粒子移动控制芯片的拓扑结构图；
图2为本发明支路通道与主通道的连接图；
图3为本发明支路通道荧光粒子移动控制芯片的系统图；
图4为本发明支路通道荧光粒子移动控制芯片参数辨识部分的实验流程；
图5为本发明荧光粒子移动耦合模型的建模流程图；
图6为本发明荧光粒子移动耦合模型的参数辨识流程图；
图7为本发明荧光粒子移动耦合模型的控制设计流程图；
图8位本发明控制效果图和溶液浓度控制曲面。
具体实施方式
下面结合附图和具体实施方式对本发明进一步说明。
实施例一
如图1至图8所示，终端封闭微通道中粒子移动耦合模型控制方法，包括模型通道和控制装置；所述模型通道包括主微通道和支路通道，所述主微通道与支路通道的一端连接固定，支路通道的另一端封闭，该主微通道的入口设置微泵控制溶液进入，支路通道为真空状态；所述溶液包括带荧光粒子的nacl溶液、nacl溶液和油滴；其中油滴的注入是为了防止带荧光粒子的nacl溶液和nacl溶液混合；所述控制装置包括微泵、移动控制芯片和高速相机，微泵、高速相机与移动控制芯片连接，具体移动芯片处理步骤如下：
101)获得支路通道的灰度图步骤：通过高速相机获取支路通道的灰度图，荧光粒子的实际流动速度，该实际流动速度包括荧光粒子扩散电泳速度和流体运动速度两者之和；
102)获得荧光粒子浓度步骤：根据高速相机获取终端封闭微通道中的荧光粒子移动情况的灰度图，因灰度图的光强度正比于荧光粒子浓度值，因此形成i(x,z,t)∝n(x,z,t)，i为像素点的灰度值，n为高速相机拍摄时刻的荧光粒子浓度值，与x方向即x轴线方向和z方向即高度方向有关，t为时间变量；
103)获取荧光粒子溶液浓度变化步骤：非均匀的溶液浓度分布产生溶液浓度梯度，从而驱动荧光粒子移动，产生荧光粒子移动速度；在终端封闭微通道中，采用一维的扩散模型来进行表示，具体公式如下：
其中dc为nacl双极性溶液的扩散系数，c(x,t)为终端封闭微通道中的nacl溶液浓度，初始条件c(x,0)＝ci，其中ci为一常数；边界条件
由公式(1)可得到溶液浓度变化公式(2)如下：
其中an＝2(1-cosλn)/λn和λn＝(2n+1)π/2，其将溶液浓度变化在时间上和空间上分离；
104)获得荧光粒子扩散电泳速度步骤：由于荧光粒子移动是二维，而溶液浓度为一维的因此在荧光粒子移动二维模型中，令c(x,z,t)＝c(x,t)，此c(x,z,t)与z方向无关，即在特定的x方向上，所有z方向的溶液浓度都相等；从而得到荧光粒子在溶液浓度梯度驱动下移动，产生的荧光粒子扩散电泳速度公式(3)如下：
其中，γp为荧光粒子扩散电泳的系数，vp(x,z,t)为荧光粒子扩散电泳速度，此速度是在溶液浓度梯度驱动下产生的；
105)获得流体速度步骤：因高速相机拍摄的是荧光粒子扩散电泳速度和流体运动速度两者之和；而流体所产生的速度是与荧光粒子扩散电泳速度相反，阻碍粒子进入终端封闭微通道中，而流体所产生的速度为
其中p(x,z,t)为流体的压力，第一项为斯托克斯方程，第二项表示流体压力跟z方向无关，且公式(4)满足公式(5)的边界条件；
其中γw为壁面滑移系数，在终端封闭管道中，流体产生的净流量为零，即
从而得到：
106)获得荧光粒子移动速度步骤：从上述公式(3)、公式(7)和公式(8)可得粒子移动的实际速度v(x,z,t)＝vp(x,z,t)+vf(x,z,t)，并令
v(x,z,t)＝[vx(x,z,t),vz(x,z,t)]·；从而建立封闭终端微通道中二维荧光粒子移动公式如下：
其中dp是粒子的扩散系数，该公式(9)的边界条件为：
其中，n为粒子在x＝0处的浓度值，为一常数；公式(10)的后三项表明，在x＝l,z＝h,z＝0处满足自然边界条件；
107)最优策略步骤：根据步骤103)获得的参数，并设置粒子在一定的时间内要达到的分布量，而采用近似梯度优化方法和形式梯度优化法来制定最佳溶液边界控制方法。
实施例二
如图1至图8所示，支路通道中荧光粒子移动耦合模型控制方法，包括模型通道、控制装置和观测装置；所述模型通道包括主微通道、支路通道和溶液进入通道，所述控制装置包括微泵、移动控制芯片，所述观测装置包括显微镜物镜、滤色片、高速相机、激光器。移动控制芯片设置在支路通道上，并且所述激光器设置的位置通过激发的光源，产生的光路通过显微镜物镜到达移动控制芯片，高速相机监测支路通道中二维荧光粒子浓度信息，从而建立微控模型。
具体包括如下步骤：
101)建立模型步骤：所述主微通道与支路通道中未封闭的一端连接固定，该主微通道的一端与溶液进入通道连接，并在进入通道处设置微泵。所述支路通道为真空状态。
102)注入溶液步骤：对微控模型通过微泵注入带荧光粒子的nacl溶液、nacl溶液和油滴，其中油滴的注入是为了防止带荧光粒子的nacl溶液和nacl溶液混合。
具体以上两步骤如图1所示为支路通道荧光粒子的拓扑结构，依次为荧光粒子1，带荧光粒子的nacl溶液2，主微通道3、微泵4、油5、油滴6、支路通道7、nacl溶液8等组成。荧光粒子为直径为1.01μm的荧光聚苯乙烯乳胶荧光粒子。移动控制芯片材料为聚二甲基硅氧烷。支路通道尺寸为200μm×10μm×10μm。两个微泵型号为tj-3a/w0109-1b，适用注射器5μl～60ml，线速度范围7.94μm/min-79.4mm/min。带荧光粒子的nacl溶液2中溶液浓度为0.02mmol/l(任意合理浓度皆可)，荧光粒子浓度为1×10-3μm-3(任意合理浓度皆可)。nacl溶液8中溶液浓度为2mmol/l(任意合理浓度皆可)。油滴的注入是为了防止带荧光粒子的nacl溶液2和nacl溶液8混合。如此设计可以使不同浓度的nacl溶液进入主微通道，从而为参数辨识和控制的实际输入做准备。支路通道的一端封闭，另一端与主通道相连接。部分主微通道和支路通道要注满2mmol/lnacl溶液，因此支路通道中初始情况必须为真空。首先将部分主微通道和支路通道抽成真空，然后将2mmol/lnacl溶液注入真空的部分主微通道和支路通道中。当然，一般在支路通道荧光粒子移动控制芯片的下端设置实验容器，用于储存实验过程中流过主通道中的溶液和荧光粒子。支路通道与主通道的连接如图2所示。
具体如图3所示：移动控制芯片的系统图，主要包括支路通道荧光粒子移动控制芯片9，显微镜物镜11，滤色片12，高速相机13，光路14，激光器15。激光器激发光源，产生光路14，通过显微镜物镜11到达支路通道荧光粒子移动控制芯片9。高速相机13测量得到支路通道中的二维荧光粒子浓度信息(聚苯乙烯荧光粒子在激光下会产生荧光，通过对荧光强度的获取，测量二维荧光粒子浓度信息)，为荧光粒子移动耦合模型的参数辨识以及控制设计提供必要的数据。
终端微通道中荧光粒子移动耦合模型包括二维荧光粒子移动时空演化模型和一维溶液浓度时空演化模型。即二维荧光粒子移动时空演化模型描述荧光粒子在支路通道中的移动情况，一维溶液浓度时空演化模型描述溶液浓度在支路通道中的演化情况。
103)参数辨识步骤：通过高速相机获取微控模型中的支路通道中的荧光粒子移动情况，从而获得相应溶液和粒子的运动参数；具体步骤如下：
301)获得荧光粒子浓度步骤：根据高速相机获取支路通道中的荧光粒子移动情况的灰度图，因灰度图的光强度正比于荧光粒子浓度值，因此形成i(x,z,t)∝n(x,z,t)，i为像素点的灰度值，n为高速相机拍摄时刻的荧光粒子浓度值，与x方向即x轴线方向和z方向即高度方向有关，t为时间变量。
即在具体试验中在激光器等实验设备的辅助下，利用高速相机拍摄支路通道中的荧光粒子移动情况，得到时间间隔非常短的两帧灰度图。设拍到的灰度图像素点为100*16。根据所拍的微通道范围和高速相机的性能有关，此实验采用100*16，这个设定可以根据实际所需进行更改。从中因灰度图的光强度近似正比于荧光粒子浓度值，即i(x,z,t)∝n(x,z,t)。i为像素点的灰度值，n为高速相机拍摄时刻的荧光粒子浓度值，与x方向(x为轴线方向)和z方向有关(z为高度方向)，t为时间变量，并考虑到支路通道尺寸为高度2h＝10μm，长度l＝200μm。由于支路通道中高度方向的荧光粒子移动是中心对称的，已z轴中线为界，因此考虑0μm到5μm之间的区域即可。
302)获取荧光粒子移动速度步骤：非均匀的溶液浓度分布产生溶液浓度梯度，从而驱动荧光粒子移动，产生荧光粒子移动速度；在支路通道中，采用一维的扩散模型来进行运算，是因为为了建立终端微通道中荧光粒子移动耦合模型，首先考虑溶液浓度时空演化模型，而在支路通道中，溶液浓度在z方向的几乎没有变化，而荧光粒子的移动并不会对溶液浓度产生影响。因此具体公式如下：
其中dc为nacl双极性溶液的扩散系数，c(x,t)为支路通道中的nacl溶液浓度，初始条件c(x,0)＝ci，其中ci为一常数；边界条件
由公式(1)可得到溶液浓度变化公式(2)如下(采用10维来逼近真实解)：
其中an＝2(1-cosλn)/λn和λn＝(2n+1)π/2，其将溶液浓度变化在时间上和空间上分离，从而大大简化后续参数辨识和控制器设计的难度。
303)获得荧光粒子扩散电泳速度步骤：由于荧光粒子移动是二维，而溶液浓度为一维的因此在荧光粒子移动二维模型中，令c(x,z,t)＝c(x,t)，此c(x,z,t)与z方向无关，即在特定的x方向上，所有z方向的溶液浓度都相等；从而得到荧光粒子在溶液浓度梯度驱动下移动，产生的荧光粒子扩散电泳速度公式(3)如下：
其中，γp为荧光粒子扩散电泳的系数，vp(x,z,t)为荧光粒子扩散电泳速度，此速度是在溶液浓度梯度驱动下产生的；因高速相机拍摄的是荧光粒子扩散电泳速度和流体运动速度两者之和；而流体所产生的速度是与荧光粒子扩散电泳速度相反，阻碍荧光粒子进入支路通道中，而流体所产生的速度为
其中p(x,z,t)为流体的压力，第一项为斯托克斯方程，第二项表示流体压力跟z方向无关，且公式(4)满足公式(5)的边界条件；
其中γw为壁面滑移系数，在终端封闭管道中，流体产生的净流量为零，即
从而得到：
304)获得荧光粒子移动速度步骤：从上述公式(3)、公式(7)和公式(8)可得荧光粒子移动的实际速度v(x,z,t)＝vp(x,z,t)+vf(x,z,t)，并令v(x,z,t)＝[vx(x,z,t),vz(x,z,t)]·；从而建立封闭终端微通道中二维荧光粒子移动公式如下：
其中dp是荧光粒子的扩散系数，该公式(9)的边界条件为：
其中n为荧光粒子在x＝0处的浓度值，为一常数；公式(10)的后三项表明，在x＝l,z＝h,z＝0处满足自然边界条件。
具体在步骤103)的应用中采用如下步骤进行具体化参数辨识，结合实际调试设置的值，设两帧之间的时间间隔为δτ，在高速相机快速地采集过程中，两帧时间间隔非常短，因此两帧之间的荧光粒子移动速度v(x,z,t)＝[vx(x,z,t),vz(x,z,t)]·可认为不变，记为
式中，ω为的已知量为100*16的像素点区域，n(x,z,t)为高速相机拍摄的前后两帧灰度值，其中前一张的荧光粒子浓度灰度值可以另为n(x,z,p)，后一张的荧光粒子浓度灰度值可以另为n(x,z,a)。目标函数(11)的物理含义是拍摄到的两帧的图像灰度值满足二维荧光粒子移动模型。易得，如果非常满足时，此目标函数等于零，即达到最优状态。dp和
像素点100*16上的灰度值为一个矩阵，令第i行第j列的像素点(xi,zj)表示为(i,j)，i＝2,3,5,4...,99，j＝2,3,4,...,15(其中，第100行和第1行，第1列和第16列可以通过边界条件得到)。公式(10)的荧光粒子浓度相关的值可以用两帧之间的灰度值来得到
其中，δx和δz为像素点在x方向和z方向上的间距。公式(12)-(14)可以通过两帧的灰度值直接得到，继而离散的目标函数(11)变为
离散的目标函数(14)对未知数
式中，κ为一个设定的常数，
根据公式(3)(7)(8)，我们得到荧光粒子移动的实际运动速度求解公式
式中，t＝0为荧光粒子进入微通道的时刻。τ为相对于t＝0时刻的拍摄荧光粒子移动的时间。
完成荧光粒子移动模型参数辨识需要两步骤，由于nacl双极性溶液的扩散系数dc存在指数函数中，而vx，vz中又含有对数函数，因此对此的参数辨识需要采用智能优化算法。考虑到vx跟荧光粒子扩散电泳系数γp，壁面滑移系数γw和nacl双极性溶液的扩散系数dc都有关。而vz只与壁面滑移系数γw和nacl双极性溶液的扩散系数dc相关。我们首先求解z方向上的目标函数
采用智能优化算法求解得到壁面滑移系数γw和nacl双极性溶液的扩散系数dc。然后求解x方向上的目标函数
用最小二乘方法直接得到得到荧光粒子扩散电泳系数γp。
104)最优策略步骤：根据步骤103)获得的参数，并设置荧光粒子在一定的时间内要达到的分布量，而采用近似梯度优化方法和形式梯度优化法来制定最佳溶液边界控制方法。即控制的目的是希望荧光粒子能够在给定的时间内尽可能达到一定的分布，(因为荧光粒子是大量的，希望得到一定分布的荧光粒子)。其中将采用近似梯度优化方法和形式梯度优化法，考虑到荧光粒子移动模型是二维模型，引入feincsproject，结合有限元方法，而得到荧光粒子移动模型，并且采用近似梯度优化法和形式梯度优化法，能有效的避免找到优化控制结果陷入局部最优点。
1、设置控制器形式，采用边界溶液浓度作为实际的控制量，且此方法在移动控制芯片上较容易实现。
根据公式(1)给出时空演化过程中的控制形式如下：
其中，初始的溶液浓度是个已知量ci，而边界溶液浓度c0是未知量。
控制目标为：
目标函数的意义就是通过边界溶液浓度调节，使得在终端时刻t时，荧光粒子浓度分布达到期望的nd(x,z)分布。
2、采用变分方法，得到公式(8)的伴随模型，其中变分量为n(x,z,t)和c0。
式中μ(x,z,t)为拉格朗日函数，满足边界条件和末端条件。
终端时刻条件为：
μ(x,z,t)＝nd(x,z)-n(x,z,t)(24)
形式梯度形式为：
公式(25)中c(x,z,t),μ(x,z,t)简写为c,μ，且
3、对于二维的荧光粒子移动模型，荧光粒子的对流项系数要远大于荧光粒子的扩散系数，荧光粒子移动耦合模型极易不稳定。因此引入fenicsproject，并结合有限元方法，得到荧光粒子移动模型。
此时，荧光粒子移动模型的fenics形式：
式中m＝1,2,3,...,m，m为[0,t]之间的离散数，δt＝t/m为时间间隔量，w(x,z)为fenicsproject设定的测试函数，满足第一(dirichlet)边界条件。此时伴随模型的fenics形式如下：
4、模型计算过程如下，伴随模型类似：
1)给定荧光粒子移动模型的fenics形式的时间和空间离散数，给定空间的网格形式，给定荧光粒子移动模型的fenics形式的边界条件。
2)定义测试函数w(x,z)的形式(比如：一次lagrange插值多项式形式)；
3)设定初始时间t＝0，求解荧光粒子移动模型的fenics形式；
4)存储t(x,z,tm)的数值，判断tm是否小于等于t，返回第3步，否则停止。
5、因此，基于fenicsproject的形式梯度优化方法步骤如下：
1)选择一个初始猜测值c0；
2)求解荧光粒子移动模型及其伴随模型；
3)采用fenicsproject自带的积分工具获得形式梯度；
4)执行最优测试.如果当前点是最优的,那就停止；否则执行第5步；
5)根据形式梯度信息得到下一步的搜索方向,确定最优步长,获取新的c0
具体为，为了防止形式梯度优化方法陷入局部最优解，使用近似梯度方法先得到一个粗略的最优解(次优解)，然后在通过形式梯度优化方法得到最优解。具体近似梯度步骤如下：
1)将需要优化的量c0在猜测区间中分成n个(我们取10等份)，那么就有10组优化猜测值；
2)加一个扰动项系数β，(取[0,0.05]之间的数，这里我们一般取0.01，根据经验而定)；
3)计算每个目标函数之差
4)在第k+1步，采用平移变换算子
5)计算
6)最后，将n个最优值中选取最小值，并记录此时的优化量c0。
如图8是δt＝0.1s，t取15s时，基于fenicsproject的形式梯度优化方法和近似梯度优化方法结合的控制效果图和溶液浓度控制曲面。其中横坐标是x方向，纵坐标是z方向，深度图为荧光粒子移动浓度图/nacl溶液浓度图。左边最后一个图的粒子分布是我们期望的分布，通过右边的浓度变化，左边第一个的粒子原始分布图最终达到了期望的左边最后一个图的期望分布。
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员，在不脱离本发明构思的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明保护范围内。