Inclination Effect on the Periodic Response of a Symmetrical MEMS Gyroscope

The inclination effect caused by fabrication errors on the periodic response of a symmetric MEMS gyroscope is investigated. The dynamic equation is established considering the inclination effect on support stiffness and electrostatic forces. The periodic response is obtained by the averaging method. The two-variable singularity theory is employed to study the bifurcation characteristics and give transition sets on the DC-AC voltage plane, which divide the plane into four persistent regions. The amplitude-frequency curves demonstrate that only the two persistent regions with low voltages are feasible for the gyroscope. Both over-etching and under-etching reduce the feasible region. The effect of parameters on the performance is present. The mechanical sensitivity and nonlinearity increase with the voltages. With the increase in the inclination angle, the mechanical sensitivity and nonlinearity decrease first and then increase. The full temperature stability of the mechanical sensitivity is also considered. The variation in mechanical sensitivity with temperature is small at a large voltage and negative inclination angle. Under-etching, which leads to small nonlinearity and good temperature stability, is more beneficial to the performance of the gyroscope than over-etching.


Introduction
MEMS gyroscopes are widely used for their small size, light weight, low manufacturing cost, and high reliability [1]. However, nonlinearities, including nonlinear stiffness and nonlinear electrostatic force, usually exist in MEMS gyroscopes. Nonlinearities can lead to multi-solutions and softening/hardening characteristics of the response, which affect the stability, sensitivity, and reliability of MEMS gyroscopes. Therefore, the nonlinear dynamics of MEMS gyroscopes have attracted significant attention from researchers.
Braghin et al. [2] obtained a high deformation range of the support beams of a tuning-fork MEMS gyroscope and investigated the nonlinear vibrations both experimentally and numerically. Tsai et al. [3] studied the stability and resonance of a micromachined gyroscope considering the nonlinearities resulting from flexure springs and electrostatic force. The external resonance and non-resonant hard excitation were considered. Martynenko et al. [4] investigated the nonlinear dynamics of a MEMS tuning fork gyroscope with elastic rods on the moving base using the asymptotic two-scales method. Lajimi et al. [5] examined the effect of the DC bias voltages, the AC drive voltage, drive frequency, and the quality factors on the nonlinear dynamics and the sensitivity of a beam-rigid body gyroscope based on a computational shooting method. Shang et al. [6] demonstrated the nonlinear vibration of a parametric exciting gyroscope considering the principal parameter resonance and 1:1 internal resonance. Hao et al. [7] discussed the influence of DC bias voltage and comb spacing on the nonlinearity of electrostatic force and the dynamic response of a 4-DOF electrostatically driven gyroscope. Liang et al. [8] analyzed the nonlinear mechanism of MEMS vibratory ring gyroscopes in which the gyroscopic coupling and geometrically and structurally nonlinear couplings are all considered. Zhang et al. [9] investigated the effect of driving and sensing voltages on the bifurcation of the periodic response of a micro gyroscope with comb driving and parallel-plate sensing.
Currently, most MEMS gyroscopes are fabricated by the deep reactive ion etching (DRIE) process. Inevitably, over-etching or under-etching will lead to the non-parallelism of the working planes of the driving and sensing combs, as well as the change in the cross-section of the support beams, which will affect the performance of the gyroscope.
To clarify the influence of fabricating error, Tay et al. [10] investigated the effect of non-parallel plates on the capacitance, sensitivity, electrostatic force, and electrostatic spring constant of a micro accelerometer using a theoretical approach and finite-element analyses. Dong et al. [11] considered the effect of non-parallel combs on the reliable operation range for step and pulse signals in different capacitive configurations and gave the analytical expressions of the critical step acceleration and the critical time duration of the pulse acceleration of a capacitive sensor whose comb electrodes had a small angle of the decline. Guo et al. [12] introduced a new theoretical model of plate capacitor for a MEMS comb actuator. The bulk fabrication process and edge effect on the capacitance are both considered using the integration method and conformal transformation theory. Zhong et al. [13] studied the inclination effect of the fingers on the capacitance, driving electrostatic force, and electrostatic spring constant and discussed the nonlinear dynamics of a MEMS resonator. Yan et al. [14] analyzed the electrostatic tuning mechanism of electrical stiffness produced by the electrostatic force and gave the dependence of resonance frequency variation on the inclination angle and DC bias voltage of an inclined disk resonator array. Zhang et al. [15] studied the influence of inclination angle on the resonance frequency, stiffness hardening effect, and sensitivity of a 4-DOF MEMS gyroscope.
However, existing studies have mainly focused on the influence of the inclination angle on the electrostatic force. For actual gyroscopes, the inclination angle can also vary the cross-section of the support beams, thereby changing the support stiffness and resonant frequency.
In this paper, we propose an investigation of the nonlinear dynamics of a symmetrical MEMS gyroscope considering the inclination effect on both the electrostatic force and support stiffness of the gyroscope. The structure of the paper is organized as follows. In Section 2, the dynamical model of the MEMS gyroscope is established considering the influence of the inclination angle. In Section 3, the periodic response is obtained based on the averaging method. In Section 4, we study the effect of parameters on the bifurcation characteristics, amplitude-frequency response, and the performance of the gyroscope. Finally, the conclusions are summarized. Figure 1a is a schematic diagram of the symmetrical gyroscope [16]. The proof mass is connected to the outer frame by internal elastic beams and then supported on the anchor by external elastic beams. Figure 1b shows the corresponding spring-mass model. The blue and green springs represent the external and internal elastic beams. x and y are the driving and sensing directions, and z is the direction of the vertical paper facing outward. m p is the proof mass, m f d and m f s are the masses of the driving and sensing outer frames, and m f d = m f s = m i . The proof mass can move in two directions, while the driving and sensing outer frame can move in only one direction. During the gyroscope operation, the DC voltage V D is applied to the movable structure, and the AC voltages V A1 and V A2 are applied to the fixed combs in the driving direction. The electrostatic forces drive the proof mass and the driving outer frame to vibrate in the driving direction. When the carrier of the gyroscope has an angular velocity around the z direction, the proof mass vibrates in the sensing direction due to the Coriolis effect. Then the vibration is transmitted to the sensing outer frame. The angular velocity can be determined by measuring the vibration in the sensing direction.  Figure 2 shows the relation between the displacement of the gyroscope coordinate system and the inertial coordinate system, where x-o-y is the gyroscope coordinate system and x a -o-y a is the inertial coordinate system. The displacement of the proof mass and the driving and sensing frames in the inertial coordinate system can be written as x aA = x A cos ω z t − y A sin ω z t, y aA = x A sin ω z t + y A cos ω z t, x aB = x A cos ω z t, y aB = x A sin ω z t, x aC = −y A sin ω z t, y aC = y A cos ω z t

The Dynamical Model of the MEMS Gyroscope
The velocities can be obtained by the derivative of Equation (1) to time t.
Thus, the kinetic energy of the system is We can also give the potential energy and Rayleigh dissipation function as where k o , k i , α o , and α i are the linear and nonlinear stiffness of the external and internal support beams and c is the damping coefficient. According to the Lagrange equation, the dynamic equation of the system can be obtained as follows: Considering that (m p + 2m i )ω z 2 4k o + 2k i ,ẏ ẋ, Equation (6) can be simplified as where where F ex and F ey denote the electrostatic forces in the driving and sensing directions, respectively. Their expressions will be given below. The relation between the displacements of the gyroscope coordinate system and the inertial coordinate system.

The Electrostatic Forces in the Driving and Sensing Directions with Inclination Effect
As shown in Figure 3, two kinds of capacitances exist between the combs, the longitudinal capacitance C L and the transverse capacitance C T . For the non-ideal state, the working faces of the combs are not parallel, which consequently changes the capacitances. Assuming the inclination angle is small, the longitudinal capacitance considering the edge effect and the inclination angle is [12] C L = εNhl hθ+d + εNh π 1 + ln 1 + πl hθ+d + ln 1 + πl hθ+d + εNl π 1 + ln 1 + πh hθ+d + ln 1 + πh hθ+d (9) where θ is the inclination angle, l is the overlap length of combs, h is the structural thickness of the gyroscope, d is the spacing between the fingers, ε is the dielectric constant, and N is the number of fingers on a single side of the combs. Since the separation between the combs is large, the inclination effect on the transverse capacitance can be ignored. The transverse capacitance is given as where w is the finger's width, and x 0 is the separation between the combs. According to Ref. [17], the electrostatic force can be obtained by where C and V are the capacitance and voltage between the electrodes, respectively. Thus, the electrostatic forces in the driving and sensing direction can be given as where C l and C r are the capacitances on the left and right as the movable structure has a displacement x. C t and C b are the capacitances on the top and bottom with the displacement y.
The non-ideal state (Over-etching)

The Inclination Effect on the Restoring Force
The inclination angle also changes the cross-section of the support beams and the stiffness of the gyroscope, as shown in Figure 4. However, it is difficult to represent the inclination effect on the stiffness analytically. Thus, a finite element analysis is carried out. We also consider the effect of the environmental temperature on Young's modulus of the material.
where T 0 = 300 K, E 0 = 170 GPa is Young's modulus of monocrystalline silicon at 300 K, and υ = 5 × 10 −5 . The FEA model is presented in Figure 5. The combs are ignored since they do not affect the calculation results of the restoring force. The locally refined hexahedron mesh generated by the sweep method is used. The minimum element edge length is 6.21 µm, and the total number of elements is 15,956. The structures of the driving and sensing directions are wholly symmetrical, so only the restoring force of the x direction is considered. The fixed constraint is applied to the surface A, B, C, and D. The environmental temperature is applied to the structure as a thermal load. By giving a displacement in the x direction at the center E of the proof mass, the restoring force of the fixed constraint can be obtained by the software ANSYS Workbench 14.5.

Anchor
Anchor Proof mass The non-ideal state (Over-etching) By changing the inclination angle from −2 • to −2 • in steps of 0.08 • , the environmental temperature from 240 K to 360 K with a step size of 10 K, and the displacement of the x direction from 0 µm to 30 µm with steps of 2 µm, we can calculate the restoring force with different parameters. Figure 6 shows the relation between the restoring force and the displacement corresponding to different inclination angles and the environmental temperature, indicating that the variation is nonlinear. For convenience, we fit the restoring force cubically as The fitting results are also shown in Figure 6. Figure 7 gives the influence of the inclination angle and the environmental temperature on K and Γ. We also fit them as polynomials about θ and T.
It was demonstrated that the influence of the inclination angle on the linear stiffness is cubic, and the effect on the nonlinear stiffness is linear. The effect of the environmental temperature is very small, so the room temperature (T = 300 K) was used for this study, except in Section 4.3.

The Static Pull-In Analysis in the Driving Direction
By setting X = x/x 0 , Y = y/x 0 , τ = ω x t, Equation (7) can be rewritten in the nondimensional form as where The non-dimensional electrostatic forces are where The amplitude in the driving direction is large, which can easily induce pull-in. Hence, the static pull-in analysis of the driving direction should be conducted. We consider a conservative system in the driving direction as follows.
The system's potential energy is obtained by setting the equilibrium position as the point with zero potential energy. Figure 8 shows the potential energy corresponding to different DC voltages and inclination angles. The parameter values are chosen in Table 1 in the following calculations if there are no special instructions. The system may have a stable periodic motion in the potential well near the equilibrium position X = 0. Pull-in is triggered when the proof mass displacement crosses the potential barrier near X = 1. The barrier height increases with the inclination angle but decreases with the DC voltage. The displacement of the peak point of the barrier can be regarded as the boundary amplitude of the pull-in.

Parameters Values
The We investigated the effect of the DC voltage and the inclination angle on the pullin boundary amplitude, as shown in Figure 9. It was found that the pull-in amplitude decreases with the DC voltage but increases with the inclination angle, which means pull-in is more likely to occur at a large DC voltage or a small inclination angle. Physically, the electrostatic force in the driving direction that causes the system to pull in increases with the DC voltage and decreases with the inclination angle. However, reducing the inclination angle leads to a rapid increase in the stiffness and the restoring force, which increases the pull-in boundary amplitude. Ignoring the inclination effect on the stiffness may bring about an opposite result. Therefore, it is necessary to consider the inclination effect on both the electrostatic forces and the stiffness.

The Periodic Solution and Its Stability of the Gyroscope
The periodic response is significant for MEMS gyroscopes. In this section, we investigate the periodic response of the gyroscope and its stability. Note that f (x) contains logarithmic terms that cannot be calculated analytically. f (x) is expanded into the Taylor series as follows. Figure 10 gives the relative error of the Taylor expansion of f (X). The maximum error is less than 3%. Substituting Equation (23) into Equation (19) leads to The solutions to Equation (17) can be assumed as where A 1 , A 2 , θ 1 , and θ 2 are the amplitudes and phase angles of X and Y. Considering the primary resonance yields Ω 2 = 1 + σ, where σ is the detuning parameter. According to the method of averaging [18], we have The integrals in Equation (26) can be solved as where the integral of the fractional term in F EXT , F EYT can be obtained by the residue theorem [9,19], and the expressions of F X11 , F X12 , and F Y12 are as shown below.
By setting A 1 = A 2 = θ 1 = θ 2 = 0, the nonlinear equations about the steady-state solutions to the periodic response can be obtained. The amplitudes and phase angles can be calculated using the homotopy continuation method. The stability of the steady-state solutions can be determined by the eigenvalues of the Jacobian matrix of Equation (27) [18].  Figure 11 shows the amplitude-frequency curves of the system. By using the fourthorder Runge-Kutta method, the solutions of Equations (18) and (20) are also given to verify the analytical solutions. The results are in good agreement. Note that the amplitudefrequency curves exhibit different dynamic characteristics. In Figure 11a, A 1 and A 2 have unique solutions corresponding to any value of Ω. In Figure 11b, a jump happens with the change in Ω, while the jump leads to pull-in in Figure 11c. Thus, it is necessary to analyze the effect of the parameters on the bifurcation characteristics.

The Singularity Analysis in the Driving Direction
In this section, the singularity theory is employed to study the influence of the parameters on the bifurcation characteristics in the driving direction. It is impossible to eliminate θ 1 to obtain the bifurcation equation about A 1 from Equation (27). Thus, the two-variable singularity theory is used in this section. By setting A 1 = u, θ 1 = v, σ = λ, the following bifurcation equations are considered.
The transition sets of Equation (29) are defined as D = B ∪ H ∪ DL, where B, H, and DL denote the bifurcation set, the hysteresis set, and the double limit point set, respectively. Their expressions are given as follows [20,21]. Bifurcation set: Hysteresis set: Double limit set: Figure 12 shows the transition sets of the system at θ = 0. The transition sets divide the parameter plane into four persistent regions. Figure 13 shows the projections of the bifurcation curves corresponding to different parameter regions on the planes of λ − u and λ − v. To illustrate the difference between the bifurcation curves, we define three key points as a, b, and c, which are the turning points of the curves. As can be seen, the system has one solution curve corresponding to the corresponding parameters in regions A and B. A single solution for u and v corresponds to any value of λ in region A, while in region B, between points a and b, u and v have three solutions. The system has two solution curves corresponding to the parameter regions C and D. In region C, point b is to the left of a, and there are five solutions between a and b. In region D, point b is to the right of point a, and a single solution exists between a and b.  Figure 14 shows the amplitude-frequency curves in the driving direction corresponding to different parameter regions. For clarity, the physical amplitudes are used instead of the non-dimensional amplitudes, where A d = A 1 x 0 and A s = A 2 x 0 . The numerical results of the Runge-Kutta method are also given. They indicate that the jump phenomenon does not occur in region A. For region B, the response jumps down at point a as the frequency increases and jumps up at point b when the frequency decreases. In region C, the response pulls in after it jumps up with the frequency increase. Furthermore, the response jumps to another solution curve when the frequency is reduced. In region D, the jump of the response at both the frequency's increasing and its decreasing leads to pull-in. The open-loop gyroscope cannot work stably in multi-solution regions because it may result in an incorrect output as the system is disturbed. Note that the peak parts of the amplitude-frequency curves of regions C and D are included in the multi-solution region. Furthermore, pull-in may occur with the change in the frequency. Thus, only regions A and B are feasible for the gyroscope operation.  Figure 14. The amplitude-frequency curves in the driving direction correspond to different persistent regions in Figure 12, where SS is the stable solution, US is the unstable solution, BP denotes the bifurcation point, A p is the pull-in boundary amplitude, NSFI is the numerical solution as the frequency is increased, NSFD is the numerical solution as the frequency is decreased; (A): Figure 15 presents the transition sets on the V D − V A plane corresponding to different inclination angles. Compared with θ = 0, both the bifurcation set and the hysteresis set shift to the bottom left at θ = 1 • and 2 • . The gyroscope's feasible region decreases with the inclination angle's magnitude. As the inclination angle is increased, the double limit point set gradually shifts to the top right.

The Influence of Parameters on the Response of the System
The amplitude of the response is meaningful for the operation of the gyroscope. This section focuses on the effect of the parameters in the feasible region (regions A and B) on the amplitude of the response. Figure 16 presents the amplitude-frequency curves corresponding to different DC and AC voltages. We notice that the amplitude-frequency curves of cases 2 to 4 are almost the same. Since there is V D V A = 100V 2 corresponding to these three cases, we suppose that the product of V D and V A is the main factor affecting the response. To verify the conjecture, Figure 17 gives the amplitude-frequency curves with different values of V D V A . The corresponding parameters are shown in Table 2. As can be seen, the amplitude-frequency curves are approximate, corresponding to the same inclination angle and V D V A values. Therefore, V D = V A = V E is used instead of other cases with the same product of V A and V D in the following.  Figure 16.
The amplitude-frequency curves of the driving direction, (b) The amplitude-frequency curves of the sensing direction.  Figure 18 shows the amplitude-frequency curves corresponding to different voltages and inclination angles. It can be seen that as the voltage V E is increased, the response amplitude increases and exhibits a stronger hardening characteristic. The natural frequency of the gyroscope decreases with the inclination angle due to the reduction in the size of the cross-section of the elastic beams as well as the bending stiffness. The amplitude of the response curves decreases first and then increases with the inclination angle. The maximum amplitude in the driving direction occurs at θ = −2 • , while the maximum amplitude in the sensitive direction appears at θ = 2 • .
We are also concerned about the mechanical sensitivity and nonlinearity of the gyroscope. These can be obtained from the linear least-squares fitting of the relation curve of the sensing amplitude and the angular velocity [9,22], as shown in Figure 19. We set the operating range of the gyroscope as ω z =0 to 50rad/s. The fitting curve is given as where S denotes the mechanical sensitivity. The nonlinearity can be given as

The Effect of the Environmental Temperature on the Mechanical Sensitivity
The change in environmental temperature may affect the output of the gyroscope. Figure 22 gives the response curve of the sensing direction corresponding to different environmental temperatures. We also set the linear natural frequency of 300 K as the working frequency. The temperature has a very small effect on the response curve, and the variation in the peak frequency of the curves is less than 1%. However, the temperature significantly affects the output amplitude of the gyroscope (shown by the red dot in Figure 22). Compared with the output of the room temperature (300 K), the output amplitude increased by 10% at 360 K and decreased by 21.5% at 240 K. Figure 23 presents the variation in the mechanical sensitivity of the gyroscope in the whole temperature range. To quantify the full temperature stability of the mechanical sensitivity of the gyroscope, we define where S T denotes the mechanical sensitivity of the full temperature range of the gyroscope and S 300 is the mechanical sensitivity at T = 300 K. Figure 24 gives the effect of the inclination angle and voltage on the full temperature stability of the mechanical sensitivity. κ increases at first and then decreases with the voltage and the inclination angle. The mechanical sensitivity changes are smaller in the whole temperature range when the voltage is large or the inclination angle is negative.

Conclusions
This paper studies the nonlinear dynamics of a symmetric MEMS gyroscope. The inclination effects on both the support stiffness and electrostatic forces are considered. The dynamic equation is established using the Lagrangian equation. The linear and nonlinear stiffness varied with the inclination angle, and environmental temperature is studied by finite element analysis. Furthermore, a static pull-in analysis in the driving direction was performed. It was found that the influence of the inclination angle on the linear stiffness is cubic, and the effect on the nonlinear stiffness is linear. The pull-in boundary amplitude decreases with the inclination angle. Ignoring the inclination effect on the stiffness leads to the opposite result.
The steady-state periodic response was obtained using the averaging method. The two-variable singularity theory was employed to study the influence of parameters on the bifurcation characteristics in the driving direction and gave the transition sets on the DC-AC voltage plane, which divided the parameter plane into four persistent regions. The corresponding bifurcation diagrams and amplitude-frequency curves are presented. The effect of the inclination angle on the transition sets is analyzed. The results show that the peak part of the amplitude-frequency curve is contained in the multi-solution region corresponding to the two regions with large voltage, and the pull-in phenomenon may occur with the change in the driving frequency. Only the two parameter regions with low voltages meet the operational requirements of the gyroscope. Both over-and under-etching may lead to the narrowing of the feasible region for the gyroscope.
The effect of parameters on the amplitude-frequency curve and the mechanical sensitivity and nonlinearity are studied in the feasible region. It demonstrates that the amplitudefrequency curves with the same product of the DC and AC voltages are almost the same. The amplitude of the response, mechanical sensitivity, and nonlinearity increase with the voltages. The response shows stronger hardening characteristics with the rise in amplitude. The increase in the inclination angle significantly reduces the natural frequency of the gyroscope. The mechanical sensitivity and nonlinearity decrease first and then increase with the inclination angle. The full temperature stability of the mechanical sensitivity is also considered. The variation in mechanical sensitivity is small at a large voltage and small inclination angle. Although the mechanical sensitivity of under-etching and over-etching is approximate, under-etching is more beneficial to the performance of the gyroscope for the small nonlinearity and good temperature stability.