`Mathematical Problems in EngineeringVolume 2011, Article ID 173459, 14 pageshttp://dx.doi.org/10.1155/2011/173459`
Research Article

## Homotopy Analysis Method for Nonlinear Dynamical System of an Electrostatically Actuated Microcantilever

1State Key Laboratory of Mechanical System and Vibration, Shanghai Jiao Tong University, 800 Dongchuan Road, Shanghai 200240, China
2Department of Mechanics, Sun Yat-Sen University, Guangzhou 510275, China

Received 21 May 2010; Revised 8 November 2010; Accepted 30 November 2010

Academic Editor: Marcelo Messias

Copyright © 2011 Y. M. Chen et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The homotopy analysis method (HAM) is employed to propose an approach for solving the nonlinear dynamical system of an electrostatically actuated micro-cantilever in MEMS. There are two relative merits of the presented HAM compared with some usual procedures of the HAM. First, a new auxiliary linear operator is constructed. This operator makes it unnecessary to eliminate any secular terms. Furthermore, all the deformation equations are purely linear. Numerical examples show the excellent agreement of the attained solutions with numerical ones. The respective effects of applied voltage, cubic nonlinear stiffness, gap distance, and squeeze film damping on vibration responses are analyzed detailedly.

#### 1. Introduction

Over the last couple of decades, Liao  described a nonlinear analytical technique which does not require small parameters and thus can be applied to solve nonlinear problems without small or large parameters. This technique is based on homotopy, which is an important part of topology, called the homotopy analysis method (HAM). Its main idea is to construct a class of homotopy in a rather general form by introducing an auxiliary parameter. This parameter can provide us with a convenient way to control the convergence of approximation series and adjust convergence rate and region when necessary. A systematical description of this method was presented in . Liao  also studied the convergence properties of the HAM and proved that as long as an HAM series is convergent, it must converge to one solution of the considered problem.

Many microelectro-mechanical systems (MEMSs) inherently contain nonlinearities, such as intrinsic and exterior nonlinearities arising from coupling of different domains . Also, there exist mechanical nonlinearities, for example, large deformations, surface contact, creep phenomena, time-dependent masses and nonlinear damping effects, and so forth [9, 10]. It seems fair to say that nonlinear dynamic analysis becomes an increasingly important task in MEMS research and manufacturing.

The nonlinear dynamical behaviors of microcantilever-based instrument in MEMS under various loading conditions have stimulated the curiosities and interests of many researchers . For example, the oscillation of an electrostatically actuated microcantilever-based device in MEMS was investigated through a simplified mass-spring-damping model subjected to nonlinear electrostatic force . Complex nonlinear terms arising from electrostatic force and from squeeze film damping make it difficult to analyze the system directly using some routine techniques for nonlinear vibrations. For this reason, Zhang and Meng  suggested an approximate treatment by expanding the nonlinearities into Taylor series and retaining only the first two terms. The harmonic balance method was then applied to solve the approximate system. Note that the attained results have not been compared with any numerical solutions. The purpose of this study is to seek highly accurate solutions of the aforementioned system on the basis of the HAM.

#### 2. Dynamical System

The electrostatically actuated microcantilever in MEMS as shown in Figure 1 is in dimensions . The governing equation describing the motion of the microstructure is [16, 20] where the superscript denotes the differentiation with respect to time the vertical displacement of the microcantilever relative to the origin of the fixed plate, the mass, and the effective spring stiffness and damping coefficient of the simplified system, respectively. According to the parallel plate theory, if the fringe effects at the edges of the plates are ignored , the electrostatic force () generated by applying a voltage between the capacitor plates (the fixed plate and the movable plate) can be expressed by where is the absolute dielectric constant of vacuum,  N/m, the overlapping area between the two plates, and (μm) is the gap between them. The other parameter values are given as  kg, and  N/m,  kg/s and  m2 , while the value of varies.

Figure 1: A simplified dynamical model of the microcantilever in MEMS.

Introducing the following dimensionless variables then one can rewrite (2.1) and (2.2) as where , , and the superscript () refers to the differentiation with respect to .

When the applied voltage, , includes alternating current (ac) voltage, , and polarization voltage, , can be considered for analyzing the nonlinearities of the system. In addition, MEMS devices are often characterized by structures that are a few microns in size, separated by micron-sized gaps. At these sizes, air damping dominates over other dissipation. Squeeze film damping may be used to represent the air damping experienced by the moving plates . When considering the effect of squeeze film damping and a cubic nonlinear spring stiffness , system (2.4) is given as where is a nondimensional parameter, and is the squeeze film damping (Pa·μm4·s).

#### 3. Mathematical Formulation

In this section, the HAM is applied to solve (2.5). First of all, introducing a new time scale , and multiplying (2.5) by , one can rewrite (2.5) as where the superscript denotes the differentiation with respect to . Compared with (2.5), it is very easy to expand (3.1) once is expressed as an HAM series.

##### 3.1. Basic Idea of Homotopy Analysis Method

In the frame of the HAM, the first procedure is to choose an auxiliary linear operator and a nonlinear operator, respectively. A new equation, usually called as the zeroth-order deformation equation, is constructed as where denotes the chosen auxiliary operator, the nonlinear one, the embedding parameter and a nonzero auxiliary parameter. The HAM is based on the continuous variation of . When the embedding parameter increases from 0 to 1, varies from an initial guess to one exact solution of the considered problem According to (3.1), it can be given as When , (3.2) becomes a linear one, whose exact solution is evident. When , it is the same as (3.1) provided that . Expanding into a Taylor series Note that is dependent upon . In fact, the auxiliary parameter is introduced to make it convenient to control and adjust the convergence of the HAM series (i.e., (3.5)). As long as series (3.5) is convergent at , one can obtain the th-order HAM solution as Substituting (3.4) and (3.5) into (3.2) and equating the coefficients of to 0, one obtains where when , , otherwise ; corresponds to the expansion of , that is, , with and . Equations (3.7) and (3.8) are always linear for every positive integer . Importantly, the nonlinear equation has so far been successfully transformed into a series of linear ones. These linear equations can be solved step by step. That is core idea and also the most outstanding advantage of the HAM.

##### 3.2. A Widely Used Linear Operator

One of key procedures of the HAM is to choose an appropriate linear operator. According to Liao , it should be chosen on the basis of the so-called rule of solution expression. Periodic/limit cycle solutions can be expressed as Fourier series in Note that the HAM has been widely used to obtain periodic and/or limit cycle solutions of nonlinear oscillators . Particularly, this technique has found applications in microsystems  and microbeams subjected applied voltages , respectively. To the best of out knowledge, most researchers  adopted such a following linear operator as If (3.11) is employed in this study, (3.7) becomes one homogenous equation in . This equation has a series of solutions that has to be determined using the first-order deformation equation, that is, (3.8) with . To this end, an initial guess for is usually given as where and are prior unknown constants.

Due to the rule of solution expression, can be given as truncated Fourier series. Substituting into , one has where refers to the highest harmonic. Considering (3.7) and (3.13), and can result in the so-called secular terms and , respectively. These terms do not comply with the rule of solution expression. One usual procedure to eliminate the secular terms is equating both and to 0, that is, When , these equations are all linear. When , however, (3.14) is nonlinear. Actually, it is essentially the first-order harmonic balance equations [25, 26] governing and . Well-known, it is a major shortcoming of the harmonic balance method that nonlinear algebraic equations have to be solved.

##### 3.3. A New Linear Operator

Considering that (3.1) contains powers like , there will be and in (3.14) with . These nonlinearities can make it rather difficult to determine the initial solution. It is worthy to find another linear auxiliary operator to avoid this obstacle. Consider the following linear operator Substitution of (3.15) into (3.7) and (3.8) results into where , otherwise when . Equations (3.16) and (3.17) provide us with a purely linear algorithm, without an additional task of eliminating the secular terms as long as . Since (3.16) is an inhomogeneous equation whose exact solution can be determined analytically, there is also no need to choose an initial guess for .

##### 3.4. A Means for Saving Computational Effort

The number of harmonics in the HAM approximations increases acceleratedly, because of the existence of the powers like in (3.1). Denote the number of harmonics in as . The starting solution () provided by (3.16) contains the first two harmonics, that is, . In the equation governing , the highest powers is , hence one knows that . As increases further, and can be obtained. In most cases, periodic or limit cycle oscillations are dominated by several lower harmonics. High harmonics contribute to the improvement of the precision quite limitedly, compared with the lower ones. Therefore, it seems to be unworthy to improve the accuracy slightly with such a vast increasing of computational cost.

Therefore, the high-order HAM approximations are obtained with a given number () of harmonics, while neglecting all harmonics higher than the th one. Denote the th-order HAM approximation as where increases with increasing. Once , retain the first to the th harmonics in , and rewrite it as When seeking higher-order HAM approximations, (3.19) but not (3.18) is employed. Through this means, an HAM solution with harmonics can be obtained. Interestingly, as long as this solution converges, it must converge to one harmonic balance solution, as revealed in .

#### 4. Results and Discussion

Theoretically, analytical solutions of can be obtained by solving (3.16) and (3.17) step by step, since all equations are linear. However, their expressions can rapidly become complicated as increasing. It is rather easy to obtain semi-analytical solutions when all parameters are given. The HAM solutions are compared with numerical results. All the numerical solutions are obtained via integrating (2.5) using the RK method, with the initial values as and .

Let denote the CPU time needed to obtain the th-order HAM approximation () with harmonics, and is the counterpart when no harmonic is neglected. In addition, let be the maximum value (amplitude) of the th-order HAM solution (i.e., ) with harmonics retained, and the counterpart. From Figure 2, one can see the convergence of and to the numerical one as increasing. It implies that the accuracy will not be lost by retaining only the first 10 harmonics (). Also shown in Figure 2 is the ratio between and . It rises rapidly when increases. Therefore, a lot of computational effort can be saved without losing substantial accuracy, through neglecting some high harmonics. Note that, all the following HAM solutions are obtained by employing (3.19) rather than (3.18).

Figure 2: (a) the HAM solutions; (b) the ratio between and , , and . Parameter values are given as , , , , , and .

A preliminary procedure to obtain an HAM solution consists in the choosing of the auxiliary parameter . Liao  suggested an -curve approach for choosing a proper value of , which can guarantee both the convergence of HAM series and a relatively rapid convergence rate. Figure 3 shows the relative errors versus varying . These errors refer to the percentage-wise differences between the amplitudes () of the 20th-order HAM solutions with 10 harmonics and the numerical results. According to the figure, an appropriate range for can be chosen from −2 to −0.5 or so. When the absolute value of is too large, for instance more than 2, it can probably make the HAM series divergent. On the other hand, as long as approaches 0 enough, the HAM series converge by a relatively slow convergence rate.

Figure 3: Relative errors versus between 20th-order HAM solutions and numerical ones. Parameter values are , , , and , respectively.

Figure 4 shows a phase plane of system (2.4). The solutions provided by the HAM approach the numerical one when increases. The 20th-order HAM solution is nearly the same as the numerical one. In fact, solutions can be obtained to any desired accuracy without additional difficulty.

Figure 4: Comparison between numerical solutions and HAM results , where heavy dots denote the former and solid lines the latter. Parameter values are given as , , , , and .

When a nonlinear spring and/or squeezefilm damping are taken into account, the proposed algorithm can still provide very accurate solutions, as shown in Figures 5 and 6. Both the figures show clearly the convergence of the HAM solutions to the numerical ones. It is worthy of emphasizing that, the vibrations described by these phase planes are nonharmonic, namely the solutions cannot be approximated by the first harmonic alone. That is probably because high harmonics contribute to the responses to an nonnegligible extent. The phase plane shown in Figure 6 even has a loop. Even though, the presented HAM can still track them with excellent precision. These illustrative examples demonstrate the feasibility and efficiency of the presented approach in analyzing nonlinear dynamical behaviors.

Figure 5: Comparison between numerical solutions and HAM solutions , where heavy dots denote the former and solid lines the latter. Parameter values are given as , , , , , and .
Figure 6: Comparison between numerical solutions and HAM solutions , where heavy dots denote the former and solid lines the latter. Parameter values are given as , , , , , and .

Frequency-response curve is of substantial importance for dynamic analysis, when external excitations are included. Figures 7, 8, and 9 show the curves of the amplitudes of periodic solutions versus varying , for different values of , and , respectively. The 20th-order HAM solutions are in excellent agreement with the numerical ones. The amplitudes decrease with increasing for all the three cases. According to (2.5), and correspond to the magnitude of the external excitation, respectively. Also, is connected with the external excitation through , which implies that increases when decreases. That is why the amplitude increases when () increases or decreases, as shown in these figures.

Figure 7: The vibrations amplitudes versus , where heavy dots denote numerical results and solid lines the 20th-order HAM solutions with , , , , , and .
Figure 8: The solutions amplitudes versus , where heavy dots denote numerical results and solid lines the 20th-order HAM solutions with , , , , , and .
Figure 9: The solutions amplitudes versus , where heavy dots denote numerical results and solid lines the 20th-order HAM solutions with , , , , , and .

Figures 10 and 11 are presented to investigate the respective effects on the periodic responses of the cubic stiffness coefficient and the squeeze film damping. Nice agreement of the 20th-order HAM solution to the numerical results can also be observed. As they show, respectively, the larger or is, the smaller the amplitudes become. According to (2.5), is the nonlinear restoring force, hence can be considered as an equivalent linear stiffness coefficient. Its value keeps nonnegative though varying. The increasing of implies an increasing of the stiffness, so that the magnitude of response decreases. As concerning , it denotes the coefficient of a positive nonlinear damping. Generally speaking, positive damping can suppress the vibration to some extent.

Figure 10: The solutions amplitudes versus , where heavy dots denote numerical results and solid lines the 20th-order HAM solutions with , , , , , and .
Figure 11: The solutions amplitudes versus where heavy dots denote numerical results and solid lines the 20th-order HAM solutions with , , , , , and .

It is worthy to point out, the influence of squeeze film damping on vibration amplitudes reduces rapidly when increases. When is relative large (), the decreasing of makes a very little difference of the amplitudes.

#### 5. Conclusions and Remarks

Based on the homotopy analysis method, we have proposed an approach to solve the nonlinear dynamical system of an electrostatically actuated microcantilever. This approach was further improved by restricting the number of harmonics as a given number . It turns out that the improved approach can save a large amount of computational effort without reducing accuracy. Numerical examples validate the proposed approach. Approximations for periodic solutions are obtained very precisely. Using the presented approach, the frequency-response curves are obtained and discussed in details. The feasibility and efficiency, illustrated by the numerical examples, imply that we could expect the presented approach be applicable in more nonlinear vibration problems, especially those arising in MEMS.

As mentioned above, we have presented an auxiliary operator (3.15) that is different from another usually-adopted one (3.11). It is advantageous to use (3.15) rather than (3.11) for the considered problem. First, there is no need to choose an initial guess for the zeroth-order HAM approximation if one employs linear operator (3.15). Second, if (3.11) is adopted, a series of algebraic equations have to be solved from time to time to eliminate secular terms when seeking every order HAM approximation. In many occasions, the first equation is nonlinear. In our procedures, there is no requirement of doing so. Last but not least, the proposed procedures are purely linear.

#### Acknowledgments

This paper is supported by the National Natural Science Foundation of China (10102023, 10972241, 11002088 and 10732060), Doctoral Program Foundation of Ministry of Education of China (20090171110042), China Postdoctoral Science Foundation (20090460628), and Shanghai Postdoctoral Science Foundation (10R21413600). Project supported by the Research Fund of State Key Lab of MSV, China (Grant no. MSV-MS-2010-9).

#### References

1. S. J. Liao, “Second-order approximate analytical solution of a simple pendulum by the process analysis method,” Journal of Applied Mechanics, vol. 59, no. 4, pp. 970–975, 1992.
2. S. J. Liao, “A kind of approximate solution technique which does not depend upon small parameters-II: an application in fluid mechanics,” International Journal of Non-linear Mechanics, vol. 32, no. 5, pp. 815–822, 1997.
3. S. Liao, “On the homotopy analysis method for nonlinear problems,” Applied Mathematics and Computation, vol. 147, no. 2, pp. 499–513, 2004.
4. S. Liao, “Notes on the homotopy analysis method: some definitions and theorems,” Communications in Nonlinear Science and Numerical Simulation, vol. 14, no. 4, pp. 983–997, 2009.
5. S. J. Liao, “An explicit, totally analytic approximate solution for Blasius' viscous flow problems,” International Journal of Non-Linear Mechanics, vol. 34, no. 4, pp. 759–778, 1999.
6. S. D. Senturia, “Simulation and design of microsystems: a 10-year perspective,” Sensors and Actuators, A, vol. 67, no. 1, pp. 1–7, 1998.
7. S. E. Lyshevski, “Nonlinear microelectromechanic systems (MEMS) analysis and design via the Lyapunov stability theory,” in Proceeding of the 40th IEEE Conference on Decision and Control, pp. 4681–4686, Orlando, Fla, USA, 1997.
8. M. Moghimi Zand and M. T. Ahmadian, “Dynamic pull-in instability of electrostatically actuated beams incorporating Casimir and van der Waals forces,” Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science, vol. 224, no. 9, pp. 2037–2047, 2010.
9. S. K. De and N. R. Aluru, “Complex nonlinear oscillations in electrostatically actuated microstructures,” Journal of Microelectromechanical Systems, vol. 15, no. 2, pp. 355–369, 2006.
10. M. Ashhab, M. V. Salapaka, M. Dahleh, and I. Mezić, “Dynamical analysis and control of microcantilevers,” Automatica, vol. 35, no. 10, pp. 1663–1670, 1999.
11. A. Passian, G. Muralidharan, A. Mehta, H. Simpson, T. L. Ferrell, and T. Thundat, “Manipulation of microcantilever oscillations,” Ultramicroscopy, vol. 97, no. 1–4, pp. 391–399, 2003.
12. Y. Ahn, H. Guckel, and J. D. Zook, “Capacitive microbeam resonator design,” Journal of Micromechanics and Microengineering, vol. 11, no. 1, pp. 70–80, 2001.
13. E. K. Chan and R. W. Dutton, “Electrostatic micromechanical actuator with extended range of travel,” Journal of Microelectromechanical Systems, vol. 9, no. 3, pp. 321–328, 2000.
14. S. N. Mahmoodi and N. Jalili, “Piezoelectrically actuated microcantilevers: an experimental nonlinear vibration analysis,” Sensors and Actuators A, vol. 150, no. 1, pp. 131–136, 2009.
15. Y. M. Fu and J. Zhang, “Nonlinear static and dynamic responses of an electrically actuated viscoelastic microbeam,” Acta Mechanica Sinica, vol. 25, no. 2, pp. 211–218, 2009.
16. S. Liu, A. Davidson, and Q. Liu, “Simulating nonlinear dynamics and chaos in a MEMS cantilever using Poincare mapping,” in Proceeding of the IEEE 12th International Conference on Solid State Sensors Actuators and Microsystems (Transducers '03), pp. 1092–1095, Boston, Mass, USA, June 2003.
17. S. Liu, A. Davidson, and Q. Lin, “Simulation studies on nonlinear dynamics and chaos in a MEMS cantilever control system,” Journal of Micromechanics and Microengineering, vol. 14, no. 7, pp. 1064–1073, 2004.
18. W. M. Zhang, G. Meng, and DI. Chen, “Stability, nonlinearity and reliability of electrostatically actuated MEMS devices,” Sensors, vol. 7, no. 5, pp. 760–796, 2007.
19. G. Meng, W. M. Zhang, H. Huang, H. G. Li, and DI. Chen, “Micro-rotor dynamics for micro-electro-mechanical systems (MEMS),” Chaos, Solitons and Fractals, vol. 40, no. 2, pp. 538–562, 2009.
20. W. Zhang and G. Meng, “Nonlinear dynamical system of micro-cantilever under combined parametric and forcing excitations in MEMS,” Sensors and Actuators A, vol. 119, no. 2, pp. 291–299, 2005.
21. R. H. Price, J. E. Wood, and S. C. Jacobsen, “Modelling considerations for electrostatic forces in electrostatic microactuators,” Sensors and Actuators, vol. 20, no. 1-2, pp. 107–114, 1989.
22. T. Veijola, H. Kuisma, and J. Lahdenperä, “The influence of gas-surface interaction on gas-film damping in a silicon accelerometer,” Sensors and Actuators A, vol. 66, no. 1–3, pp. 83–92, 1998.
23. S. J. Liao, “An analytic approximate approach for free oscillations of self-excited systems,” International Journal of Non-Linear Mechanics, vol. 39, no. 2, pp. 271–280, 2004.
24. S. Liao and Y. Tan, “A general approach to obtain series solutions of nonlinear differential equations,” Studies in Applied Mathematics, vol. 119, no. 4, pp. 297–354, 2007.
25. Y. M. Chen and J. K. Liu, “Homotopy analysis method for limit cycle flutter of airfoils,” Applied Mathematics and Computation, vol. 203, no. 2, pp. 854–863, 2008.
26. Y. M. Chen and J. K. Liu, “A study of homotopy analysis method for limit cycle of van der Pol equation,” Communications in Nonlinear Science and Numerical Simulation, vol. 14, no. 5, pp. 1816–1821, 2009.
27. LI. N. Zhang and JI. H. He, “Homotopy perturbation method for the solution of the electrostatic potential differential equation,” Mathematical Problems in Engineering, vol. 2006, Article ID 83878, 6 pages, 2006.
28. M. Moghimi Zand and M. T. Ahmadian, “Application of homotopy analysis method in studying dynamic pull-in instability of microsystems,” Mechanics Research Communications, vol. 36, no. 7, pp. 851–858, 2009.
29. M. Moghimi Zand, M. T. Ahmadian, and B. Rashidian, “Semi-analytic solutions to nonlinear vibrations of microbeams under suddenly applied voltages,” Journal of Sound and Vibration, vol. 325, no. 1-2, pp. 382–396, 2009.
30. Y. M. Chen, J. K. Liu, and G. Meng, “Relationship between the homotopy analysis method and harmonic balance method,” Communications in Nonlinear Science and Numerical Simulation, vol. 15, no. 8, pp. 2017–2025, 2010.