环形微机电陀螺作为一种MEMS(micro-electromechanical system)器件，利用Coriolis加速度效应实现角速度的测量。由于采用微机械加工工艺制作，它具有体积微小、易于与电路集成及成本低等特点。与采用质量块和梁的结构形式的微机电陀螺相比，环形特别是由多层环构成的微机电陀螺结构具有高度的对称性^[1-3]。

环形微机电陀螺的性能特点与它的对称性紧密联系。理想的环形结构可以实现陀螺驱动与检测模态的谐振质量、刚度系数和阻尼系数等参数的一致，易于实现驱动与检测振动频率调谐。通过好的结构设计可以降低结构的热弹阻尼，完美加工的环结构在采用真空封装的条件下可以达到很高的振动品质因数。在驱动与检测振动频率调谐且具有高振动品质因数的条件下，环形微机电陀螺可以具有很高的角速度灵敏度^[4-6]。

环形微机电陀螺的支承点可以设计在其中心处，其工作振动模态具有对称性，因而可以有效地抑制环境振动产生的虚假敏感信号^[1]；并且材料的热膨胀也很难导致振动结构产生内部应力，所以结构振动特性稳定，对温度变化不敏感。环形微机电陀螺既可以具有高灵敏度，又对环境条件变化不敏感，因而具有达到高稳定性的潜力。

但是实际的MEMS加工工艺并不能保证多环结构的完全对称性。多环结构误差的来源主要有密度、Young's模量、结构厚度和梁宽度等结构参数的不均匀性^[7-10]。这些缺陷的存在使得环形完美的对称性被破环，工作模态的频率不再一致，而且本征轴(固有振动模态下振幅达到极大值的轴向)也不再是环上的任意方向，而是与缺陷的等效位置相关，这会严重降低陀螺的稳定性。

结构参数不均匀性中对谐振子影响最严重的是Fourier展开式的4次谐波^[11-12]。为了确定本征轴位置，李巍等^[13]研究了通过推导振动位移的表达式并测试最佳检测时间来确定本征轴的方法。本文主要研究质量分布的4次谐波对环结构振动模态的影响，并提出一种基于谐频响应曲线来确定有加工缺陷的环形微机电陀螺本征轴方向的方法。该方法可以确定结构缺陷等效于质量分布4次谐波分量的位置，可为后续的陀螺的机械修调提供依据。

1 环形微机电陀螺的敏感结构 1.1 敏感结构与工作原理

本文所分析的环形微机电陀螺敏感结构如图 1所示，采用硅为结构材料，玻璃为基片，以静电键合和体硅深刻蚀技术加工。环形结构的中心支承锚点与玻璃基底键合连接。在环的外侧沿圆周均匀设置16个驱动和检测电极，这些电极与环之间形成电容。在驱动电极上施加静电电压，驱动陀螺以工作频率振动；当在垂直于振动平面内输入角速度时，因为Coriolis力的作用，检测模态被激励；由检测电极来检测电容的变化，经电路处理得到角速度信息。

1.2 环结构的工作振动模态

环形微机电陀螺的工作振动模态是在平面内2个本征轴成45°夹角的2节点弯曲振动模态，如图 2所示。这2个振动模态具有相等的自然频率^{[11-12, 14-15]}。

在没有结构缺陷的理想情况下，本征轴可以处在任意方向。以环结构设计时确定的电极所在方向为参考方向，令本征轴与参考方向重合，其上各点的位移变化可以表示为^[11]

$ \left\{\begin{array}{l} {u(t)=2 q_{1}(t) \cos 2 \varphi+2 q_{2}(t) \sin 2 \varphi} , \\ {v(t)=-q_{1}(t) \sin 2 \varphi+q_{2}(t) \cos 2 \varphi}. \end{array}\right. $

其中：φ为环上质点的空间方位角，u为其径向位移，v为其切向位移，q₁(t)和q₂(t)分别为2节点弯曲模态中第1和第2模态的变形量。2个模态的空间方程是sin2φ和cos2φ的组合。

2 环结构质量分布对工作模态的影响

实际加工的环结构存在尺寸的系统性和随机性误差，材料也可能存在密度分布的不均匀，这些因素都会破坏环的对称性。结构缺陷会影响2个工作振动模态的本征轴方向，2个自然角频率也不再相等。以下称其中自然角频率较低的本征轴为“重本征轴”，自然角频率较高的轴为“轻本征轴”^[11]。

结构尺寸误差对环的振动模态的影响效果可以等效为沿环的周向的密度不均匀分布，也即存在密度缺陷。设重本征轴的参考角度为φ₀，把环形谐振子的密度分布以Fourier级数展开，有

$ \rho=\rho_{0}+\sum\limits_{n} \Delta \rho \cos n\left(\varphi-\varphi_{0}\right). $

其中：Δρ为密度缺陷相对值，n为谐波的次数。显然重本征轴的位置直接对应于缺陷所在的角度。

使用多物理场仿真软件COMSOL对单层环结构进行动力学数值分析，检查不同的谐波分量对结构频率分裂的影响。所分析的环的直径为5 940 μm，高度为80 μm，壁厚为60 μm，分析的谐波次数为1~8。为了降低数值计算误差的影响，Δρ的取值范围较大，为0~200 kg/m³。数值分析的结果如图 3所示。

所分析结构的工作模态自然频率为12 928 Hz。从数值分析结果可见，1、3、5、6、7、8次谐波造成的频率分裂均小于1 Hz，2次谐波在分析范围内造成的频率分裂达到几Hz，而4次谐波可造成几百Hz的频率分裂。因此，对频率分裂影响最大的是密度缺陷的4次谐波分量，在同等条件下它们比2次谐波高出2个数量级，对频率分裂起着决定性的作用；其他次谐波分量对频率分裂的影响则可忽略。

3 本征轴方向测算原理

在频率分裂比较大的情况下，需要对结构进行机械修调，使得2个工作模态的自然频率接近于相等。机械修调的前提条件是需要找到本征轴的位置。但是由于无法直接确定结构的加工误差三维分布情况，实际结构本征轴的方向是未知的，需要通过间接的实验测试来推算本征轴的方向。

测算的基本方法是对环施加激励，在各电极处检测环振动的频率响应特性，再根据频率响应数据计算本征轴的方向。在驱动电极上施加激励，静电驱动力按谐波分解作用到整个环的圆周上，环的2个工作振动模态被激励起来^[13-16]；通过检测电极敏感振动信号，得到环的频率响应特性。

3.1 理论分析

环形振动陀螺的驱动模态和检测模态均为二阶系统，不考虑耦合量作用，两模态的传递函数为^{[2, 15]}

$ \left\{\begin{array}{l} {G_{1}(s)=\frac{1}{m_{1}} \cdot \frac{1}{s^{2}+2 \zeta_{1} \omega_{1} s+\omega_{1}{^{2}}}}, \\ {G_{2}(s)=\frac{1}{m_{2}} \cdot \frac{1}{s^{2}+2 \zeta_{2} \omega_{2} s+\omega_{2}{^{2}}}}. \end{array}\right. $

(1)

其中：m₁和m₂分别为2个模态的谐振质量，ω₁和ω₂分别为2个本征轴的自然角频率，ζ₁和ζ₂分别为阻尼比。

假设轻本征轴的角度为φ₀(也即4次谐波密度缺陷位置)，重本征轴的角度为φ₀±45°，且φ₀处于±45°之间。在0°电极上施加简谐力，同时引起两个模态的振动，在检测电极处的振动则是这2个模态振动叠加的结果。

设施加的简谐力幅值为1，分解为作用在2个振动模态上的分布力N₁和N₂：

$ \left\{\begin{array}{l} {N_{1}=\cos 2 \varphi_{0}}, \\ {N_{2}=\sin 2 \varphi_{0}}. \end{array}\right. $

(2)

在检测电极处，检测幅值是2个固有振动分量叠加的结果。例如在0°和22.5°处，频率特性分别为:

$ y_{0}(\mathrm{j} \omega)=N_{1} \cos 2 \varphi_{0} G_{1}(\mathrm{j} \omega)+N_{2} \sin 2 \varphi_{0} G_{2}(\mathrm{j} \omega), $

(3)

$ \begin{array}{c} y_{22.5}({\rm j} \omega) =N_{1} \cos 2\left(22.5^{\circ}+\varphi_{0}\right) G_{1}({\rm j} \omega)+\\ N_{2} \sin 2\left(22.5^{\circ}+\varphi_{0}\right) G_{2}(\mathrm{j} \omega). \end{array} $

(4)

对0°处的检测信号，把式(1)和(2)代入式(3)，有

$ \begin{aligned} y_{0}(\mathrm{j} \omega) &=\cos ^{2} 2 \varphi_{0} \cdot \frac{1 / m_{1}}{s^{2}+2 \zeta_{1} \omega_{1} s+\omega_{1}^{2}}+\\ & \sin ^{2} 2 \varphi_{0} \cdot \frac{1 / m_{1}}{s^{2}+2 \zeta_{1} \omega_{1} s+\omega_{1}^{2}}. \end{aligned} $

(5)

其实部、虚部分别为：

$ \begin{array}{c} \operatorname{Re}\left[y_{0}({\rm j} \omega)\right]= \frac{1}{m_{1}} \cdot \frac{\left(-\omega^{2}+\omega_{1}^{2}\right) \cos ^{2} 2 \varphi_{0}}{\left(-\omega^{2}+\omega_{1}^{2}\right)^{2}+\left(2 \zeta_{1} \omega_{1} \omega\right)^{2}}+\\ \frac{1}{m_{2}} \cdot \frac{\left(-\omega^{2}+\omega_{2}^{2}\right) \sin ^{2} 2 \varphi_{0}}{\left(-\omega^{2}+\omega_{2}^{2}\right)^{2}+\left(2 \zeta_{2} \omega_{2} \omega\right)^{2}}, \\ \operatorname{Im}\left[y_{0}({\rm j} \omega)\right]=\frac{1}{m_{1}} \cdot \frac{-2 \zeta_{1} \omega_{1} \cos ^{2} 2 \varphi_{0}}{\left(-\omega^{2}+\omega_{1}^{2}\right)^{2}+\left(2 \zeta_{1} \omega_{1} \omega\right)^{2}}+\\ \frac{1}{m_{2}} \cdot \frac{-2 \zeta_{2} \omega_{2} \sin ^{2} 2 \varphi_{0}}{\left(-\omega^{2}+\omega_{2}^{2}\right)^{2}+\left(2 \zeta_{2} \omega_{2} \omega\right)^{2}}. \end{array} $

(6)

在相位为90°奇数倍的角频率ω^*处，实部应为0，也即式(6)中Re[y₀(jω^*)]=0，可得

$ \begin{array}{c} {\frac{m_{2}}{m_{1}} \cdot \frac{\left(-\omega^{* 2}+\omega_{1}^{2}\right) \cos ^{2} 2 \varphi_{0}}{\left(-\omega^{* 2}+\omega_{1}^{2}\right)^{2}+\left(2 \zeta_{1} \omega_{1} \omega^{*}\right)^{2}}+} \\ {\frac{\left(-\omega^{* 2}+\omega_{2}^{2}\right) \sin ^{2} 2 \varphi_{0}}{\left(-\omega^{* 2}+\omega_{2}^{2}\right)^{2}+\left(2 \zeta_{2} \omega_{2} \omega^{*}\right)^{2}}=0}. \end{array} $

(7)

由式(7)推导得

$ \tan ^{2} 2 \varphi_{0}=\frac{m_{2}}{m_{1}} \frac{\omega^{* 2}-\omega_{1}^{2}}{\omega_{2}^{2}-\omega^{* 2}} \cdot \frac{\left(\omega_{2}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{1} \omega_{2}\right)^{2}}{\left(\omega_{1}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{2} \omega_{1}\right)^{2}}. $

(8)

进一步有：

$ \begin{array}{c} \tan 2 \varphi_{0} =\sqrt{1-\frac{m_{1}-m_{2}}{m_{1}}}\cdot \\ \sqrt{\frac{\omega^{* 2}-\omega_{1}^{2}}{\omega_{2}^{2}-\omega^{* 2}} \cdot \frac{\left(\omega_{2}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{1} \omega_{2}\right)^{2}}{\left(\omega_{1}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{2} \omega_{1}\right)^{2}}}，\\ \sqrt{1-\frac{m_{2}}{m_{1}}} \approx \sqrt{1-\frac{\Delta m}{m}}. \end{array} $

Δm/m₁表示了谐振质量的不一致性，由于$f=\sqrt{\frac{k}{m}} $，假设刚度k不变，则有$\frac{\Delta f}{f}=-\frac{1}{2} \frac{\Delta m}{m} $，也即当有1%的频率分裂时，谐振质量不一致的量约为2%，此时计算本征轴的角度时，由质量不一致计算的本征轴角度与忽略质量不一致所计算的结果相差1%左右。由于环结构设计的对称性，实际计算时，可假设m₂/m₁≈1，从而简化计算过程为

$ \tan ^{2} 2 \varphi_{0}=\frac{\omega^{* 2}-\omega_{1}^{2}}{\omega_{2}^{2}-\omega^{* 2}} \cdot \frac{\left(\omega_{2}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{1} \omega_{2}\right)^{2}}{\left(\omega_{1}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{2} \omega_{1}\right)^{2}}. $

(9)

其中Q₁=1/(2ζ₁)，Q₂=1/(2ζ₂)，分别为2个振动模态的品质因数。

同理，在其他电极处测试也可得到本征轴位置的计算关系。例如在22.5°电极处检测，有

$ \begin{array}{c} \frac{\tan ^{2} 2 \varphi_{0}+\tan 2 \varphi_{0}}{1-\tan 2 \varphi_{0}}=\frac{\omega^{* 2}-\omega_{1}^{2}}{\omega_{2}^{2}-\omega^{* 2}}·\\ \frac{\left(\omega_{2}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{1} \omega_{2}\right)^{2}}{\left(\omega_{1}^{2}-\omega^{* 2}\right)^{2}+\left(Q_{2} \omega_{1}\right)^{2}}. \end{array} $

(10)

从0°电极测试频率响应特性，获得2个振动模态的自然角频率、Q值和ω^*，再用式(9)即可计算得到环结构本征轴角度的2个可能值，但是尚不能唯一确定；利用式(10)计算也可以得到2个可能取值；两组解中相同者即为真解。

本征轴方向不同，在各电极上会显示出不同的频率特性，如图 4和5所示，重本征轴方向在10.0°时，从0°和22.5°电极处测得的相频曲线穿越±90.0°的角频率在两自然角频率之间；重本征轴方向在30.0°时，从0°和67.5°电极处测得地相频曲线穿越±90.0°的角频率在两自然角频率之间。反过来说，可根据实测频率响应曲线的这个特征确定重本征轴方向是在0°~22.5°还是在22.5°~45.0°。

3.2 数值分析检验

用有限元分析方法对具有密度缺陷的环结构进行谐响应分析，将所获得的谐响应特性参数代入本征轴方向的理论结果，与实际密度缺陷位置对比，检验上述测算方法的准确性。

在有限元分析软件ANSYS中建立有四点缺陷的模型^{[4, 6, 16]}用来近似密度缺陷的4次谐波分量。设置品质因数为Q₁=Q₂=1 000，在0°和180°位置施加方向相反的简谐力，在0°位置检测振动的谐响应特性。

仿真中，设置密度缺陷的角度φ₀^′分别为8.00°、10.00°、11.25°、15.00°和20.00°共5种情况，与通过测算方法得到的φ₀值进行对比。仿真模型的2个工作模态的自然频率分别为f₁^′=6 802 Hz，f₂^′=6 896 Hz。扫频频率步长设置为1 Hz，分析得到5种密度缺陷位置下的谐频响应曲线如图 6所示。从图 6中获取ω₁、ω₂和ω^*对应的频率值f₁、f₂和f^*，并通过式(9)可以得到缺陷角度φ₀，结果如表 1所示。

由表 1数据可见，本文所提出的本征轴角度测算方法所确定的结果与模型设置值的相对误差在2%以内。此误差中包含了数值分析过程的误差，而有限元分析的误差一般约1%。另外理论上存在的一个误差源是基于一个假设：2个振动模态的等效质量相等，当然实际上是不相等的。总的来说，2%的相对误差在工程实际中是可以接受的。

4 结论

本文通过在不同角度的电极处测试环结构的频率响应特性，推算环结构的本征轴角度，即可确定等效为密度分布4次谐波分量的结构缺陷所在的位置。且2个工作模态的频率分裂程度与密度缺陷的大小呈线性关系。本文提出的本征轴向推算方法，可以用于测算环形微机电陀螺敏感结构缺陷的位置，为采用飞秒激光加工等方法对结构误差进行修调、提高敏感结构的对称性、提高环形微机电陀螺的性能提供支撑。