Gap Dependent Bifurcation Behavior of a Nano-Beam Subjected to a Nonlinear Electrostatic Pressure

This paper presents a study on the gap dependent bifurcation behavior of an electro statically-actuated nano-beam. The sizedependent behavior of the beam was taken into account by applying the couple stress theory. Two small and large gap distance regimes have been considered in which the intermolecular vdW and Casimir forces are dominant, respectively. It has been shown that changing the gap size can affect the fundamental frequency of the beam. The bifurcation diagrams for small gap distance revealed that by changing the gap size, the number and type of the fixed points can change. However, for large gap regime, where the Casimir force is the dominant intermolecular force, changing the gap size does not affect the quality of the bifurcation behavior.


INTRODUCTION
Nowadays, with the rapid development in nano technology, the nano electro-mechanical systems (NEMS) have become one of the hottest research topics.High speed, accuracy and performance and low energy consumption have increased the possibility of substituting the nano technology with micro technology.
NEMS devices such as random access memory, nano tweezers, super sensitive sensors, resonators, etc are widely designed, analyzed, fabricated and used (Adu et al., 2001;Collins et al., 2000;Kim and Lieber, 1999;Mobki et al., 2013;Rueckes et al., 2000;Tahami et al., 2009).There are several actuation mechanisms used in micro and nano electromechanical systems (MEMS & NEMS) such as electrostatic, piezoelectric, thermal, etc (Rezazadeh et al., 2011).Electrostatically actuated devices form a broad class of MEMS& NEMS devices due to their simplicity, as they require few mechanical components and small voltage levels for actuation (Senturia,   Mohammad Fathalilou a  Morteza Sadeghi b Ghader Rezazadeh c * Latin American Journal of Solids and Structures 11 (2014) 2426-2443   2011).In such circumstances, a conductive flexible beam or plate is suspended over a ground plate and a potential difference is applied between them.As the microstructure is balanced between electrostatic attractive force and mechanical (elastic) restoring force, both electrostatic and elastic restoring forces are increased by the upsurge of the electrostatic voltage.When the voltage reaches the critical value, pull-in instability will occur.Pull-in is the point at which the elastic restoring force can no longer balance the electrostatic force.Further increasing the voltage, gives rise to the dramatic displacement jump in the structure, causing structure collapse and fail.Pull-in instability is a snap-through like behavior and a saddle-node bifurcation type of instability (Zhang et al., 2006).
In NEMS devices, by decreasing the geometric dimension, intermolecular surface forces such as van der Waals force and Casimir force can affect the structure behavior.vdW force and Casimir force can both be connected with the existence of zero-point vacuum oscillations of the electromagnetic field (Batra et al., 2007;Bordag et al., 2001;Klimchitskaya, 2000;Lifshitz, 1956;Moeenfard et al., 2012).
The microscopic approach to the modeling of both vdW and Casimir forces can be formulated in a unified way using Quantum Field Theory (Bordag et al., 2001;Klimchitskaya, 2000;Lamoreaux, 2005;Lifshitz, 1956).Batra et al. (2008) have reported that the Casimir force is generally effective at larger separation distances between the bodies than the vdW force.Whereas the Casimir force between semi-infinite parallel plates is inversely proportional to the fourth power of the gap, vdW force is inversely proportional to the third power of the gap.The dependence of these forces on the dielectric properties of the plates and the filling medium is studied in detail by Bordag et al. (2001).It is important to note that vdW and Casimir forces cannot in general be considered to simultaneously act in MEMS, since they describe the same physical phenomenon at two different length scales (Batra et al., 2008).
Although these forces have been explored for decades but they are rarely considered in analyzing the mechanical behavior and calculating the Pull-in voltage of the NEMS systems (Tahami et al., 2009).Lin et al. (2005) have studied the pull-in phenomena and calculated the pull-in voltage for a nano-electromechanical switch with the assumption of one degree of freedom for a nano switch without taking into account the effect of vdW force.Dequesnes et al. (2002) have calculated the pull-in voltage for carbon-nanotube-based nano-electromechanical switchs.They have considered the vdW force in their model and introduced an analytical formula for pull-in voltage of a lumped model of the switch.They haven't considered the Casimir force in their model.Moeenfard et al. (2012) have investigated the static behavior of nano and micro-mirrors under the effects of Casimir force.Palasantzas and DeHosson (2005) have explored the influence of self-affine roughness on the pull-in parameters for nanoelectromechanical switches in the presence of the Casimir force.Tahami et al. (2009) studied static and dynamic pull-in phenomena of a capacitive nano-switch considering both Casimir and vdW effects.They have considered both vdW and Casimir forces in all gap distances between two parallel plates.But, as said before, many researchers such as Batra et al. (2007), Palasantzas et al . (2008), Boström et al. (2012) and Lambrecht and Reynaud (2000) have clarified that at smaller distances dominant force is the vdW force whereas for larger distances the Casimir force is dominant and there is a regime in which the transition from vdW to Casimir force Latin American Journal of Solids and Structures 11 (2014) 2426-2443 takes place.Palasantzas et al. (2008) have shown that the transition from vdW to Casimir regime is rather weak.They have obtained for two parallel gold surfacesabout 18 nm gap distance for the crossover from vdW to Casimir regime, whereas Lambrecht and Reynaud (2000) have predicted theoretically about 13 nm for the transition.
On the other hand, many researches showed that in micro and nano scales the materials have strong size dependence in deformation behavior (Papargyri-Beskou et al., 2003;Lazopoulos and Lazopoulos, 2010;Asghari et al., 2010;Park and Gao, 2006;Fu and Zhang, 2011;Kong et al., 2008Kong et al., , 2009;;Fathalilou et al., 2011).Size-dependent behavior is an inherent property of materials, which appears for a beam when the characteristic size such as thickness or diameter is close to the internal length-scale parameter of materials (Fathalilou et al., 2011).However, in some materials, remarkable discrepancies are observed between the experimental results and those obtained using classical elasticity theory.These discrepancies are mainly due to the dominance of atomic structures of the material neglected in classical elasticity (singh, 2008).So, the static or dynamic behavior of these structures cannot be correctly described by the classical linear elasticity and higher order mechanical theories must be applied (Papargyri-Beskou and Beskos, 2008).
Voigt was the first who tried to correct these shortcomings of the classical elasticity by taking into account of the assumption that interaction between the two parts through an area element inside the body is transmitted not only by a force vector but also by a moment vector giving rise to a 'couple stress theory' (Voigt, 1887).The complete theory of asymmetric elasticity was developed by Cosserat and Cosserat (1909), which was non-linear from the beginning.They assumed that each material point of a three dimensional continuum is associated with a 'rigid triad' and during the process of deformation; it can rotate independently, in addition to the displacement (singh, 2008).After a gap of about fifty years, Cosserats theory drew attention of researchers and several Cosserat-type theories were developed independently, e.g., Gunther (1958), Grioli (1960), Ra-jagopal (1960), Palmov (1964), Aero and Kuvshinskii(1960), Mindlin and Tiersten (1962), Toupin (1962), Eringen (1962), Koiter (1964), Nowacki (1974) among several others.Later, the general Cosserat continuum theory acquired the name of 'micropolar continuum theory' following Eringen (1966), in which the micro-rotation vector is taken independent of the displacement vector.Eringen and Suhubi (1964) and Suhubi and Eringen (1964) developed a non-linear theory for 'micro-elasticity', in which intrinsic motions of the microelements were taken into account.The general theory of Mindlin (1964) includes three equivalent forms which are defined on the basis of three different expressions for the strain energy density (Papargyri-Beskou and Beskos, 2008).The first expression involves gradients of displacements, the second one the gradients of strain, and the third one the gradients of rotation.The couple stress theory is based on the third expression of the strain energy density while second form gives rise to the gradient elastic theory.
Recently, some researchers have applied the non-classic theories of elasticity to study the mechanical behavior of electrostatically actuated micro and nano beams Asghari et al., 2010;Fu and Zhang, 2011;Kong et al., 2008Kong et al., , 2009;;Fathalilou et al., 2011).
In spite of the mentioned works about the electrostatically actuated nano-beams, there is no comprehensive study about their stability from bifurcation view point.In this paper, distributed as well as the lumped models of the nano-beam with electrostatic actuation are introduced, considering the couple stress theory of elasticity.The vdW and Casimir forces are taken into account for small Latin American Journal of Solids and Structures 11 (2014) 2426-2443 and large gap distances, respectively.It is obtained that the fundamental frequency of the nanobeam can be gap dependent.Also, it is understood from bifurcation analysis that changing the initial gap size can changes the number and type of bifurcation points.

Distributed model
Figure 1 shows an electrostatically actuated fixed-fixed Euler-Bernoulli nano-beam.The device consists of a beam, suspended over a dielectric film deposited on top of the center conductor and fixed both ends to the ground conductor.When a direct voltage (DC voltage) is applied between the upper and lower electrodes, the upper deformable beam is pulled down due to the electrical force.The rectangular cross-section beam is considered with length L, thickness h, and width b, and initial gap of 0 g .The nano-beam is assumed to be isotropic with Young modulus E and density  .
Together with the attractive electrostatic force, the vdW and Casimir forces pull the nano-beam down towards the substrate.Minimum length of the nano-beam in which vdW and Casimir forces lead the nano-beam to collapse on the substrate (in lack of the electrostatic force) is known as the detachment length (Tahami et al., 2009).
In the linear couple stress theory the strain energy of the deformed body is assumed to depend upon the strain  and the rotation gradient  , so that the associated stress quantities are the sym- metric Cauchy stress tensor  and the deviatoric couple stress tensor  (Mindlin, 1964).It follows that the strain energy s E in a deformed isotropic linear elastic material occupying region V will be: Also, the kinetic energy of the body motion is as following: The constitutive equations of the theory are given considering the following expression for the strain energy density s e (Mindlin, 1964): where  and  are two Lame's constants of the classical elasticity whereas  and are two non- classic Lame-type material constants which introduce the couple stress effects.Eq. ( 3) leads to the following constitutive equations: It is time now to introduce the governing equation for the beam of Figure 1 considering couple stress effects.Primarily, it is mentioned that i, j and k indices vary from 1 to 3; representing the variables in x, y and z directions in Cartesian coordinates, respectively.Using the coordinate system (x,z) shown in Figure 1, where x-axis coincides with the centroidal axis of the un-deformed beam and z-axis is the symmetry axis, the displacement components of an Euler-Bernoulli beam, neglecting the mid-point displacement in x direction can be represented by (Park and Gao, 2006): where u, v and w are, respectively, the (x, y) and z components of the displacement vector.Considering Eq. ( 5), the components of the symmetric strain tensor for the plane stress and plane strain conditions are as following, respectively: Plane stress: The components of the rotation vector can be written as: , 0 Now, the components of the asymmetric rotation-gradient tensor are as following: Substituting Eqs. ( 6) and ( 8), into Eq.( 4) the following relation is obtained for the strain energy density of the plane strain condition: where  is the shear modulus and ( ) ( ) , is the elasticity modulus of the material in plane strain condition.As seen in this relation, the strain energy in beam model does not depend on , and one non-classical material constant only appears in this model which is defined as: where l is the material length scale parameter (Anthoine, 2000;Mindlin, 1964).Considering mentioned relations and applying the Hamilton principle the governing equation with couple stress is obtained as following: with the following boundary conditions: Clearly, when the couple stress effect is suppressed by letting l=0, the present model defined by Eq. ( 11) will reduce to the classical Euler-Bernoulli beam model.In Eq. ( 11), ext q is considered as the sum of the electrostatic, vdW and Casimir forces as following (Tahami et al., 2009): where ( ) bV Ab b hc q q q g w g w g w In these equations,  et al., 2009).For convenience, Eq. ( 11) may be cast into a non-dimensional form; in particular, both the transverse displacement, wand the spatial coordinate, x are normalized by characteristic lengths of the system and the gap size, and beam length, respectively, according to: is the classic characteristic fundamental period of the system.Substituting these parameters into Eq.( 11), the following non-dimensional equation is obtained: Latin American Journal of Solids and Structures 11 (2014) 2426-2443 The parameters α, β, γ and η appeared in Eq. ( 15) are: (1 Due to the nonlinearity of the derived static equation, the solution is complicated and time consuming.Direct application of either Galerkin or finite difference methods will produce a set of nonlinear algebraic equation.To avoid this problem, we adopted the step-by-step linearization method (SSLM) (Rezazadeh et al., 2011), followed by Galerkin method to solve the obtained linear set of algebraic equations.For using SSLM, it is assumed that ˆk w is the displacement of the beam due to the applied voltage k V .The vdW and Casimir forces only depend on the non-dimensional gap size and their value cannot be controlled manually.In order to solve this problem, it is assumed that these forces being applied gradually in a virtual manner.For this purpose, a virtual variablel is introduced which varies from zero to one.Multiplying this variable by vdW and Casimir forces, the assumption of the step by step change of the load can be satisfied (Zarei and Rezazadeh, 2008).

Now, ˆk
w is defined as the non-dimensional displacement of the nanostructure when subjected to k V under the virtual nano scaled vdW and Casimir forces.Then by increasing the voltage and consequent virtual force variable l , the deflection of the (k+1) th step can be obtained as: where V d and dl are the voltage and virtual force variable, changing between two successive steps, respectively.By considering a small value of V d and dl the y will be small enough so that we can use the first term of the Taylor expansion in each step instead of exact main excitation function.The equation for k th step is: and for (k+1) th step: Subtracting Eq. ( 19) from ( 20) and considering Eq. ( 18), we have: Latin American Journal of Solids and Structures 11 (2014) 2426-2443 The obtained linear differential equation is solved by Galerkin method.ˆ ( ) x y having basis in function space, may be expressed as: The ( ) j x f is jth shape function which satisfies the given boundary conditions.In this paper, the shape functions are selected as the linear undamped natural mode shapes of the beam.The unknown ˆ ( ) x y is approximated by truncating the summation series to a finite number, n: Substituting Eq. ( 23) into ( 21), multiplying by as a weight function in Galerkin method and integrating the outcome from ˆ0 to 1 x = a set of linear algebraic equation is generated as: where Dynamic loading response may be obtained using Galerkin-based reduced order model (Nayfeh and Mook, 1979).Due to the non-linearity of the external forces, direct application of the reduced order model to the dynamic equation (Eq.( 15)) leads to generation of n nonlinear-coupled ordinary differential equation.To solve this difficulty, the forcing terms in Eq. ( 15) is considered as a constant in each integration time step, which takes the value of the previous step.Selecting time steps small enough, leads to satisfactory results.To achieve a reduced order model, ( ) ˆˆ, ŵ x t may be approximated as: By substituting Eq. ( 26) into Eq.( 15) and multiplying by ( ) i x f as a weight function in Galerkin method and integrating the outcome from ˆ0 to 1 x  , the Galerkin based reduced order model is generated as: Latin American Journal of Solids and Structures 11 (2014) 2426-2443   1 where M and K are mass and mechanical stiffness matrices, respectively.Also F introduces the forcing vector.The mentioned matrices and vector are given by: Now, Eq. ( 27) can be integrated over time by various integration methods such as Runge-Kutta method where ext q in each time step of integration takes the value of previous step.

Lumped model
A lumped model can be helpful to gain a rough quantitative estimation for the response of a wide range of electrostatically-actuated micro and nano-structures (Rezazadeh et al., 2011).The lumped model shown in Figure 2 is utilized to represent a NEMS device employing electrostatic actuation.
The device has a movable nano-beam of mass m, which forms one side of a variable capacitor.A spring with coefficient k is used to model the effective stiffness of the nano-beam, which is due to the elastic restoring force (Rezazadeh et al., 2011).Considering the couple stress theory, the equivalent mass (m) and stiffness (k) for a lumped model of the fixed-fixed nano-beam are considered as: deh et al., 2011).
Using the equivalent mass and stiffness, the equation of motion for lumped mass-spring model of the nano-beam can be written as following: Latin American Journal of Solids and Structures 11 (2014) 2426-2443 where z introduces displacement of the mass defined to be positive downward.
For adjusting the lumped model with distributed one, a model corrective coefficient 1.11 a = can be multiplied to k to obtain a modified stiffness (Rezazadeh et al., 2011).For convenience in analysis, Eq. ( 29) can take a nondimensional form with nondimensional variables as following: where * t is the characteristic time and equal to / m k .
Using the nondimensional parameters, the governing equation is rewritten as following: (1 ) where 1 2 , D D and 3 D are nondimensional parameters and defined as following: The static equation of deflection can be obtained by dropping time from Eq. ( 31) as following: (1 ) Using this equation, one can obtain the fixed points of the system by changing voltage as a parameter control.For analyzing the stability of the fixed points, a phase diagram for each voltage can be plotted.To this end, Eq. ( 31) can be solved by various methods such as Runge-Kutta method.

RESULTS AND DISCUSSION
In order to validate our results, the obtained pull-in voltages are compared to those obtained by Tahami et al. (2009) for a fixed-fixed silicon nano-beam with 0 10 and 6 g n m h n m = = .They have considered classic elasticity theory to obtain their results.They have obtained 6.78 V for pull-in voltage of this beam by considering the vdW force when the length is taken as , they obtained 4.11 V. Our results are 6.77V and 4.11 V for these two cases, respectively.As seen, there is a good agreement among the results.
For analyzing the bifurcation behavior, a capacitive gold nano-beam is considered with geometries and material properties as shown in Table 1.The length scale parameter is taken as 0.5 h .It is considered two gap distance regimes which in small distances, 0 10 g n m < , the only intermolecular force is vdW force whereas in large distances, 0 15 g n m > , the Casimir force is the only dominant intermolecular force (Palasantzas et al., 2008).As mentioned the vdW and Casimir forces pull down the movable beam.So it is expected that considering these forces changes the fundamental frequency of the beam.Figures 3a and b show the gap dependent fundamental frequency (  ) of the nano-beam without an electrostatic actuation for small and large gap regimes, respectively.As illustrated, for small gaps the rate of increasing the fundamental frequency with increasing the gap is higher than large distances.Also, for the beam with mentioned specifications, after the gap of 30 nm the frequency converges to a constant value.Figure 4 shows the bifurcation behavior of the nano-beam in state-control space for 0 6 g nm = .
It can be seen that for applied voltages smaller than pull-in voltage, there are two fixed points.For voltages between 0.7 1.7 V -, no fixed points are observed in this space.For voltages higher than 1.7 V two fixed points are observed again, but for 1 z > , which is physically impossible.Figure 5 illustrates the stability of these fixed points.Figures 5a, b and c include phase diagram of the nano-beam for 0 , 1 and 2 , respectively.In Figure 5a, it can be found that the first equilibrium position is a stable centre point (CP) and the second is an unstable saddle-node (SN).As shown in this Figure, there are two basins of attraction of stable centers and a basin of repulsion of unstable saddle node.The first basin of attraction of the first stable center is bounded by a closed orbit and the second basin of attraction of the second stable is an unbounded region.Depending on the location of the initial condition the system can be stable or unstable.Fig- ure 5b, confirms that there is not any types of fixed points at above or under the substrate for . As shown in Figure 5c, for , there are mathematically stable center and saddle node under the substrate which are physically impossible.
In addition as shown in Figure 4, in the state-control space, the stable and unstable branches of the fixed points with changing applied voltage meet together at a saddle-node bifurcation point.The voltage corresponding to the saddle-node bifurcation pointon the upper side of the substrate is a critical value, which is known as static pull-in voltage in the MEMS and NEMS Literature.In other words, when the applied voltage equal to the static pull-in voltage there is no any basin of stable attractors on the upper side of the substrate and the micro-beam is unstable for every initial condition.Now, with attention to Figure 5, the stable and unstable branches (S.B and U.S.B) of the bifurcation diagram can be determined.Also, it must be mentioned that for a given voltage there is a singular point (SP) at z=1.It must be mentioned that as seen in Figure 4, for V=0V, the stable fixed point is not displaced in not displaced position (z=0).This is due to the existence of the intermolecular force, which has displaced the mass initially.In Figure 6 the state-control space is plotted when 0 8 g nm = .
As shown, besides increasing the pull-in voltage, the voltage ranges in which no fixed points are observed have became smaller.It can be understood that by increasing the gap distance the voltage range with no fixed points tends to zero and even for some ranges four fixed points can be seen in state-control space.As an instant this problem can be seen in Figure 7, where the gap size is 0 10 g n m = .
Figures 8a, b and c illustrate the stability of these fixed points for 0 , 2 and 2.5 , respectively.As shown, for the stable center point and saddle node appear only on upper side of the substrate, whereas for  , respectively.It is understood that for all gaps above 15 nm there are three fixed points before pull-in voltage.One of them is under the substrate which is physically impossible, but other two points are above the substrate.On the other hand, for voltages higher than pull-in voltage there is only one fixed point which is under the substrate.To investigate the stability of these fixed points, the phase diagram for 0 15 g n m = is plotted in Figure 11.As shown, for V=0V, the first fixed point is a stable center point, the second is a saddle node which are on upper side of the substrate and the third one is a mathematically stable center which is physically impossible.

CONCLUSIONS
In presented work, gap dependent bifurcation behavior of an electrostatically-actuated gold nanobeam was studied.Both distributed and lumped models were introduced to explain the nano-beam deformation considering the couple stress theory.The results showed that fundamental frequency of the nano-beam can be gap dependent and the rate of the change in the frequency was higher for small gap regimes than large distances.In the bifurcation analysis, the following results were obtained: 1.In the both small and large gap regimes, for the voltages lower than the pull-in voltage, two stable and unstable fixed points appear on upper side of the beam.2. For the small gap regime, a voltage range can be found in which no fixed point appears, whereas for large gaps there is not such a range.3.For the large gap distances, for all voltages, there is one mathematically stable fixed point under substrate plate which is physically impossible, whereas for the small gaps we have two voltage ranges; in first range, there is not any fixed pint under the substrate and in second one two mathematically stable and unstable branches meet together in a saddle node.In this case the distance between two saddle nodes on the upper and lower sides of the substrate plate varies with changing gap size.
can be derived by dropping the time dependent terms from the dynamic equation of motion:

Figure 2 :
Figure 2: A lumped model of the nano-beam.

Figure 3 :
Figure 3: Variation of the nondimensional fundamental frequency versus Center gap; a) small distance gap, b) large distance gap.

Figure 4 :Figure 5 :
Figure 4: Variation of the center gap of the nano-beam with applied voltage for 0 6 g nm =

Figure 6 :
Figure 6: Variation of the center gap of the nano-beam with applied voltage for 0 8 g nm = .
center point and saddle node exist under the substrate which as mentioned is physically impossible.But, for 2 V V = as illustrated in Figure 8b, Latin American Journal of Solids and Structures 11 (2014) 2426-2443 both side of the substrate include a center point and a saddle node which is in agreement with Figure 7.

Figure 7 :Figure 8 :
Figure 7: Variation of the center gap of the nano-beam with applied voltage for 0 10 g n m =.

Figure 9 :
Figure 9: Variation of the center gap of the nano-beam with applied voltage for 0 15 g n m = .

Figure 10 :
Figure 10: Variation of the center gap of the nano-beam with applied voltage for 0 25 g n m = .

Figure 11 :
Figure 11: Phase diagram of the nano-beam for various initial conditions for 0 15 g n m = and V=0V.
non-dimensional form is given as:

Table 1 :
Geometrical and material properties of the gold nano-beam.