图1 静电电机结构示意图
Fig.1 Structure of electrostatic motor
摘要 静电电机不使用磁铁、无需线圈、结构简单、效率高，在微机电系统中有着潜在的应用前景，但目前基于电容可变原理的直线型静电电机还缺乏理论模型。该文以一种双激励直线型静电电机为例，通过求解多层介质情况下的泊松方程和拉普拉斯方程，获得电机内部电场分布情况；继而通过静电能量的虚位移法获得电机推力特性，由此建立起完整的电机数学模型，并计算得到电机动子的推力曲线和最大电场强度数值，该推力曲线可用于指导电机换相系统设计；随后建立电机的二维有限元模型，仿真结果验证了所提数学模型的正确性；最后根据电机各个参数与最大推力的关系并结合Sobol灵敏度分析结论对电机进行优化设计。结果表明，输入电压，介电常数和定、动子导体间距离对电机推力的影响最大，电极宽度对推力的影响最小。据此可指导该类电机的设计。
关键词：直线型静电电机 数学模型 多介质电场分析 灵敏度分析
静电电机是一种利用介电弛豫原理或者电容可变原理运行的电机[1]。由于它不使用磁铁、无需线圈、结构简单、效率高、容易实现微型化，特别是微型静电电机的平面结构和制造工艺与微机电系统（Micro-Electro-Mechanical System, MEMS）器件高度兼容[2]，因此研究超微型旋转型静电电机较多，如顶驱动（top-drive）、侧驱动（side-drive）、晃动驱动（wobble-drive）、中心钉驱动（center-pin-drive）、轴瓦、法兰盘（flange）和快门静电型（shuffle）等[3]。
相较于超微型旋转型静电电机，宏观尺寸的直线型静电电机研究较少[4-8]。这方面有代表性的是1995年由东京大学提出的一种名为双激励多相静电驱动（Dual Excitation Multiphase Electrostatic Drive, DEMED）的新型静电直线电机[9]。该电机定子和动子都装有电极，在定子和动子上施加电压可以让动子做直线运动。动子电极为三相，定子电极为三相，采用矩形波电压驱动；通过调节电压大小和频率可以调节动子的运动速度。该电机只有cm级大小，但是推力可达到数十N且其移动速度超过1m/s。鉴于上述特点，该类型直线电机不仅可应用在微型机械领域，常规的设备如位移器和机器人伺服驱动等同样可以应用静电电机[10-11]。在此基础上，日本学者A. Yamamoto等进一步对电极的形状进行改进[12]，为减小施加电压过高而导致的推力波纹，他们提出将电极设计成倾斜状或者V形，这种设计可以使动子平滑移动，有效地降低了电机运行时的噪声及振动，提升了电机的运行效率。M. Bahoura等将静电薄膜电机应用到仿生鱼的驱动中，利用静电薄膜电机和一个小型灵活的电力传输系统制造出一种能够通过自身振动实现前进运动的仿生鱼[13-14]。
为了对电机性能进行分析，A. Yamamoto利用电容网络矩阵对已制成电机的推力特性进行了分析[13]。由于DEMED一共有六相，所以电机可以看作是一个6端口网络，6个端口之间互相以电容相连，这些电容可以用电容矩阵来表示。动子板在运动时，电容矩阵内各个容值都随动子位置的改变而改变。电容的变化反映静电能的变化，进而反映推力的变化。因此基于电容可变原理的静电电机推力特性可以用电容网络来进行描述，通过测量电机动子板处于不同位置处的电容矩阵可以推算该位置处的推力大小。由于电容的大小不受外部驱动电路的影响，通过这种方法可以对未连接外部驱动的电机的推力特性进行分析[13]。
但上述方法是一种后验式测量方法。目前针对该类型直线静电电机缺乏理论模型以预估其机械性能。虽然用有限元方法可以较准确地获得电机性能参数，但其存在计算时间长的固有缺点，不易得到电机性能与电机参数的关系，不适用于电机优化设计，因而获得一套用以预估和优化电机性能的理论模型是十分必要的。本文以一种双激励直线型静电电机为例，通过求解静电方程组，获得电机内部电场分布情况，从而得到电机的推力曲线和最大电场强度，并用有限元模型进行了验证。最后分析了电机各个参数对输出力的灵敏度，以期指导该类电机的设计。
本文以定子三相6根导体、动子单相4根导体构成一个空间周期的静电电机为例，对静电电机进行理论研究和性能仿真分析。
静电电机主要由定子板、动子板和绝缘液构成，其结构示意图如图1所示。电极采用电路板工艺蚀刻在极板表面。由于静电电机的驱动电压较高，因此需要在两板之间注入绝缘液来避免高压造成的击穿，同时减小运动时两块板子之间的摩擦。
图1 静电电机结构示意图
Fig.1 Structure of electrostatic motor
静电电机电压变换示意图如图2所示。将动子电极极距记作
。图2所示黑色矩形代表带负极性电压的电极，白色矩形代表带正极性电压的电极。在运动过程中动子板电极通电电压Us1与Us2始终保持不变，定子板电极的电压Ua、Ub和Uc随动子板位置变化而变化，具体变化如图2a所示。时态1处为电机一个运动周期的初始位置，动子板受到向右的静电力移动
到时态2处，此处动子板的横向合力为零，将该位置称为平衡点位置。令定子板电压按时态2上标注的方式换相，动子板获得继续向右的静电力，再次运动到时态3处的平衡位置，此时定子再次换相。
图2 静电电机电压变换示意图
Fig.2 The schematic diagram of voltage transformation of electrostatic motor
以此类推，通过在平衡位置处为定子板电压换相，动子板可以获得持续向右的力。图2b为定子板三相电压的换相时序。从图中看出，定子板在一个周期内共换相6次，将一个周期划分成6个时态，每个时态内动子板都经历一次水平推力从固定值减小到零的过程。
首先做如下两点假设：
（1）由于电极导体厚度（mm级）与电机厚度（mm级）相比非常小，电机导体视为平面导体。
（2）电极导体电荷均匀分布，即面电荷密度为常数。忽略边缘效应。
为简化分析，假设电机定子介质、动子介质和绝缘液的介电常数相同，此时电机模型如图3所示。采用图示坐标系后，每根导体Ci的中心点坐标为（xi, 0, zi），带电量为Qi，电势为ji（以无穷远处为零电势），施加电压为Ui，i=1,…,10。
图3 电机模型
Fig.3 The model diagram of motor
计算静电场中导体受力一般有两种方法：①通过计算导体处电场强度得到导体受力；②通过静电系统能量的虚位移法计算受力。由于电场强度计算涉及到矢量运算比较复杂，本文采用第②种方法。图3所示静电系统的静电能为
（1）当动子在x方向有一微小位移Dx时，各个电极的带电量Qi和电势ji都会产生变化，由此系统静电能会有相应变化DW，则动子的受力为
（2）
因此需计算在固定电压下，动子处于不同位置时导体的带电量Qi和其电势ji。
第Ci根导体的电势为
（3）式中，
为第j根导体在第i根导体处产生的电势，称为互电势。
在多层介质存在情况下，导体的电势不仅与其自身带电量有关，还与介质界面上的极化电荷有关。因此本文首先计算单一介质情况下的电势，然后计算三层介质下的电势。
由于导体是等势体，其自电势可认为是导体上电荷在导体中心点处产生的电势，单一介质下自电势如图4所示，表示为
图4 单一介质下自电势
Fig.4 Self potential in a single medium
（4）式中，jii为第i根导体在其自身位置产生的电势，称为自电势；e 为介质的介电常数；a为导体宽度；L为导体长度；电势的上标0表明该电势为原导体产生的电势，以区别于后面的镜像导体。
单一介质中的互电势可认为是导体j上电荷在导体i中心点处产生的电势，单一介质互电势如图5所示，表示为
（5）
图5 单一介质互电势
Fig.5 Mutual potential in a single medium
为计算本文两种材料三层介质下静电体系的电势分布，首先计算一般性的三种材料三层介质中点电荷的电势，各层介质的介电常数分别为e1、e2和e3，介质界面相距l，电荷距第Ⅰ层介质距离为d，三层介质点电荷系统如图6所示。采用柱坐标系，点电荷q位于原点，用静电场的电势傅里叶积分变换求解边界条件下的Ⅱ区泊松方程和Ⅰ、Ⅲ区拉普拉斯方程，各区的电势[15]可写为
（6）
图6 三层介质点电荷系统
Fig.6 Point charge system in three-layer dielectric
（7）
（8）
式中，A(k)、B1(k)、B2(k)和C(k)为待定系数；J0(kr)为零阶贝塞尔函数。
边界条件为
本文关心的是点电荷所在II区域的电场分布j2。由边界条件和方程组可以求出
（9）
（10）其中
利用贝塞尔函数的积分
可得
由此可将式（7）改写为
由此可知在三层介质情况下，Ⅱ区中的电势j2是原来点电荷q产生的电场及z向无穷多镜像电荷产生电场的叠加（全空间为介电常数为e2的介质）。各电荷的位置及大小见表1。
表1 电荷位置及大小
Tab.1 The position and amount of charge
电荷位置（z轴）电荷大小 0
三层介质导体自电势如图7所示。对应于本文的情况，由于导体都位于中间介质中，Ci导体产生的自电势jii（见图7a）即该导体在Ci处产生的自电势
（见图7b）以及所有镜像导体在Ci处产生的自电势
的叠加（见图7c）。
（11）式中，的计算方法同互电势（见式（5））。
图7 三层介质导体自电势
Fig.7 Self potential in 3 layers medium
由于和都与Qi有关，因此可将自电势写为
（12）式中，ki,i为自电势系数。
同理，互电势可写为
（13）式中，ki,j为互电势系数。
由式（4）和式（5）可知，自电势系数和互电势系数仅与材料属性、电机几何尺寸以及导体相对位置有关，与电学参数无关。不同的动子位置对应不同的电势系数。
至此，静电电机各导体的电势可写为
（14）其中
由于
（15）
因此，将式（14）通过式（15）的运算后能够得到9个方程。同时
（16）将式（15）和式（16）联立，最终得到关于Qi的10个方程，求解得到Qi，代入式（14）得到ji，然后由式（1）获得动子在特定位置处的系统静电能量。将动子在各个位置处的静电能量由式（2）进行差分，最终获得静电电机动子在各位置处的受力。
一般而言，电机施加的电压越高，电机输出力越大。但可施加的最大电压受限于电机材料的击穿场强，因此需校核特定电压下电机的最大电场强度。
Ci导体在空间某点处产生的电势可参考式（13）得出，记为jn，因此全部导体在该处产生的电势为
由于电场强度E为
（17）
其中
因此电场强度大小为
（18）
据此可以得出确定的电压下动子在不同位置处静电电机全空间的电场强度大小，由此校核电机 参数。
电机的计算参数见表2。
表2 电机计算参数
Tab.2 The parameters of motor for caculation
参 数数 值 电压/V1 000 电极长度/mm60 电极宽度/mm0.135 定子、动子电极间距/mm0.1 定子极距/mm0.32 动子极距/mm0.48 板子厚度/mm0.5 介质相对介电常数2
设置Ua、Ub和Us1的电压为1kV，Uc和Us2的电压为0V。由于定子3根导体构成一个空间周期，因此计算过程中动子运动范围为3根定子导体极距（计算步数为150步），计算在该电压下动子运动过程中的受力过程和电机中的最大电场强度。电机长度为一个定子空间周期，约2.5mm。电机动子受力曲线如图8所示。
图8 电机动子受力曲线
Fig.8 The thrust curve of motor slider
从图8看出，电机动子板在第一根定子导体处受力为零。当它移动到接近第二根定子导体处时，获得最大推力，随后减小。移动到定子第二和第三根导体中间时受力为零，随后受力方向反向，并在经过第三根导体后达到反向最大力。电机最大推力约为38mN。同时可得在该运动周期内的最大电场强度为
V/m。该推力结果与文献[4]的实验结果一致。
利用同样方法获得图2a所示的6个时态下动子板的推力曲线合并如图9所示。每个时态下推力曲线的波形都类似，且后一个时态的波形比前一个时态右移半个定子极距。同时，利用该图可以对定子板矩形波电压的换相点进行设计，如在推力最大、推力为零或平均推力最大处换相，以实现电机在一个方向上的持续运动。
图9 静电电机6个时态推力曲线合并
Fig.9 Thrust curves combination under 6 steps of electrostatic motor
为了验证上述理论模型的准确性，利用Maxwell 2D模型对该静电电机进行有限元仿真，建立静电电机的有限元模型如图10所示。
图10 静电电机的有限元模型
Fig.10 FEM of electrostatic motor
选择材料时，将定子板、动子板和绝缘液都选择为相对介电常数为2的绝缘材料Polyimide Quartz，外层空气介质选择Air，电极材料选择为Copper。边界条件选为气球边界，半径为1cm。图11为动子在一个周期内的电机受力仿真曲线。
图11 电机受力仿真曲线
Fig.11 The simulation curve of thrust of the motor
对比图8和图11，两者的大小和变化趋势基本一致，从而验证了所提数学模型的正确性。
在用Matlab搭建的数学模型计算电场强度时，一个周期内的最大电场强度为
V/m；同等电机参数下采用有限元模型计算得到的最大电场强度为
V/m。两者的差别是由于理论模型推导时做了诸多简化造成的。在实际电机设计时，可先由理论模型给出一组合适的电机参数，再由有限元进行校验。
由于电机参数很多，因此有必要了解各参数分别对电机最大推力的影响程度，以期指导电机设计。
在前述静电电机数学模型的基础上，利用Sobol灵敏度分析法对各参数的灵敏度进行分析。数学模型中有8个主要变量：输入电压、相对介电常数、电极长度、电极宽度、板间距、定子板厚度、动子板厚度和电极间距，输出量设定为动子板一个周期内推力的最大值。此时的样本数量为8，设置采样的数量为4。用Sobol Sequence对样本进行取样，生成4行16列矩阵。将矩阵的前8列设置为矩阵A，后8列设置为矩阵B，用矩阵B中的第i列替换矩阵A中的第i列，得到矩阵ABi，由此得到10组矩阵。灵敏度分析各参数的取值范围见表3。
表3 灵敏度分析各参数的取值范围
Tab.3 The value range of parameters for sensitivity analysis
参 数取值范围 输入电压U/V50～5 000 相对介电常数1～80 电极长度L/mm20～100 两板间距dSR/mm0.05～1 定子板厚度dSt/mm0.08～0.5 动子板厚度dRt/mm0.08～0.5 电极宽度a/mm0.05～0.5 电极间距b/mm0.02～0.2
利用Sobol法灵敏度分析计算公式，得到各个自变量对动子最大推力灵敏度的影响指数，静电电机各参数对推力的灵敏度如图12所示。
图12 静电电机各参数对推力的灵敏度
Fig.12 The sensitivity of each parameter of electrostatic motor to thrust
从图12中可以看出，输入电压U和相对介电常数er对电机推力的影响程度最大，其次为两板间距dSR，电极宽度b对推力的影响最小。据此可优化电机本体来增大电机的输出力。
由电机的理论模型可以获得电机最大推力和最大电场强度与电机参数之间的关系，如图13所示。
结合各参数对受力输出的影响结果，对电机参数进行了优化，优化后电机参数见表4。
对优化后的模型进行理论计算，得到优化后的电机推力曲线如图14所示。与优化前数据（见图8）相比，最大推力提高至约4倍。由有限元计算获得的电场强度分布如图15所示，可以看出，优化后电机的受力显著增大，最大电场强度为
V/cm，可据此选择电机材料。
图13 最大推力和最大电场强度与某些电机参数的关系
Fig.13 The relationship of maximum thrust and field strength to some parameters of the motor
表4 优化后电机参数
Tab.4 The optimised motor parameters
参 数数 值 输入电压/V1 000 相对介电常数6.8 电极长度/mm60 两板间距/mm0.05 定子板厚度/mm0.12 动子板厚度/mm0.24 电极宽度/mm0.135 电极间距/mm0.185
图14 优化后的电机推力曲线
Fig.14 The thrust curve with optimised motor parameters
图15 优化后的电机电场强度
Fig.15 The electric field intensity with optimised motor parameters
本文针对目前基于电容可变原理的直线型静电电机缺乏理论模型的现状，以一种双激励直线型静电电机为例，通过求解多层介质情况下的泊松方程和拉普拉斯方程，获得电机内部电场分布情况。继而用静电能量的虚位移法获得电机推力特性，由此建立起完整的电机数学模型并计算得到电机动子运动一个机械周期的推力曲线和最大电场强度曲线。然后由二维有限元模型验证了所提数学模型的正确性。最后结合数学模型，利用Sobol灵敏度分析法获得电机各个参数对电机最大推力的灵敏度，以指导该类电机的设计。
参考文献
[1] 王冰, 尤彩红, 王文襄. 微型及超微型压力传感器的开发[J]. 仪表技术与传感器, 2009(增刊1): 180-181.
Wang Bing, You Caihong, Wang Wenxiang. Deve- lopment of micro and ultra-micro pressure sensors[J]. Instrument Technology and Sensors, 2009(S1): 180-181.
[2] 齐臣杰, 刘理天, 谭智敏. 微电机和微动力MEMS[J]. 电子科技导报, 1999(12): 21-23.
Qi Chenjie, Liu Litian, Tan Zhimin. Micro motor and micro power MEMS[J]. Electronic Technology, 1999(12): 21-23.
[3] 唐苏亚. 非电磁原理微电机的应用前景展望[J]. 微电机, 2007, 35(2): 70-74.
Tang Suya. Application prospect of non-electromagnetic micromotors[J]. Micromotors, 2007, 35(2): 70-74.
[4] 姬相超, 赵希梅. 永磁直线同步电动机的自适应时滞控制[J]. 电工技术学报, 2020, 35(6): 1231-1238.
Ji Xiangchao, Zhao Ximei. Adaptive time delay control of permanent magnet linear synchronous motor[J]. Transactions of China Electrotechnical Society, 2020, 35(6): 1231-1238.
[5] 凌志健, 赵文祥, 吉敬华. 高推力永磁直线作动器及其关键技术综述[J]. 电工技术学报, 2020, 35(5): 1022-1035.
Ling Zhijian, Zhao Wenxiang, Ji Jinghua. Overview of high force density permanent magnet linear actuator and its key technology[J]. Transactions of China Electrotechnical Society, 2020, 35(5): 1022- 1035.
[6] 张博, 齐蓉, 林辉. 激光切割永磁直线伺服系统的反演滑模控制[J]. 电工技术学报, 2018, 33(3): 642-651.
Zhang Bo, Qi Rong, Lin Hui. Back-stepping sliding mode control of laser cutting permanent magnet linear servo control system[J]. Transactions of China Electrotechnical Society, 2018, 33(3): 642-651.
[7] 武志涛, 杨永辉. 一种永磁直线电机驱动X-Y平台精密轮廓跟踪控制策略[J]. 电工技术学报, 2018, 33(17): 4037-4043.
Wu Zhitao, Yang Yonghui. A precise contour tracking control method for X-Y table driven by permanent magnet linear motors[J]. Transactions of China Electrotechnical Society, 2018, 33(17): 4037- 4043.
[8] 孙宜标, 仲原, 刘春芳. 基于LMI的直线伺服滑模位移跟踪控制[J]. 电工技术学报, 2019, 34(1): 33-40.
Sun Yibiao, Zhong Yuan, Liu Chunfang. LMI-based sliding mode displacement tracking control for per- manent magnet linear synchronous motor[J]. Transa- ctions of China Electrotechnical Society, 2019, 34(1): 33-40.
[9] Niino T, Higuchi T, Egawa S. Dual excitation multiphase electrostatic drive[C]//Proceedings of the Industry Applications Society Annual Meeting, Orlando, FL, USA, 1995: 1318-1325.
[10] 赵正平. 典型MEMS和可穿戴传感技术的新发展[J]. 微纳电子技术, 2015, 52(1): 1-13.
Zhao Zhengping. New developments in typical MEMS and wearable sensor technology[J]. Micro- nano Electronics, 2015, 52(1): 1-13.
[11] 王欣利, 程树康. 静电电动机及其研究发展状况[J]. 微特电机, 2001, 29(6): 12-15.
Wang Xinli, Cheng Shukang. Electrostatic motor and its development[J]. Micro Motor, 2001, 29(6): 12-15.
[12] Yamamoto A, Niino T, Higuchi T. Modeling and identification of an electrostatic motor[J]. Precision Engineering, 2006, 30(1): 104-113.
[13] Bahoura, Mohammed. FPGA implementation of high- speed neural network for power amplifier behavioral modeling[J]. Analog Integrated Circuits & Signal Processing, 2014, 79(3): 507-527.
[14] Kimura F, Yamamoto A, Higuchi T. FPGA imple- mentation of a signal synthesizer for driving a high- power electrostatic motor[C]//IEEE International Symposium on Industrial Electronics, 2011, DOI: 10.1109/ISIE.2011.5984345.
[15] 谢宝昌. 电磁能量[M]. 北京: 机械工业出版社, 2016.
Theoretical Research and Performance Simulation of Electrostatic Motor
Abstract Electrostatic motor has the merits of no magnets, no coils, simple structure and high efficiency, which has potential application prospects in microelectromechanical systems (MEMS). However, there is no theoretical model for the linear electrostatic motor based on variable capacitance. Taking a double-excited linear electrostatic motor as an example, the distribution of electric field inside the motor was obtained by solving the Poisson equation and Laplace equation in the case of multi-layer media. Then, the thrust characteristics of the motor were obtained through the virtual displacement method, thereby establishing a complete mathematical model of the motor and calculating the thrust curve and the maximum field strength of the electric motor. The thrust curve can be used to guide the design of the motor commutation system. Then, a two-dimensional finite element model of the motor was established, and the simulation results have verified the correctness of the mathematical model. Finally, according to the relationship between each parameter of the motor and the maximum thrust, combined with the Sobol sensitivity analysis, the motor was optimized. The results show that the input voltage, the permittivity and the distance between the stator and the slider have the greatest influence on the motor thrust, while the electrode width has the least influence. Accordingly, the design of such motors can be guided.
keywords：Linear electrostatic motor, mathematical model, analysis of multi-medium electric field, sensitivity analysis
中图分类号：TM359.4
DOI: 10.19595/j.cnki.1000-6753.tces.200598
国家重点基础研究发展计划资助项目（2015CB057503）。
收稿日期 2020-06-03
改稿日期 2020-09-30
李华峰 男，1974年生，教授，博士生导师，研究方向为新型电机的运行理论及控制。E-mail: lihuaf@nuaa.edu.cn（通信作者）
王伊凡 女，1994年生，硕士，研究方向为电机本体设计。E-mail: wangyifanwy@163.com
（编辑 陈 诚）