The dynamic behavior of MEMS arch resonators actuated electrically

ARTICLE IN PRESS International Journal of Non-Linear Mechanics 45 (2010) 704–713

Contents lists available at ScienceDirect

International Journal of Non-Linear Mechanics journal homepage: www.elsevier.com/locate/nlm

The dynamic behavior of MEMS arch resonators actuated electrically Hassen M. Ouakad, Mohammad I. Younis  Department of Mechanical Engineering, State University of New York at Binghamton, Binghamton, NY 13902, USA

a r t i c l e in f o

a b s t r a c t

Article history: Received 8 June 2009 Received in revised form 16 April 2010 Accepted 20 April 2010

In this paper, we investigate the dynamic behavior of clamped–clamped micromachined arches when actuated by a small DC electrostatic load superimposed to an AC harmonic load. A Galerkin-based reduced-order model is derived and utilized to simulate the static behavior and the eigenvalue problem under the DC load actuation. The natural frequencies and mode shapes of the arch are calculated for various values of DC voltages and initial rises. In addition, the dynamic behavior of the arch under the actuation of a DC load superimposed to an AC harmonic load is investigated. A perturbation method, the method of multiple scales, is used to obtain analytically the forced vibration response of the arch due to DC and small AC loads. Results of the perturbation method are compared with those obtained by numerically integrating the reduced-order model equations. The non-linear resonance frequency and the effective non-linearity of the arch are calculated as a function of the initial rise and the DC and AC loads. The results show locally softening-type behavior for the resonance frequency for all DC and AC loads as well as the initial rise of the arch. & 2010 Elsevier Ltd. All rights reserved.

Keywords: Shallow arch MEMS Electrostatic force Vibrations Non-linear dynamics Perturbation

1. Introduction Micro-electromechanical systems (MEMS) resonators have been investigated thoroughly in the literature for their potential to build high sensitive sensors [1–11]. Examples of these include mass sensors to detect chemical and biological substances [1], force and acceleration sensors [2], and temperature sensors [3]. MEMS resonators have several advantages that are related to their fabrication technology allowing them to be compatible with the complementary metal oxide semiconductor (CMOS) processes [4,5]. This results into lower cost, low power consumption [6], increased reliability and manufacturability, and enabling single chip solutions. Papers dating back to the sixties by Nathanson et al. [7,8] have described the utilization of microbeams as resonators. Resonators made of clamped–clamped beam have been modeled and measured by various research groups [9–11]. Initially curved beams have been under increasing interest in the research community in recent years. Curved beams refer here to beams that are fabricated intentionally to be curved (arches) or made curved by buckling straight beams through compressive axial loads (buckled beams). Many groups [12–30] studied the bistability behavior of initially curved microbeams, which were found to be suitable for applications such as micro-shutter positioning, microvalves, and electrical microrelays. These beams have been proposed also as switches and actuators based on their snap-through motion. Most of the MEMS literature has been

 Corresponding author. Tel.: + 1 607 777 4983; fax: + 1 607 777 4620.

E-mail address: [email protected] (M.I. Younis). 0020-7462/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijnonlinmec.2010.04.005

focused on utilizing snap-through as a static phenomenon due to the actuation of static forces. Those forces can be mechanical [21], magnetic [22–24], thermal [25,26], and electrostatic [27–30]. When actuating a curved beam by a parallel-plate electrostatic force, the stability behavior of the arch becomes more interesting. Studies of electrostatically actuated curved beams have shown that they may exhibit snap-through buckling or pull-in instability as well as bi-stable behavior depending on the interaction between mechanical and electrostatic non-linearities [27–30]. Casals-Terre´ and Shkel [27] studied theoretically and experimentally the possibility of triggering the snap-through motion of a bi-stable electrically actuated beam driven dynamically by means of mechanical resonance. Zhang et al. [28] and Krylov et al. [29] conducted theoretical and experimental investigations of initially curved clamped–clamped microbeams actuated by DC loads. Their simulations were based on the Galerkin method and they have shown good agreement among their theoretical and experimental results. Das and Batra [30] conducted a transient analysis of curved microbeams using coupled finite-element and boundary-element methods. They have shown the softening effect of the MEMS arch may be dominant before it experiences its snap-through motion. Few groups studied the dynamic behavior of curved microbeams, for example [31–34]. Buchaillot et al. [31] explored theoretically and experimentally the dynamic snap-through phenomenon of MEMS clamped–clamped buckled beams when they are subjected to large vibration amplitudes. Cabal and Ross [32] devised a methodology to describe theoretically the snapthrough phenomena of a bilayer micromachined curved beam. Poon et al. [33] studied the dynamic buckling response to

ARTICLE IN PRESS H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

705

sinusoidal excitation of a clamped–clamped curved beam using the Runge–Kutta numerical integration method. Their work predicted various characteristic features, such as softening and hardening behaviors and chaotic motion of intermittent snapthrough. Noijen et al. [34] presented a higher-order model for fast and accurate evaluation of the non-linear dynamics of a buckled pinned–pinned beam. The dynamic behavior of macro-scale shallow arches has been longly investigated thoroughly in the structural and non-linear mechanics literature [35–41]. Humphreys [35] studied experimentally and theoretically the dynamics of a circular arch under uniform dynamic pressure loading. Hsu et al. [36] investigated the problem of clamped–clamped shallow arch under step loads and derived conditions for the stability of the considered structure. Hung [37] studied, using energy methods, the dynamic buckling of some elastic hinged–hinged shallow arches under sinusoidal load. Bi and Dai [38] investigated the dynamical behavior of a shallow arch subjected to periodic excitation. They showed that this kind of structure exhibits internal resonance and period doubling cascade bifurcations leading to chaos. Internal resonances have been also addressed for shallow arches including one-to-one, two-to-one and three-to-one internal resonances [39–41]. This paper aims to investigate the dynamics of a MEMS resonator based on an initially curved (arch) clamped–clamped microbeam, which is common in MEMS, under the effect of small electric loads. By reviewing the state of the art, we can see that the dynamic behavior of clamped–clamped shallow arches under harmonic electrostatic forces has not been investigated. In previous works [11,42,43], we used analytical and computational methods (perturbation technique, reduced-order models) to investigate the static and dynamic behavior of straight microbeams under electrostatic loads. In this work, we develop similar approaches to investigate the static and dynamic behaviors of electrostatically actuated clamped–clamped shallow arches. The motivation of this investigation is the fact that deep understanding of the dynamic behavior of these structures is needed to enable successful utilization of them as sensors and actuators and to reveal their unique dynamical characteristics. The organization of this paper is as follows. In Section 2 we present the model for an electrically actuated clamped–clamped shallow arch. In Section 3, a reduced-order model is derived. We then study the arch static problem under the effect of DC electrostatic forces in Section 4. In Section 5, the linear vibration problem (natural frequencies and mode shapes) of the arch under a DC load is solved. Section 6 deals with the perturbation analysis of the shallow arch to study its dynamics near its first natural frequency when excited by a small DC and AC loads. Finally, Section 7 provides a summary and some conclusions of the paper.

where

2. Problem formulation

w0 ðxÞ ¼

In this section, the problem governing the static and dynamic behavior of a MEMS shallow arch is formulated. We consider a

The non-dimensional parameters used in the above equations are presented in Table 1.

ˆz

clamped–clamped shallow arch, Fig. 1, of initial shape ^ ¼ b0 ½1cosð2px=LÞ=2, ^ ^ 0 ðxÞ where b0 is the initial rise, actuated w by an electrostatic force applied by a parallel-plate electrode underneath it with a gap width d using a DC load of amplitude VDC superimposed to an AC harmonic load of amplitude VAC and ~ . In this investigation, we assume a shallow arch, in frequency O ^ 0 0 {1. Hence, when actuated by electrostatic forces, the which w parallel-plate assumption can be considered valid. In another word, the axial component of the electrostatic force, due to the upper deformed electrode (the arch), is assumed negligible. This assumption however may not be valid for deep arches. Therefore, assuming an Euler–Bernoulli beam model, the non-linear equation ^ of the arch of ^ x, ^ tÞ of motion governing the transverse deflection wð width b, thickness h, and length L is expressed as [29,44] EI

^ @4 w

^ @2 w

^ @w @t^ #" Z ( 2  ) # ^ ^ dw ^0 ^ ^ 0 EA L @2 w d2 w @w @w dx þ þ2 ¼ 2 2 2L 0 @x^ @x^ dx^ dx^ @x^ eb½VDC þ VAC cosðO~ tÞ2  ^ 0 þ wÞ ^ 2 2ðd þ w

4 @x^ "

þ rA

2 @t^

þ c~

where E is its effective Young’s modulus, A¼bh its cross-sectional area, I its moment of inertia, r is the material density, c~ is the viscous damping coefficient (here assumed to be linear since the gap width is large, which reduces any squeeze-film damping effect), and e0 is the dielectric constant of the air. The boundary conditions for a clamped–clamped arch are: ^ wð0, t^ Þ ¼ 0,

^ @w ^ ¼ 0, ð0, tÞ @x^

^ ¼ 0, ^ tÞ wðL,

^ @w ^ ¼0 ðL, tÞ @x^

w¼

^ w , d

w0 ¼

^0 w , d

x¼

x^ , L

t¼

t^ T

ð3Þ

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where T is a time constant defined by T ¼ rAL4 =EI. Therefore, the non-dimensional equations of motion and associated boundary conditions can be written as " #"Z (   ) # 1 @4 w @2 w @w @2 w d2 w0 @w 2 @w dw0 ¼ a1 dx a2 Fe þ 2 þc þ þ2 4 2 2 @t @x @x dx dx @x @t @x 0

wð0,tÞ ¼ 0,

@w ð0,tÞ ¼ 0, @x

b0 ½1cosð2pxÞ, 2d

wð1,tÞ ¼ 0,

Fe ¼

@w ð1,tÞ ¼ 0 @x

½VDC þ VAC cosðOtÞ2 ð1 þw0 þ wÞ2

L b0

ˆx L

VDC d

ð2Þ

For convenience, we introduce the following non-dimensional variables:

ˆˆ ˆ (x,t) w

ˆ ˆ 0(x) w

ð1Þ

b

h

VAC

Fig. 1. (a) Schematic of an electrically actuated clamped–clamped arch and (b) a 3-D schematic picture of the arch.

ð4Þ ð5Þ

ð6Þ

ARTICLE IN PRESS 706

H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

Table 1 Summary of the non-dimensional variables used in the derivation. Parameter

Definition

a1 ¼ 6ðd=hÞ2 4 0 bL a2 ¼ e2EId 3 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi O ¼ O~ =on , where on ¼ EI=rAL4 ,

The stretching parameter

c ¼ c~ L4 =EIT

The damping parameter

The electric force parameter The excitation frequency parameter

To simulate the dynamic behavior, Eq. (10) can be integrated with time using Runge–Kutta technique. To simulate the static response, all time dependent terms in the differential equations, Eq. (10), are set equal to zero, then the modal coordinates ui(t) are replaced by unknown constant coefficients ai (i ¼1, 2, y, n). This results in a system of non-linear algebraic equations in terms of the coefficients ai. The system is then solved numerically using the Newton–Raphson method to obtain ai and hence the static deflection of the arch.

3. Reduced-order model To simulate the response of the shallow arch, Eqs. (4) and (5) are discretized using the Galerkin procedure to yield a reducedorder model (ROM) [43,44]. Hence, the deflection of the shallow arch is approximated as wðx,tÞ ¼

n X

ui ðtÞfi ðxÞ

ð7Þ

i¼1

where fi(x) (i¼1, 2, y, n) are the trial or admissible functions satisfying the boundary conditions, Eq. (5). Here, we investigate using either the normalized linear undamped mode shapes of the straight unactuated microbeam (w0 ¼0) or the exact mode shapes of the unactuated arch microbeam (w0≠0). The variables ui(t) (i¼1, 2, y, n) are the non-dimensional modal coordinates. The mode shapes of a clamped–clamped straight beam can be expressed as

fi ðxÞ ¼ cosh bi xcos bi x þ li ðsin bi xsinh bi xÞ

ð8Þ

where bi ¼ ðð2iþ 1ÞpÞ=2, and l1 ¼ 0:9825, li  1 for ig1:The linear undamped mode shapes of a clamped–clamped arch can be expressed as [45]

fi ðxÞ ¼ cosh bi xcos bi x þ li ðsin bi xsinh bi xÞ þ l cos 2px

ð9Þ

where the coefficients li and l are constants depending on the initial rise value b0, and b1 ¼6.54251, b2 ¼7.8532, b3 ¼11.3149, b4 ¼14.1372, and b5 ¼17.306. To obtain the ROM, we first multiply Eq. (4) by (1+w0 + w)2. This reduces the computational cost since the electrostatic force term in the discretized equation will not require complicated numerical integration (integrating a numerator term over a denominator term numerically is computationally expensive) [43]. Then, substituting Eq. (7) into the resulting equation, multiplying by fi(x), and then integrating the outcome from 0 to 1, yield the following differential equations in terms of the modal coordinates ui(t): n X

Mij u€ i ðtÞ þ

i¼1

n X

Cij u_ i ðtÞ þ

i¼1

n X

Kij ui ðtÞ ¼ Fj ðtÞ,

8j ¼ 1, . . . ,n

ð10Þ

i¼1

where Mij ¼

Z

1

Z

0 1

½fi ðxÞfj ðxÞf1 þ w0 ðxÞ þ wðxÞg2  dx,

8i,j ¼ 1, . . . ,n

Cij ¼ c ½fi ðxÞfj ðxÞf1 þw0 ðxÞ þ wðxÞg2  dx, 8i,j ¼ 1, . . . ,n 0 Z 1 iv ½fi ðxÞfj ðxÞf1 þ w0 ðxÞ þ wðxÞg2  dx, 8i,j ¼ 1, . . . ,n Kij ¼ 0 (Z Z 1 n 2 1 X 0 fj ðxÞ dxþ a1 G ui ðtÞf i ðxÞ dx Fj ðtÞ ¼ a2 V 2 0

Z

1

þ2

G¼

Z 1( 0

0

0

n X

dw0 ui ðtÞf i ðxÞ dx i¼1 n X

0

!i ¼ 1) dx , !

8j ¼ 1, . . . ,n

d2 w0 ð1 þw0 þ wÞ2 fj ðxÞ ui ðtÞf i ðxÞ dx2 i¼1

V ¼ VDC þ VAC cosðOtÞ

00

) dx ð11Þ

4. The static response As a case study, we consider the fabricated clamped–clamped shallow arch made of silicon of Krylov et al. [29] of L¼1000 mm, h¼2.4 mm, b¼ 30 mm, d¼10.1 mm, and initial rise b0 ¼3.5 mm. First, we investigate using the mode shapes of a straight beam, Eq. (8), as basis functions in the Galerkin procedure, Eq. (7). We examined the contribution from the anti-symmetric modes in the ROM and found it negligible for the considered arch (unlike the case cited in [46] for a deep arch loaded by a static point load). Hence, for the static analysis, we retain only the symmetric mode shapes. Next, we examine the convergence of the ROM. Fig. 2a shows the maximum static deflection of the shallow arch (wmax ¼ w(x¼0.5)) when using one up to six symmetric mode shapes while varying the DC load. Fig. 2 and all subsequent figures show the absolute value of the mud-deflection. This for convenience, since the actual deflection toward the lower electrode should have a negative sign. It follows from the figure that using five symmetric modes yields acceptable converged results. As seen in the figure, the shallow arch undergoes a snap-through motion near VDC ¼88 V and then a pull-in instability near VDC ¼106 V. Next, we investigate using the exact mode shapes of the shallow arch, Eq. (9), in the reduced-order model. Also here, only symmetric mode shapes are used. In Fig. 2b, we compare the obtained results using the first five exact mode shapes of the arch to those using the first five mode shapes of a straight beam. The figure shows excellent agreement between both approaches. It is worth to mention that the five exact mode shapes of the arch are the minimum needed for acceptable convergence, like the case of straight beam mode shapes. The conclusion here is that using the mode shapes of a straight beam is sufficient to yield accurate results for arches and this confirms the approach used in [29] to study the static behavior of MEMS arches under electrostatic loads. From a computational point of view however, straight beam mode shapes are much favorable. Therefore, in this paper, we will use only the mode shapes of a straight beam as basis functions. In the same figure, we validate the results using the five modes of the ROM for the clamped–clamped arch by comparing them with the experimental data of Krylov et al. [29]. It is clear that the experimental data and the ROM results are in good agreement. In the previous results, we neglected the influence of fringing of the electrostatic fields. This approximation is shown to be valid for devices of tiny aspect ratio (gap/length51) [47]. Note here that the initial gap and rise of the considered shallow arch are comparable with its width. Hence, the influence of the fringing fields might be important [48]. To investigate this, we use the Mejis–Fokkema formula [48] in Eq. (4) for the electrostatic force: (   2 VDC 1 þw0 þ w 3=4 1 þ 0:265 Fe ¼ ð1 þ w0 þ wÞ2 b^  1=2 ) h 1þ w0 þ w ð12Þ þ 0:53 h b^ where b^ ¼ b=d.

ARTICLE IN PRESS H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

10

12

Without fringing fields effect

9

With fringing fields effect

10

8 7

6

wmax (m)

8 wmax (m)

707

One mode Two modes

4

Three modes Five modes

1

0 20

4

2

Six modes

0

5

3

Four modes

2

6

40

60 VDC (Volt)

80

100

0

120

0

20

40

60

80

100

120

VDC (Volt)

14 Fig. 3. Variation of the static deflection of the shallow arch with the DC voltage with and without including the effect of the fringing field of the electrostatic force.

Five mode shapes of the straight beam

12

Five mode shapes of the shallow arch Experimental data [29]

wmax (m)

10

12

8

b0 = 0m

10

b0 = 2m

6

b0 = 3m

8 wmax (m)

4 2

b0 = 4m

6

4

0 0

20

40

60 VDC (Volt)

80

100

120

Fig. 2. (a) Variation of the static deflection of the shallow arch with the DC voltage for various number of symmetric mode shapes of a straight beam in the ROM. (b) Comparison between the obtained static deflection of the shallow arch with the DC voltage using five mode shapes of a straight beam and using the exact mode shapes of the arch with experimental data.

Fig. 3 compares between the static deflection of the shallow arch with and without considering the fringing effect. The figure shows negligible effect of the fringe fields of the electric load on the response of the microbeam. Hence, the initial assumption of neglecting the fringe effect is shown to be acceptable for this case. In Fig. 4, the effect of the initial rise on the static deflection of the arch is shown. For the cases of b0 ¼ 3 and 4 mm, the arch snapsthrough first and pulls-in while increasing VDC. However, for b0 ¼4 mm the arch undergoes immediate pull-in after snapthrough. The figure shows that the snap-through voltage increases and the pull-in voltage decreases when increasing the initial rise value of the shallow arch. This indicates that the stiffness of the shallow arch increases before snap-through and then decreases in the buckled position with the increase of b0.

5. Natural frequencies and mode shapes under a DC load In this section, we investigate the variation of the natural frequencies and mode shapes of the shallow arch with various rise levels under the actuation of the DC voltage. Toward this, we

2

0 0

20

40

60 VDC (Volt)

80

100

120

Fig. 4. Variation of the static deflection of the shallow arch with the DC voltage for various values of initial rise b0.

consider the ROM obtained in Eq. (10), with the mode shapes of a straight beam, which can be represented in a matrix form as MðuÞu_ ¼ RðuÞ

ð13Þ

where u ¼ ½u1 , u2 , u3 , . . . , un , u_ 1 , u_ 2 , u_ 3 , . . . , u_ n 

ð14Þ

is the modal amplitudes vector, M(u) is a non-linear matrix _ and R(u) is a right-hand side representing the coefficients of u, vector representing the forcing, stiffness, and damping coefficients. Both M(u) and R(u) are non-linear functions of the modal coordinates ui(t). Note here that we consider the symmetric and anti-symmetric mode shapes in the ROM to get all the possible natural frequencies and mode shapes of the shallow arch. Next, we split u into a static component Xs, representing the equilibrium position due to the DC actuation, and a dynamic component Z(t) representing the perturbation around the

ARTICLE IN PRESS 708

H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

equilibrium position, that is:

1.8

u ¼ Xs þ ZðtÞ

1.6

Then, substituting Eq. (15) into Eq. (13), using the Taylor series expansion assuming small Z, eliminating the quadratic terms, and using the fact that R(Xs) ¼0, we obtain the following equation:

1.4 1.2 1 (x)

ð15Þ

MðXs ÞZ_ ¼ JðXs ÞZ

1

where J(Xs) is the Jacobian matrix calculated at the equilibrium points [49]. To calculate the natural frequencies of the shallow arch for a given voltage, we substitute the stable static solution, Xs, into the matrix M  1J and then find its corresponding eigenvalues. The eigenvalues were calculated using Mathematica [50] by solving the below equation, which gives a characteristic algebraic equation for the eigenvalue l

b0 = 0m

0.8

b0 = 2m

0.6

b0 = 3m b0 = 4m

0.4

ð16Þ

0.2

detðM 1 ðXs ÞJðXs ÞlIÞ ¼ 0

0 0

2

4

6

8

10

8

10

where I is the identity matrix. Then by taking the square root of each individual eigenvalue, we obtain the natural frequencies of the system. Finally, by solving for the eigenfunctions for each eigenvalue, we get the mode shapes of the arch. Fig. 5a shows the effect of varying the initial rise of the arch investigated in Fig. 3 on its first mode shape. The mode shape is R1 normalized such that 0 F21 ðxÞ dx ¼ 1. As seen from the figure, as the arch initial rise increases, its first mode shape shows more deformation and deviation from that of a straight beam (b0 ¼0 mm), which is the first trial function used in the ROM (f1(x) in Eq. (8)). Fig. 5b shows the effect of the DC voltage on the first mode shape for values ranging from 40 to 100 V, where the arch undergoes snap-through. Fig. 6 shows the variation of the first five natural frequencies of the shallow arch for the case of b0 ¼3.5 mm. As seen in the figure, the fundamental natural frequency decreases to zero when snapthrough occurs. The reason for this is because of the softening effect of the quadratic geometric non-linearity of the arch in addition to the quadratic effect of the electrostatic force. Beyond snap-through, the cubic geometric non-linearity of the arch, which is of hardening effect, dominates both the geometric and electrostatic quadratic non-linearities for some range before it is dominated again by the electrostatic force near pull-in. Therefore,

x (µm) 1.5

1 (x)

1

VDC = 40 Volt VDC = 80 Volt

0.5

VDC = 100 Volt

0 0

2

4

6 x (µm)

Fig. 5. The first mode shape of the arch while varying (a) its initial rise b0 and (b) the DC electrostatic load VDC.

400 350 Nondimentional frequency

5

300 250 4

200

50 40

150

3

100

30 20

2

10

50

0

1

90

95

100

105

0 0

20

40

60 VDC (Volt)

80

100

120

Fig. 6. Variation of the first five natural frequencies with the DC voltage of a clamped–clamped shallow arch.

ARTICLE IN PRESS H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

where e is a bookkeeping parameter. We seek a solution in the following form:

50 45

wðx,t, eÞ ¼ ws ðxÞ þ uðx,tÞ

40 Nondimentional frequency

709

¼ ws ðxÞ þ eu1 ðx,T0 ,T2 Þ þ e2 u2 ðx,T0 ,T2 Þ þ e3 u3 ðx,T0 ,T2 Þ þ   

ð19Þ

35

where ws is the static component of the arch deflection and u is its dynamic component. For simplicity, we define Z 1 GðpðxÞ,gðxÞÞ ¼ ½p0 ðxÞg 0 ðxÞ dx ð20Þ

30 25 20

0

b0 = 0 m

15

Substituting Eqs. (17)–(20) into Eqs. (4) and (5) and equating like powers of e, we obtain

b0 = 2 m

10

b0 = 3 m

5

 Order e0: (the static equation)

b0 = 4 m

0

00

0

20

40

60 VDC (Volt)

80

100

w0 s ð0Þ ¼ w0 s ð1Þ ¼ 0

ws ð0Þ ¼ ws ð1Þ ¼ 0,

Fig. 7. Variation of the fundamental natural frequency with the DC voltage for a clamped–clamped shallow arch for various values of initial rise b0.

we observe from the zoomed view of Fig. 6 that the fundamental natural frequency increases just after snap-through while increasing the DC voltage. Then, near pull-in, it drops again to zero. Furthermore, the figure shows that the higher natural frequencies are sensitive to the variation of the DC voltage compared to those of the straight microbeams which have been shown to be insensitive to this variation [51]. Next, we show the effect of the initial rise on the fundamental natural frequency of the arch, Fig. 7. For all values of initial rise, the natural frequency drops to zero before snap-through. Beyond snap-through, some curves show increase of the natural frequency followed by a sudden drop to zero at pull-in. One curve however, b0 ¼ 4 mm, indicates that the arch undergoes immediate pull-in after snap-through if the DC load is increased continuously. In this case, the softening effect of the electrostatic force is always dominant. As seen in this figure, the snap-through voltage increases and the pull-in voltage decreases when increasing the value of the initial rise of the shallow arch.

2 a2 VDC

00

wiv s ¼ a1 ½w s þ w 0 ½Gðws ,ws Þ þ 2Gðws ,w0 Þ

120

ð1 þ w0 þws Þ2 ð21Þ

 Order e1: Lðu1 Þ ¼ 0

ð22Þ

where 00

00

Lðu1 Þ ¼ D20 u1 þ uiv 1 a1 ½Gðws ,ws Þ þ 2Gðws ,w0 Þu 1 2a1 ½w 2 2a2 VDC 00 þ w 0 ½Gðws ,u1 Þ þ Gðw0 ,u1 Þ u1 ð1 þ w0 þws Þ3

s

 Order e2: 00

00

Lðu2 Þ ¼ a1 Gðu1 ,u1 Þðw s þ w 0 Þ þ 2a1 ½Gðw0 ,u1 Þ þ Gðws ,u1 Þu 

2 2 VDC

3a

ð1 þ w0 þ ws Þ4

u21

00

1

ð23Þ

 Order e3: 00

00

Lðu3 Þ ¼ 2D0 D2 u1 cD0 u1 þ 2a1 Gðu1 ,u2 Þðw s þ w 0 Þ 00

þ 2a1 ½Gðws ,u1 Þ þ Gðw0 ,u1 Þu

2

00

þ 2a1 ½Gðws ,u2 Þ þ Gðw0 ,u2 Þu 1 þ a1 Gðu1 ,u1 Þu 

2a2 VDC VAC cosðOtÞ ð1þ w0 þ ws Þ2



2 6a2 VDC

ð1 þw0 þws Þ4

00

1 2 4a2 VDC

u1 u2 þ

ð1 þw0 þws Þ5

u31

ð24Þ 6. Response to small DC and AC loads In this section, we investigate in-depth the dynamic response of the clamped–clamped arch when subjected to small DC and AC loads as it is the case in MEMS resonators. Perturbation analysis is utilized. The motivation behind using this technique is that it represents a good analytical tool for the improvement and optimization of the arches to be used as MEMS resonators. Next, we present some details of the perturbation analysis.

Because in the absence of internal resonances and in the presence of damping all the indirectly exited modes die out, the solution of Eq. (22) is assumed to consist of only the directly excited mode, f(x). Accordingly, we express the dynamic component u1 as u1 ðx,T0 ,T2 Þ ¼ ½AðT2 ÞeioT0 þ AðT2 ÞeioT0 fðxÞ

ð25Þ

where A(T2) is a complex-valued function, the over bar denotes the complex conjugate, and o and f(x) are the natural frequency and corresponding eigenfunction of the directly excited mode, respectively. Substituting Eq. (25) into Eq. (23), we obtain

6.1. Perturbation analysis

2

Lðu2 Þ ¼ ð2AA þA2 e2ioT0 þ A e2ioT0 ÞhðxÞ Here, we apply the method of multiple scales and a direct attack of the equations of motion [11,52], Eqs. (4) and (5). To this end, we define the following variables for the time scale (Ti) and its derivative T0 ¼ t,

D0 ¼

@ , @T0

T1 ¼ et,

D1 ¼

@ , @T1

T2 ¼ e2 t,

D2 ¼

@ @T2

ð17Þ

Next, we scale the damping coefficient c and the forcing amplitude VAC so that the non-linearity balances their effects in the modulation equations (Eqs. (39) and (40)) [11,52], hence 2

c ¼ e c,

3

VAC ¼ e VAC

ð18Þ

ð26Þ

where 2 3a2 VDC

2

00

f þ 2a1 ðGðws , fÞ þ Gðw0 , fÞÞf ð1 þ w0 þws Þ4 00 00 þ a1 Gðf, fÞðw s þ w 0 Þ

hðxÞ ¼ 

ð27Þ

The solution of Eq. (26) can be expressed as follows: 2

u2 ðx,T0 ,T2 Þ ¼ c1 ðxÞA2 ðT2 Þe2ioT0 þ 2c2 ðxÞAðT2 ÞAðT2 Þ þ c1 ðxÞA ðT2 Þe2ioT0

ð28Þ

ARTICLE IN PRESS 710

H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

where c1 and c2 are the solutions of the following boundaryvalue problems: Hðci ,2od1i Þ ¼ hðxÞ,

i ¼ 1, 2

ci ¼ 0 and c0 i ¼ 0 at x ¼ 0 and x ¼ 1, i ¼ 1, 2

ð29Þ

Substituting Eqs. (25) and (28) into Eq. (19) and setting e ¼1, we obtain, to the second-order approximation, the following response of the MEMS arch to the AC harmonic excitation: wðx,tÞ ¼ ws ðxÞ þa cosðOtgÞfðxÞ þ 12a2 ½c1 ðxÞcos2ðOtgÞ þ c2 ðxÞ þ    ð41Þ

where dij is the Kronecker delta operator and the linear differential operator H is defined as 00

00

Hðf ðxÞ, oÞ ¼ f iv ðxÞo2 f ðxÞ2a1 ðw 0 þw s Þ½Gðf ðxÞ,w0 Þ þ Gðf ðxÞ,ws Þ 2 2a2 VDC 00 a1 ½Gðws ,ws Þ þ 2Gðw0 ,ws Þf ðxÞ f ðxÞ ð1 þ w0 þ ws Þ3 ð30Þ Note here that the eigenfunction f(x) is the solution of Hðf, oÞ ¼ 0

ð31Þ

In order to describe the nearness of the excitation frequency O to the fundamental natural frequencyo, we introduce a detuning parameter s defined by

O ¼ o þ e2 s

ð32Þ

Substituting Eqs. (25), (28), and (32) into Eq. (24) we obtain Lðu3 Þ ¼ ½ioð2A0 þ cAÞfðxÞ þ wðxÞA2 A þ F ðxÞeisT2 eioT0 þcc þ NST where 2a2 VDC VAC

ð34Þ

ð1þ w0 þ ws Þ2

In Eq. (33), A0 denotes the derivative of A with respect to T2, cc denotes the complex conjugate of the preceding terms, NST stands for the terms that do not produce secular terms, and w(x) is defined by 00

wðxÞ ¼ a1 f ½3Gðf, fÞ þ2Gðws , c1 Þ þ 2Gðw0 , c1 Þ þ 4Gðws , c2 Þ 00 00 þ 4Gðw0 , c2 Þ þ a1 ½w s þ w 0 ½2Gðf, c1 Þ þ 4Gðf, c2 Þ 00 00 þ 2a1 ½c 1 þ2c 2 ½Gðf,ws Þ þ Gðf,w0 Þ 2 2 2 12a2 VDC 6a2 VDC 12a2 VDC 3 þ

ð1 þ w0 þws Þ5

f 

ð1þ w0 þ ws Þ4

fc1 

ð1 þ w0 þ ws Þ4

fc2 ð35Þ

Multiplying the right-hand side of Eq. (33) by fðxÞeioT0 , where R f is normalized such that 10 f2 dx ¼ 1, integrating the result from x¼ 0 to 1, and equating the secular terms to zero [52] yield the following solvability condition: ioð2A0 þ cAÞ þSA2 A þ FeisT2 ¼ 0

ð36Þ

where Z 1 F¼ fF dx

ð37Þ

and

S¼

0

Z

1

fw dx 0

Next, we express A in the polar form A ¼ aeib =2, where a ¼a(T2) and b ¼ b(T2) are real-valued functions, representing, respectively, the amplitude and phase of the response. Substituting for A in Eq. (36) and letting g ¼ sT2  b, we obtain !    3 Z 1 Z 1 a 0 1 0 ib ioa þ oab  ioca e þ fw dx eib þ eisT2 fF dx ¼ 0 2 8 0 0 ð38Þ Separating the real and imaginary parts in Eq. (38), we obtain the following modulation equations: 1 sin g a0 ¼  ca þ F o 2 ag0 ¼ as þ

a3 cos g Sþ F o 8o

Solving Eq. (42) for the detuning parameters, we obtain vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi !ffi u u F2 a2 1 s ¼ t 2 2  c2  0 S ð43Þ 8o o a0 4 Recalling Eq. (32), setting e ¼1, and noting that the amplitude a is maximum when the expression inside the square root in Eq. (43) is zero, we obtain the following expression for the non-linear resonance frequency:

Or ¼ o þ ð33Þ

F ðxÞ ¼ 

It follows from Eq. (41) that periodic solutions correspond to constant a and g; that is, the fixed points (a0, g0) of Eqs. (39) and (40). Thus, letting g0 ¼ 0 and a0 ¼ 0 in Eqs. (39) and (40) and eliminating g0 yield the following frequency-response equation: " # 2 a20 F2 1 2 2 ð42Þ ¼ a s þ S þ c 0 8o 4 o2

ð39Þ

ð40Þ

S 2

8o5 z

F2

ð44Þ

where c¼2zo. We can see from Eq. (44) that the resonance frequency depends on the forcing coefficient F, the effective non-linearity coefficient S, the natural frequency o, and the damping ratio z. Next, we evaluate numerically the parameters o, f, c1, c2, and ws associated with Eq. (37) using the derived ROM in Section 3. Once these parameters are computed, frequency-response curves can be generated using Eq. (42). Note here that the stability of each calculated fixed point (a0, g0) is determined by examining the eigenvalues of the Jacobian matrix of Eqs. (39) and (40) evaluated at the corresponding fixed point. For a fixed point to be considered stable, all of the eigenvalues of the Jacobian matrix must lie in the left-half of the complex plane, eitherwise it will be marked as unstable. 6.2. Results To describe the dynamic response of the shallow arch, we need to determine the natural frequencyo, the excitation amplitude F, and the effective non-linearity of the system S. As a first step, Eq. (21) is numerically integrated using a five modes ROM, as done in Section 4, to determine the static deflection ws for a given DC voltage. Using the static solution ws, we solve the boundary-value problem, Eq. (31), for o and its corresponding eigenfunction f(x). Then, we solve the two boundary-value problems in Eq. (29) to evaluate the functions c1 and c2. Finally, we evaluate w, S, and F from Eqs. (35) and (37). Fig. 8a–c show different frequency-response curves for various AC and DC loads for the considered shallow arch of Fig. 3. It can be seen that the frequency-response curves bend more to the left while increasing the AC and DC values, which indicates softening behavior of the arch even for small motion. In the same figures, we compare the obtained frequency-response curves using the perturbation technique to those using a numerical time integration of the ROM differential equations using five mode shapes. We notice excellent agreement between both approaches for the cases of Fig. 8a and b. But, looking at Fig. 8c, one can see that the perturbation technique is limited to small DC and AC amplitudes. It cannot capture the snap-through motion of the arch.

ARTICLE IN PRESS H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

711

4.5 Perturbation (stable)

Perturbation

1.1

4

Long-time integration

1

Long-time integration

3.5

0.9 0.8

3 wmax (µm)

wmax (µm)

Perturbation (unstable)

0.7 0.6

2.5

0.5

2

0.4

1.5

0.3 1

0.2

0.5

0.1 25

30

35

40 

45

50

20

55

25

30

35

40

45

50

55



14 Perturbation (stable) Perturbation (unstable)

12

Long-time integration

wmax (µm)

10 8 6 4 2 0 20

25

30

35

40

45

50

55

 Fig. 8. Simulated frequency-response curves obtained using both the perturbation analysis and the long-time integration of the reduced-order model, Eq. (12). The assumed damping ratio for all cases is z ¼0.1: (a) VDC ¼20 V, VAC ¼20 V, (b) VDC ¼40 V, VAC ¼30 V, (c) VDC ¼40 V, and VAC ¼40 V.

One important advantage of the perturbation analysis is that it enables studying the variation of the effective non-linearity of the system analytically. In Fig. 9, we show the variation of the effective non-linearity coefficient S of the arch with the DC load. For the case of zero initial curvature (b0 ¼0 mm), S is positive for small values of DC voltages indicating a hardening behavior and then it switches to negative sign indicating softening-type behavior. On the other hand, the figure shows that for non-zero initial rise, S is always negative indicating a softening-type behavior of the shallow arch. These findings indicate that dynamically, the shallow arch is locally dominated by a softening-type behavior. This means that the quadratic nonlinearities coming from the initial curvature and the electrostatic force dominate the dynamic behavior of the arch. We notice also here that near the snap-through motion, the algorithm that we used to calculate the effective non-linearity coefficient experiences some numerical problems. For resonators application, it is important to predict accurately the non-linear resonance frequency of the arch. Fig. 10a and b

show the variation of the normalized non-linear resonance frequency Or =o0 for various values of VAC for two different DC loads and initial rises, respectively, as calculated from Eq. (44). We can see that the resonance frequency decreases with the AC load of both cases showing a softening-type behavior for the clamped– clamped arch confirming the results of Fig. 9. We notice also from those figures that the DC load enhances the softening behavior of the considered shallow arch, Fig. 10a, whereas the initial rise weakens it, Fig. 10b. This can be understood from Eq. (44), in which the resonance frequency is proportional to the forcing term F and hence to the AC and DC loads. Also, the resonance frequency is inversely proportional to the natural frequency, which is higher for high initial rise. Finally, we end this section by showing the effect of the damping ratio on the frequency-response curve. Fig. 11 shows that decreasing the damping ratio amplifies the effect of the nonlinearity (softening behavior) and can change the dynamic state of the arch from having a unique stable state (z ¼ 0.01) to multivalued states with possibilities of jumps and hysteresis (z ¼0.003

ARTICLE IN PRESS 712

H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

and 0.005). Again, this is expected from Eq. (44), in which the resonance frequency is inversely proportional to the damping ratio z.

7. Conclusions In this paper, an investigation into the non-linear dynamic behavior of an electrically actuated clamped–clamped shallow arch when actuated by a DC force and an AC harmonic load was presented. A perturbation method was used to study the dynamic behavior of the shallow arch for small DC and AC loads. Analytical expressions were derived, which represent good tools for engineers for the improvement and optimization of the MEMS shallow arch design. Dynamic analysis was conducted to explore the vibration of the arch near its primary resonance. The results indicate that for certain DC and AC loads, the arch resonator can

experience a dynamic snap-through motion for some band of frequency near primary resonance. MEMS designers need to be aware of these bands to ensure a reliable and stable operation of these arches as resonators. Using the perturbation technique, the non-linear resonance frequency of the arch as a function of the AC harmonic amplitude, the DC load, and the arch initial rise was calculated for small values of AC and DC loads. The results indicate that the resonance frequency is significantly lowered as the values of DC and AC loads are increased and for smaller values of initial rise. This indicates that, to be used as linear resonators, shallow arches should be designed with an optimized initial rise at low DC and AC loads. Further, the obtained results show that the dynamic response of the arch is of softening type even for small DC and AC values, as revealed through the frequency-response curves, the effective non-linearity, and the variation of the resonance frequency with the DC and AC loads.

x 105

3

1

Stable branch

0

2.5

-1

2 wmax (µm)

Effective nonlinearity - S

 = 0.003

-2 b0 = 0 m

-3

Unstable branch

 = 0.005

1.5

b0 = 2 m

1

b0 = 3 m

-4

 = 0.01

b0 = 4 m

0.5 -5

0

-6 0

20

40

60 VDC (Volt)

80

100

Fig. 9. Variation of the effective non-linearity coefficient S with the DC voltage.

1

1

0.95

0.95

0.9

0.9

0.85

0.85

0.8 VDC = 10 Volt VDC = 20 Volt

0.75

0.8

0.7

0.65

0.65

10

20

30

40 50 VAC (Volt)

60

70

80

36

37

38 

39

40

41

42

b0 = 4 m b0 = 3 m

0.75

0.7

0

35

Fig. 11. Simulated frequency-response curves obtained using the perturbation analysis for various damping ratios and for VDC ¼ 5 V, VAC ¼ 5 V, and b0 ¼ 3 mm.

r / ω0

r / ω0

34

120

0

10

20

30

40 50 VAC (Volt)

60

70

80

Fig. 10. The non-linear resonance frequency Or normalized to the linear natural frequency of the system o0 versus VAC for (a) various DC voltages (b0 ¼ 3 mm) and (b) various initial rises b0 (VDC ¼ 20 V).

ARTICLE IN PRESS H.M. Ouakad, M.I. Younis / International Journal of Non-Linear Mechanics 45 (2010) 704–713

Acknowledgment This work has been supported by the National Science Foundation through CAREER Award # 0846775. References [1] T.P. Burg, A.R. Mirza, N. Milovic, C.H. Tsau, G.A. Popescu, J.S. Foster, S.R. Manalis, Vacuum-packaged suspended microchannel resonant mass sensor for biomolecular detection, J. Microelectromech. Syst. 15 (2006) 1466–1476. [2] H.C. Kim, S. Seok, I. Kim, S.-D. Choi, K. Chun, Inertial-grade out-of-plane and inplane differential resonant silicon accelerometers (DRXLs), in: Transducers 05, Solid-State Sensors, Actuators and Microsystems, Seoul, South Korea, 2005. [3] W.T. Hsu, J.R. Clark, C.T.C. Nguyen, A resonant temperature sensor based on electrical spring softening, in: Transducers 01, Eurosensors XV, The 11th International Conference on Solid-State Sensors and Actuators, Munich, Germany, 2001. [4] B. Kim, M. Hopcroft, C.M. Jha, R. Melamud, S.A. Chandorkar, M. Agarwal, K.L. Chen, W.T. Park, R. Candler, G. Yama, A. Partridge, M. Lutz, T.W. Kenny, Using MEMS to build the device and the package, in: Transducers 07 Eurosensors, Lyon, France, 2007. [5] B. Kim, M.A. Hopcroft, R. Melamud, C.M. Jha, M. Agarwal, S.A. Chandorkar, T.W. Kenny, CMOS Compatible Wafer-scale Encapsulation with MEMS Resonators, ASME InterPACK, Vancouver, Canada, 2007. [6] C.M. Jha, M.A. Hopcroft, S.A. Chandorkar, J.C. Salvia, M. Agarwal, R.N. Candler, R. Melamud, B. Kim, T.W. Kenny, Thermal isolation of encapsulated MEMS resonators, J. Microelectromech. Syst. 17 (2008) 175–184. [7] H.C. Nathanson, R.A. Wickstrom, A resonant gate silicon surface transistor with high Q bandpass properties, IEEE Appl. Phys. Lett. 7 (1965) 84–86. [8] H.C. Nathanson, W.E. Newell, R.A. Wickstrom, J.R. Davis, The resonant gate transistor, IEEE Trans. Electr. Dev. 14 (1967) 117–133. ¨ P. Rantakari, [9] T. Mattila, O. Jaakkola, J. Kiihamaki, J. Karttunen, T. Lamminmaki, ¨ H. Kattelus, I. Tittonen, 14 MHz micromechanical oscillator, A. Oja, H. Seppa, Sensors Actuators A 97–98 (2002) 497–502. [10] R.M.C. Mestrom, R.H.B. Fey, J.T.M. van Beek, K.L. Phan, H. Nijmeijer, Modelling the dynamics of a MEMS resonator: simulations and experiments, Sensors Actuators A 142 (2008) 306–315. [11] M.I. Younis, A.H. Nayfeh, A study of the nonlinear response of a resonant microbeam to an electric actuation, J. Nonlinear Dyn. 31 (2003) 91–117. [12] M.T.A. Saif, On a tunable bistable MEMS—theory and experiment, J. Microelectromech. Syst. 9 (2000) 157–170. [13] N.D. Masters, L.L.A. Howell, Self-retracting fully compliant bistable micromechanism, J. Microelectromech. Syst. 12 (2003) 73–80. [14] J. Rossiter, B. Stoimenov, T. Mukai, A self-switching bistable artificial muscle actuator, in: Proceedings of the SICE-ICASE International Joint Conference, Busan, Korea, 2006. [15] R.A.M. Receveur, C.R. Marxer, R. Woering, V.C.M.H. Larik, N.F. de Rooij, Laterally moving bistable MEMS DC switch for biomedical applications, J. Microelectromech. Syst. 14 (2005) 89–98. [16] B. Charlot, W. Sun, K. Yamashita, H. Fujita, H. Toshiyoshi, Bistable nanowire for micromechanical memory, J. Micromech. Microeng. 18 (2008) 045005–045012. [17] M. Sulfridge, T. Saif, N. Miller, M. Meinhart, Nonlinear dynamic study of a bistable MEMS: model and experiment, J. Microelectromech. Syst. 13 (2004) 725–731. [18] M. Vangbo, An analytical analysis of a compressed bi-stable buckled beam, Sensors Actuators A 69 (1998) 212–216. [19] Y.J. Yang, C.J. Kim, Testing and characterization of a bi-stable snapping microactuator based on thermo-mechanical analysis, in: Proceedings of the IEEE Solid State Sensor and Actuator Workshop, Stockholm, Sweden, 1995. [20] B. Wagner, H.J. Quenzer, S. Hoershelmann, T. Lisec, M. Juerss, Bi-stable microvalves with pneumatically coupled membranes, in: Proceedings of the IEEE MEMS’96, San Diego, CA, 1996. [21] J. Qui, J.H. Lang, A.H.A. Slocum, Curved beam bistable mechanism, J. Microelectromech. Syst. 13 (2004) 137–146. [22] J.S. Han, J.S. Ko, Y.T. Kim, B.M. Kwak, Parametric study and optimization of a micro-optical switch with a laterally driven electromagnetic microactuator, J. Micromech. Microeng. 12 (2002) 39–47.

713

[23] P. Seunghoon, H. Dooyoung, Pre-shaped buckled beam actuators: theory and experiment, Sensors Actuators A: Phys. 148 (2008) 186–192. [24] J.S. Ko, M.G. Lee, J.S. Han, J.S. Go, B. Shin, D.S. Lee, A laterally-driven bistable electromagnetic microrelay, ETRI J. 28 (2006) 389–392. [25] A. Michael, C.Y. Kwok, Design criteria for bi-stable behavior in a buckled multi-layered MEMS bridge, J. Micromech. Microeng. 16 (2006) 34–43. [26] J. Qui, J.H. Lang, A.H. Slocum, A.C. Weber, A bulk-micromachined bistable relay with U-shaped thermal actuators, J. Microelectromech. Syst. 14 (2005) 99–109. [27] J. Casals-Terre´, A.M. Shkel, Snap-action bistable micromechanism actuated by nonlinear resonance, in: 2005 IEEE Sensors, Irvine, CA, USA, 2005. [28] Y. Zhang, Y. Wang, Z. Li, Y. Huang, D. Li, Snap-through and pull-in instabilities of an arch-shaped beam under an electrostatic loading, J. Microelectromech. Syst. 16 (2007) 684–693. [29] S. Krylov, B.R. Ilic, D. Schreiber, S. Seretensky, H. Craighead, The pull-in behavior of electrostatically actuated bistable microstructures, J. Micromech. Microeng. 18 (2008) 055026–055046. [30] K. Das, R.C. Batra, Pull-in and snap-through instabilities in transient deformations of microelectromechanical systems, J. Micromech. Microeng. 19 (2009) 035008–035027. [31] L. Buchaillot, O. Millet, E. Que´vy, D. Collard, Post-buckling dynamic behavior of self-assembled 3D microstructures, Microsyst. Technol. 14 (2007) 69–78. [32] A. Cabal, D.S. Ross, Snap-through bilayer microbeam, Nanotechnology 1 (2002) 230–233. [33] W.Y. Poon, C.F. Ng, Y.Y. Lee, Dynamic stability of a curved beam under sinusoidal loading, Proc. Inst. Mech. Eng., Part G—J. Aero. Eng. 216 (2002) 209–217. [34] S.P.M. Noijen, N.J. Mallon, R.H.B. Fey, H. Nijmeijer, G.Q. Zhang, Periodic excitation of a buckled beam using a higher order semianalytic approach, Nonlinear Dyn. 50 (2007) 325–339. [35] J. Humphreys, On dynamic snap buckling of shallow arches, AIAA J. 4 (5) (1966) 878–886. [36] C.S. Hsu, C.T. Kuo, R.H. Plaut, Dynamic stability criteria for clamped shallow arches under timewise step loads, AIAA J. 7 (1969) 1925–1931. [37] N.C. Hung, Dynamic buckling of some elastic shallow structures subjected to periodic loading with high frequency, Int. J. Solids Struct. 8 (1972) 315–326. [38] Q. Bi, H.H. Dai, Analysis of non-linear dynamics and bifurcations of a shallow arch subjected to periodic excitation with internal resonance, J. Sound Vib. 233 (2000) 557–571. [39] W.M. Tien, N.S. Namachchivaya, A.K. Bajaj, Non-linear dynamics of a shallow arch under periodic excitation-I. 1-2 internal resonance, Int. J. Non-Linear Mech. 29 (1994) 349–366. [40] W.M. Tien, N.S. Namachchivaya, N. Malhorta, Non-linear dynamics of a shallow arch under periodic excitation-II. 1-1 internal resonance, Int. J. NonLinear Mech. 29 (1994) 367–386. [41] F. El-Bassiouny, Three-to-one internal resonance in the nonlinear oscillation of shallow arch, Phys. Scr. 72 (2005) 439–450. [42] A.H. Nayfeh, M.I. Younis, E.M. Abdel-Rahman, Reduced-order models for MEMS applications, Nonlinear Dyn. 41 (2005) 211–236. [43] M.I. Younis, E.M. Abdel-Rahman, A.H. Nayfeh, A reduced-order model for electrically actuated microbeam-based MEMS, J. Microelectromech. Syst. 12 (2003) 672–680. [44] A.H. Nayfeh, Nonlinear Interactions, Wiley Interscience, New York, 2000. [45] W. Lacarbonara, G. Rega, Resonant non-linear normal modes. Part II: activation/orthogonality conditions for shallow structural systems, Int. J. Non-Linear Mech. 38 (2003) 873–887. [46] A.B. Pippard, The elastic arch and its modes of instability, Eur. J. Phys. 11 (1990) 359–365. [47] R.C. Batra, M. Porfiri, D. Spinello, Typical review: review of modeling electrostatically actuated microelectromechanical systems, J. Micromech. Microeng. 16 (2007) R23–R31. [48] R.C. Batra, M. Porfiri, D. Spinello, Capacitance estimate for electrostatically actuated narrow microbeams, Micro Nano Lett. 1 (2006) 71–73. [49] A.H. Nayfeh, B. Balachandran, Applied Nonlinear Dynamics, Wiley, New York, 1995. [50] Wolfram Research, Inc., Mathematica, Version 7.0, Champaign, IL, 2008. [51] E.M. Abdel-Rahman, M.I. Younis, A.H. Nayfeh, Characterization of the mechanical behavior of an electrically actuated microbeam, J. Micromech. Microeng. 12 (2002) 759–766. [52] A.H. Nayfeh, Introduction to Perturbation Techniques, Wiley Interscience, New York, 1981.