Nonlinear Dynamics of Electrostatic Comb-drive with Variable Gap Under Harmonic Excitation

This paper provides an extensive study of the nonlinear dynamics of a variable gap electrostatic comb-drive. The amplitude- and phase-frequency response, as well as the amplitude- and phase-force response of the comb-drive were obtained and analyzed with and without taking into account the cubic nonlinearity of the suspension. A signiﬁcant variation in the frequency and force response is demonstrated in the presence of nonlinearity of the elastic suspension. Using numerical methods of bifurcation theory, solutions are obtained that correspond to the resonance peak of the frequency response when the constant and variable components of the voltages change. The result obtained makes it possible to determine the range of excitation voltage values that provide the required vibration amplitude in the resonant mode. The inﬂuence of the second stationary electrode on the dynamics of the system is estimated. The signiﬁcant inﬂuence of this factor on the resonant-mode characteristics is revealed.


Introduction
The development of nano-and microsystem technology (NMST or, in other words, nano-and micro-electromechanical systems -NEMS / MEMS) was initially closely associated with the development of semiconductor technologies. From a manufacturing point of view, the very fabrication technology based on selective etching/deposition of silicon has evolved along with semiconductors since the 1950s, but the use of such technologies for the manufacture of microscale devices with moving sensing elements -MEMS sensors -did not begin until the mid-1970s. -1980s [Iannacci]. Manufacturing techniques that use anisotropic etching to produce various three-dimensional suspended structures from a silicon substrate have led to the implementation of miniature pressure sensors, accelerometers, gyroscopes, switches and other devices for various applications, for example, in the optical and biomedical fields. A further impetus in the development of microsystem technology is associated with the improvement of surface micromachining technology and today MEMS technologies are used in all spheres of human life from rocketry to inkjet printing technologies [Varadan, Raspopov].
Any MEMS sensor includes a sensing element, a system of elastic suspensions to ensure the mobility of the sensing element, as well as a system for exciting oscillations and picking up a signal by electrostatic, piezoelectric, magnetic or thermal interactions. Electrostatic excitation is the most commonly used method because it is easy to implement and compatible with Complementary metal-oxide-semiconductor (CMOS) circuits. This type of actuation of a micro-electromechanical system is based on the attractive force between plates with different electric charges and can be divided into two subtypes -movement perpendicular (variable gap) and parallel (variable area) to the electrode plane [Acar]. The advantage of variable gap drives over variable area ones lies in the magnitude of the generated electrostatic force -when the gap changes, it is much greater. However, the disadvantage of such a system is the restriction on the displacement of the movable electrode -when the displacement exceeds one third of the interelectrode gap, the electrodes collapse, and the system fails.
When modeling drives with variable gap, due to the nonlinearity of the dependence of the electrostatic force on the displacement of the moving electrode [Acar], it is necessary to use asymptotic methods and conduct a qualitative analysis of the system behavior.
When investigating the nonlinear dynamics of a variable gap electrostatic drive, the mathematical model usually takes into account the cubic nonlinearity of the elastic suspension [Zhang1, Ilyas, Rhoads, Hajjaj, Han] and the viscous friction (Squeeze Film Damping) [Zhang1]. For the analysis, asymptotic methods of nonlinear dynamics are used, in particular, the multiple scales method [Zhang1, Ilyas, Hajjaj, Han]. The expansion of the expression for the electrostatic force in the Taylor series is also used [Zhang1]. Sometimes for a more precise description, the expansion in the Taylor series is avoided by using the averaging theorem together with the residue theorem to obtain an approximate analytical solution [Zhang2].
As a result of the application of the above methods, it becomes possible to obtain and then analyze the amplitude-frequency response (AFR) of the system when considering the main and/or subharmonic resonances. For example, [Zhang1] describes the influence of the parameters of the constant and variable voltage components on the character of the frequency response and the value of natural frequencies.
An electrostatic drive with a variable gap can be used not only as a excitation system of the device. For example, [Shkel] demonstrates the effect of amplitude modulation on the Q-factor and noise level in the output of a micromechanical gyroscope. Electrothermal modulation allows you to separate the useful signal from the parasitic ones in the frequency domain. It is shown that the carrier signal of an alternating current is an accurate method for adjusting the Q-factor of an micromechanical gyroscope, and it is also demonstrated that with the correct setting of the parameters, the signal-to-noise ratio can be improved by a factor of 6.
The work [Pistorio] proposes a model of an LL-type dual mass gyroscope, in which the drive system is an electrostatic actuator with a variable area, but a variable gap drive has also found its application -it is used to adjust the resonant frequency along the sensitivity axis closer to the midrange along drive axis to provide the desired difference between them with temperature interference and manufacturing errors.
The purpose of this work is to study the nonlinear dynamics of an electrostatic comb drive with a variable gap to establish the dependence of the oscillation amplitude of the moving electrode on the constant and variable volt-age components, taking into account the specified geometric and mass-inertial characteristics.
Relevance of the work is due to the intensive development of nano-and microsystem technology, leading to the development and manufacture of new and more complex micromechanical sensors for various purposes, requiring the construction and study of mathematical models of MEMS elements, including electrostatic actuators [Hajjaj]. In particular, there is a need for an accurate assessment of the comb-drive excitation voltages values, depending on electromechanical characteristics of the device, at the stage of it's initial design (namely, micromechanical vibration gyroscopes and resonant, as well as mode-localized, accelerometers) [Zotov,Lukin,Lukin2] The novelty of the work lies in a detailed qualitative study of resonant operating mode of an electrostatic drive with a variable gap, including taking into account the second stationary electrode, which requires the joint finding of solutions to both the problem of nonlinear statics and nonlinear dynamics of a moving mass.
The work is divided into four sections. A mathematical model and transformation of equations for further application of asymptotic methods of nonlinear dynamics are given in Section 2. Further, in Section 3, systems of equations in slow variables obtained after applying the multiple scales method are presented. Finally, Section 4 is a presentation and analysis of the results obtained.

Mathematical model
A model of an electrostatic drive with a variable gap, consisting of a movable and two stationary electrodes with a gaps d 1 and d 2 between them, shown in the figure 1, is considered. The dynamics equation of an electrostatic drive considering the nonlinearity of an elastic suspension: where x − vertical displacement of movable electrode, m − mass of movable part of device, µ − damping coefficient, k − linear stiffness of elastic suspension, k 3 − cubic nonlinearity coefficient, N − number of electrodes in the comb-drive, F − electrostatic force.
The electrostatic force of attraction of one pair of combs, arising between the plates, in case of coincidence of the direction of action of the force and positive direction of the axis x will take the form: where V DC − amplitude of constant component of voltage, V AC and ω − amplitude and frequency of variable component of voltage, ε 0 − dielectric constant, ε r − dielectric constant of the medium in the gap, S = w · L− electrode's area, where w and L− electrodes width and length respectively, d 1 − air gap between movable and top stationary electrodes (in side of which the displacement of movable electrode is a positive), d 2 − air gap between movable and bottom stationary electrodes.
Moving on to dimensionless parameters the equation of motion is transformed to the form Unknown function of displacement u( τ ) is decomposed into static displacement of the electrode u 0 under the action of time-independent components of electrostatic force and a dynamic time-dependent term ξ( τ ): The static displacement of the electrode is determined from the equation of static equilibrium: As a result of substituting the expansion of the displacement function into the equation of motion, an equation for the dynamic term ξ( τ ) will be obtained: where() − derivative with respect to dimensionless time d d τ . Afterwards it is necessary to exclude the static part (4) from this equation.

Transformation of the equation of motion for the application of MSM
The resulting equation (5) is nonlinear and cannot be qualitatively investigated without the use of asymptotic methods of nonlinear dynamics. In this work, the multiple-scale method (MSM) is used [Nayfeh], but to apply it, it's necessary to transform the right-hand side of the equation so that it does not contain terms with an unknown function ξ( τ ) in the denominator. In this regard, further transformations of the resulting equation will be carried out in two ways: 1. multiplication of the whole equation by the denominator of right-hand side; 2. expansion of a fraction in a Taylor series up to the third order.
Also, in addition to the complete formulation of the problem, a simplified one will be also considered -without taking into account the influence of the second stationary electrode.
The equation of static equilibrium in the absence of a second stationary electrode is: and also the equation of motion (5): Further, the equation (7) will be reduced to a more analytically convenient form for applying the MSM in two ways -by multiplying by the denominator of the right-hand side and using the Taylor series expansion.

Dynamics equation: multiplication by the denominator
After multiplying the equation (7) by the denominator of the right-hand side, it becomes: From this equation it is necessary to exclude the terms representing the equation of static equilibrium (6). Also, for ease of use in the future of the multiple-scale method, another change of the time variable was made in order to obtain the dimensionless frequency of the linear part of the equation, equal to unity: After all the transformations and grouping of the terms of the equation by powers of ξ, the final equation of motion was obtained: ξ + C 1ξ + ξ − C 2 ξξ − C 3 ξξ + C 4 ξ 2ξ + C 5 ξ 2ξ + C 6 ξ 2 + C 7 ξ 3 + C 8 ξ 4 + +C 9 ξ 5 = C 10 λ 1 sin ητ − C 10 λ 2 cos 2ητ, where the coefficients C i , i = 1 : 10 are

Dynamics equation: Taylor series expansion
The expression for the electrostatic force in the absence of the influence of the second stationary electrode has the form The factor in front of the bracket can be expanded into a Taylor series in function ξ at the point ξ 0 = 0 up to the third order inclusive: After substituting this expansion into the original equation (7) and excluding the terms that make up the equation of static equilibrium (6), similarly to the previous case, the time variable was also replaced: After all the transformations and grouping of the terms of the equation by powers of ξ, the final equation of motion was obtained, suitable for applying the multiple-scale method where the coefficients C i , i = 1 : 7 are Due to the fact that when the second fixed electrode is taken into account, the right-hand side of the equation of motion contains several fractions with two different denominators, in this problem only the method with the expansion of the electrostatic force in the Taylor series will be considered.
The derivation of the equation of motion is similar to the derivation in the absence of the second fixed electrode in the model. The equation of motion has the form (14), and the dimensionless frequency after changing the time variable, and the coefficients C i , where i = 1 : 7, are

Multiplication equation by the denomanator of righthand side
Terms in equation (10) are scaled as follows: where ε − small parameter. The multi-scale decomposition [Nayfeh] will be performed up to third order inclusive. Therefore, before the main member of electrostatic force (when considering the main resonance) there is a small parameter in third degree. Also terms with the first time derivativeξ fall into the third order of decomposition. The double-frequency term of electrostatic force must fall into the first-order solution, so a small parameter is also added before it.
A three-term expansion is considered in the form where T n = ε n τ . Time derivatives are defined as: After substituting the expansion (16) into the equation (15) and balancing in powers of a small parameter, the following expressions were obtained: where Λ = C10 λ2 1−4η 2 , A = A(T 1 , T 2 ) = 1 2 a exp{iβ}, where a and β− slow amplitude and phase, and cc − complex conjugate members.
The main resonance is subject to consideration, so After carrying out the further procedure of the MSM, a system of equations in slow variables was obtained where a and γ − dimensionless amplitude and phase, γ = β − τ σ, σ− frequency offset parameter, which is determined by formula (20) for ε = 1, () ′ − derivative with respect to dimensionless time.
The established solution is then sought from the conditions a ′ = 0 and γ ′ = 0. Substituting them into equations (21) and (22), a system of algebraic equations is obtained, the solutions of which are stationary amplitude and relative phase of oscillations.

Taylor series expansion for the electrostatic force
Terms in equation (14) are scaled as follows: After substituting the expansion (16) into the equation (23) and balancing in powers of a small parameter, the following expressions were obtained: After carrying out the further procedure of the MSM, a system of equations in slow variables was obtained When considering the second stationary electrode in the model, the system of equations in slow variables has a form similar to (27), (28) with the coefficients C i , i = 1 : 7 described above.

Validation of equation systems in slow variables
The systems of equations in slow variables obtained above were validated by comparing the solutions obtained by numerically integrating the systems of equations in slow variables and the initial ones for a simplified formulation of the problem without taking into account the second fixed electrode. The assigned parameter values are shown in the table. Results of comparing solutions for systems of equations in slow variables and initial ones for different values of parameters V DC , V AC and σ are presented below.
Based on the graphs 2 and 3, we can conclude that the results of numerical integration of systems of equations in slow variables are in good agreement with the initial ones for various values of the parameters V DC , V AC and σ for a Taylor series expansion method.
For some values of the parameters of the constant and variable components of voltages and frequency detuning, the method of the denominator multiplication gives a slowly varying amplitude less or more in value than when numerically integrating the original equation. This effect is associated with the implementation of the MSM in this case, and more specifically, caused by the fact that the expansion is implemented up to 3 orders in the small parameter ε inclusive. In this regard, some terms, including the damping coefficient, do not completely fall into the equations of the considered orders, which causes some difference in the value of the amplitude. In this regard, in further calculations, only the method with the expansion of the electrostatic force in the Taylor series will be used.  On the graph 4 it can be seen that with an increase in the amplitude of the impact, the difference between the steady-state oscillation amplitudes for the two methods increases.

Static Equilibrium Diagrams and Static Displacement Function Approximation
Since the purpose of this work is to study the dynamic modes of operation of an electrostatic drive when changing the active parameters of voltages V DC , V AC and frequency tuning parameter σ, then all other parameters depending on the above should be explicitly expressed through them in order to correctly describe the influence of active parameters changes on the behavior of the system. Among the parameters that depend on active ones, one can name the dimensionless amplitudes of the electrostatic force λ, λ 1 and λ 2 , dimensionless parameter Ω, arising as a result of repeated replacement of the time variable, as well as the static displacement u 0 of the movable electrode under the action of only the time-independent component of the electrostatic force, which is included in the expression for the Ω parameter and in all the coefficients C i introduced earlier in Sect. 3. All of the above parameters have an explicit expression of dependence on active parameters, with the exception of the static displacement u 0 , which, if the second electrode is not taken into account, can be found from the equation (4): Static diagrams in dimensional and dimensionless form for various values of the parameter k 3 as a percentage of n % to the linear stiffness of the suspension at the maximum displacement value before buckling (in force equivalent): are presented in Figures 5 and 6. With a nonzero value of the cubic stiffness parameter of the elastic suspension k 3 , to find the value u 0 depending on λ, it is necessary to solve an equation of the 5th degree, for which there is no analytical solution.
In this regard, the function u 0 (k 3 , λ) was approximated by a polynomial function of two variables. The approximation results are shown in the figure 7. Taking into account the influence of the second stationary electrode, the static equilibrium equation of the electrostatic drive system is modified (equation (4) instead of (6)). In this regard, the function which approximates the static displacement of the movable electrode depending on the parameters of the static voltage and cubic nonlinearity of the elastic suspension was modified.
The figures 8 and 9 below show static diagrams in dimensional and dimensionless form for various values of the gap ratio parameter p.

Analysis of the amplitude-frequency response of the system
As mentioned above, in the obtained systems of equations in slow variables, there are three active parameters -V DC , V AC and σ. In this subsection, the continuation of the solution by the parameter σ will be considered, that is, as a result, the amplitude and phase-frequency response (AP F req R) of the system will be obtained with different variations of the remaining parameters -V DC and V AC for a system of equations in slow variables obtained by expanding the electrostatic force in a Taylor series.
Hereinafter, continuation by parameter implies continuation of a stable equilibrium position when one of the active parameters changes. Continuation by parameter is implemented using the MATCONT package [Dhooge].
a. Analyzing the above graphs, several conclusions can be drawn. First, in the absence of cubic suspension nonlinearity, which makes the system more rigid and tilts the frequency response to the right, nonlinearity in the system is present only due to the presence of a nonlinear electrostatic force. It is known that this kind of nonlinearity reduces the stiffness of the system ("spring-softening") [Zhang1, Ilyas, Sharma] and the frequency response under the influence of electrostatic forces tilts to the left, which is observed in the figures 10, 11.
When the parameter V AC is varied upward, the amplitude of movable electrode displacement increases, as well as when the parameter V DC is varied.
Also, the graphs show a clear dependence of resonant frequency of the device on voltage V DC when constant voltage parameter increases, the frequency decreases. This effect can be used for tuning the resonance frequency on required values in some applications. a.
Analyzing the dependences 12 -13 in comparision with 10 -11, it can be concluded, that their clear difference is a slope AP F req R to the right due to addition a cubic component of the elastic suspension stiffness, which increases the overall rigidy of the system. The graph 14 below shows a series of frequency response for different values of the parameter V DC at a fixed value of the parameter of the cubic nonlinearity of the suspension, demonstrating the change in the slope of the frequency response with an increase in the influence of the "soft" nonlinearity of the electric field. Figure 14: AP F req R of the system by varying the parameter V DC , when The behavior of AP F req R depending on the changing of constant and variable components of voltage is similar to the case with absence of cubic nonlinearity -the vibrations amplitude increases with growth of voltage amplitude, and frequency decreases with increase of V DC .
Below are presented the AP F req R of the system, taking into account the influence of the second stationary electrode.
a. Analyzing the figure 15 (a), we can conclude that when varying the value of the variable voltage component parameter V AC in the absence of cubic nonlinearity of the elastic suspension, the presence of a second stationary electrode does not qualitatively affect the character of the frequency response, but only decreases the amplitude value, and the greater V AC , the greater the influence of the second stationary electrode.
When varying the constant voltage component V DC in the absence of cubic nonlinearity (Figure 15 (c)), taking into account the second stationary electrode leads not only to a decrease in the amplitude at the frequency response, but also to a shift in the resonant frequency to the left (frequency become smaller).
Taking into account the cubic nonlinearity of the suspension ( Figures  15 (b) and (d)), the previously noted changes, introduced by presence of the second stationary electrode, are amplified. Amplitude and resonance frequency are reduced by greater values than in the absence of cubic nonlinearity. Figure 16 below shows the AP F req R of the system with varying the ratio of the gaps between the movable and two stationary electrodes.  Figure 16: AP F req R of the system with varying the ratio of the gaps between the electrodes , k 3 = 1% Thus, at this stage, it can be concluded that it's necessary to take into account the second stationary electrode when studying the nonlinear dynamics of an electrostatic drive with a variable gap due to significant changes in the key characteristics of the system when considering a model with two electrodes.

Analysis of the amplitude-force response of the system
In this subsection, the continuation of the solution by the parameter V AC will be considered, that is, as a result, the amplitude and phase-force response (AP F orce R) of the system will be obtained with different variations of the remaining parameters -V DC and σ. a. In the presence of cubic nonlinearity, the frequency response of the systems has a slope to the right, and the AP F orce R graph shows the change in displacement when the parameter V AC changes at a certain frequency, that is, when moving along the frequency response from bottom to top along a vertical straight line passing through the selected frequency value (parameter is frequency detuning from resonance). Moreover, the lower branch of the frequency response (on the interval where three displacement values correspond to one frequency value) is unstable and the modes it describes are not realized. In this regard, on the AP F orce R graphs for values of the frequency detuning parameter greater than 0 (in the region to the right of the resonance), an inflection in the dependence should be observed with a sharp increase in the amplitude and subsequent smooth growth with an increase in the parameter V AC , which is observed in Figures 17 (а) and 18 (a), and the closer to resonance, the lower voltage values the described effect will be observed.
In the zone before resonance, the dependence of the frequency response is unambiguous, and therefore the graphs of the AP F orce R should show just a smooth increase in the amplitude with increasing voltage, which is demonstrated in Figures 17 (b) and 18 (b).
Further, the AP F orce R of the system will be obtained, taking into account the influence of the second stationary electrode.
a. The analysis of the above graphs leads to the conclusion that outside the resonance zone (at σ = −0.005) the influence of the second electrode on the AP F orce R of the system is limited by a decrease in the amplitude. However, in the over-resonance zone (σ = 0.005), the influence of the second electrode on the slope of the ACH curve is seen.

Continuation of the resonant-mode solution
The next step is to study the dependence of the resonance point, which describes the desired operating mode of the electrostatic drive, on the values of the voltage parameters -V AC and V DC . From a mathematical point of view, we perform continuation of the bifurcation point of the fold type ("limit point") [Kuznetsov] by the frequency detuning parameter when increasing or decreasing two other active parameters in turn. Figure 20 demonstrates series of AP F req R of the system for various values of the parameter V AC . Red dashed line on this graph demonstrates changing of bifurcation point by σ with increasing V AC . It must be noted, that when moved along this line, not only value of frequency tuning parameter is changing, but also the vibrations amplitude. Hence, the graph demonstrating the dependence of resonant amplitude and frequency on active parameters V AC and V DC , will define the values of excitation voltages required to maintain the vibration at given resonant frequency with given amplitude, that is a main purpose of this work.

Resonant regime dependence on V AC
In this subsection the dependencies of resonant frequency and amplitude on parameter of the variable voltage component V AC are presented. Calculations are performed using slow-variable system obtained for Taylor series expansion method.  Figure 21 on the left shows the dependence of the resonant frequency on V AC . This dependence has a limiting value in the region near zero according to the parameter V AC , after which it changes its growth pattern. This effect is due to the fact that not for all values of active parameters, a fold-type bifurcation point will be observed on the frequency response, since at low voltage values the frequency response function is unambiguous and does not have an unstable branch. In this connection, with the continuation of the bifurcation point, no solutions are observed in this area. This effect is also observed in the figure 21 on the right, which shows the dependence of the resonance amplitude on V AC . Figure 22 below shows the frequency response for different sets of parameter values, chosen so that one of them falls into the zone of existence of the solution in the figure 21, and the other into the zone where no solutions are observed . It can be seen that at small values of the parameter V AC , the dependences obtained on the linear and nonlinear models coincide, and with an increase in the amplitude of the AC voltage, the nonlinear model makes it possible to more accurately determine the dependence of the moving electrode amplitude of oscillations on the control voltage value.
Of interest is the three-dimensional graph of the dependence of the resonant amplitude and frequency on the active parameter V AC , which is shown below in the figure 23. Similarly to the conclusion from two-dimensional plots, of the two branches in the figure 23, the lower one is realized, which is characterized by a rapid increase ξ with an increase V AC . At a given value of the parameter V DC , according to these dependencies, it is possible to determine the value of the variable voltage component to ensure a given displacement of the movable electrode.
Below are the dependences of the resonance frequency and amplitude when changing the parameter of the variable voltage component V AC , taking into account the influence of the second stationary electrode.
a. The effect of the second stationary electrode is expressed in a decrease in resonant amplitude and frequency.
Also below is a three-dimensional graph of the dependence of the amplitude and frequency on the variable voltage component.

Resonant regime dependence on V DC
In this subsection the dependencies resonant frequency and amplitude on parameter of the constant voltage component V DC are presented. Calculations are performed using slow-variable system obtained for Taylor series expansion method.
It should be noted that the parameter V DC , along which the continuation is carried out, is also included in the equation of static equilibrium (4) and (6), the solution of which is the static displacement of the electrode. Here, while continuing the dynamic stationary regime, the solution also varies (does not remain equal to the original value). The continuation is carried out due to the explicit dependence of the static solution on the constant voltage component, obtained as a result of the approximation of the electrode displacement function in Section 4.2.  Figure 26 on the left shows the dependence of the resonant frequency on V DC . This dependence has a limiting value in the region near zero according to the parameter V DC , after which it changes its character of growth. This effect is also due to the fact that not for all values of active parameters, a fold-type bifurcation point will be observed on the frequency response. This effect is also observed in the figure 26 on the right, which shows the dependence of the resonant amplitude on V DC .
Also of interest is the three-dimensional graph of the dependence of the resonant amplitude and frequency on the active parameter V DC , which is shown below in the figure 27. Similarly to the conclusion from two-dimensional plots, of the two branches in the figure 27, the lower one is realized, which is characterized by a rapid increase ξ with an increase V AC . At a given value of the parameter V AC , according to these dependencies, it is possible to determine the value of the constant voltage component to ensure a given displacement of the movable electrode.
In conclusion, similar dependences are presented when taking into account the second stationary electrode when varying the ratio of the gaps between the electrodes. a.

Conclusion
In this paper, an extensive study of the nonlinear dynamics of an electrostatic comb-drive with a variable gap is carried out. The amplitude-and phase-frequency response and the amplitude-and phase-force response of the comb-drive were obtained and analyzed with and without regard to the cubic nonlinearity of the suspension. A decrease in the stiffness of the system under the action of electrostatic forces ("spring-softening"), as well as a significant variation of the frequency and force response in the presence of nonlinearity of an elastic suspension was demonstrated. Using numerical methods of bifurcation theory, solutions are obtained that correspond to the resonance peak of the frequency response when the constant and variable components of the voltages change. This results makes it possible to determine the range of excitation voltages values that provide the required vibrations amplitude in resonant regime. The influence of the second stationary electrode on the dynamics of the system is estimated. The significant influence of this factor on the resonant-mode characteristics is revealed.

Conflict of Interest
The authors declare that they have no conflict of interest.