- Email: [email protected]
Nonlinear damping in doubly clamped beam resonators due to the attachment loss induced by the geometric nonlinearity
Journal of Sound and Vibration 372 (2016) 255–265
Contents lists available at ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
Nonlinear damping in doubly clamped beam resonators due to the attachment loss induced by the geometric nonlinearity André Gusso Departamento de Ciências Exatas-EEIMVR, Universidade Federal Fluminense, Volta Redonda, RJ 27255-125, Brazil
a r t i c l e i n f o
abstract
Article history: Received 21 September 2015 Received in revised form 2 February 2016 Accepted 25 February 2016 Handling Editor: L.N. Virgin Available online 11 March 2016
A nonlinear damping mechanism relevant for doubly clamped beam resonators vibrating transversally is proposed and investigated theoretically. The energy loss is a consequence of the axial stress induced along the beam due to the geometric nonlinearity. As the beam vibrates a time varying normal stress is induced at the attachment point which results in acoustic energy loss. Analytical expressions for the resulting amplitude dependent quality factor and the nonlinear damping parameter in a reduced order model are derived considering supports modeled as semi-spaces and as semi-infinite thin plates. The results are expected to be particularly relevant in the analysis of the nonlinear dynamics of suspended beam micro- and nanoresonators, but are not restricted to these particular devices, being valid for similar macroscopic systems. & 2016 Elsevier Ltd. All rights reserved.
Keywords: Nonlinear damping Clamping loss Geometric nonlinearity Microresonator Nanoresonator
1. Introduction The investigation of the nonlinear dynamics of micro- and nanodevices (MEMS and NEMS) has been receiving great attention because of the need to establish the dynamical range of devices required to operate in the linear regime, as well as, because of potential applications of the nonlinear phenomena in new devices [1,2]. The nonlinearity in MEMS and NEMS can originate from many sources, the most common being the electrostatic force and the geometric nonlinearity. More recently, nonlinear damping has also been considered as a relevant source of nonlinear behavior [2] and attracted special attention due to the experimental observation of strong nonlinear damping in micro and nanoresonators [3–5]. One interesting aspect of these observations is that well-known nonlinear damping mechanisms that are present in usual macroscopic engineering problems [6] are absent in these experiments. For instance, the drag force of fluids, that can be proportional to the velocity squared, did not contribute because the experiments were carried out in high vacuum. Hence, new nonlinear damping mechanisms must be responsible for the observed dissipation. We can expect that the large nonlinear dissipation observed in micro- and nanoresonators results from the simultaneous contribution of more than one nonlinear damping mechanism, in analogy to what has been observed with the dissipation in the linear regime [7]. In this last case, the linear damping can result from the contribution of several known mechanisms such as clamping loss, two-level systems (quantum dissipation), mechanical defects, thermoelastic damping [7] and noncontact friction [8,9]. It can also be expected that the observed large nonlinear damping could result from novel underlying physical process or from the same physical processes that lead to the linear damping mechanism mentioned above. There are also some recently proposed damping mechanisms that are specific to micro- and nanoresonators and that can result in nonlinear damping, but cannot explain the observed large nonlinear damping. A well-known case is that of squeeze film E-mail address: [email protected] http://dx.doi.org/10.1016/j.jsv.2016.02.044 0022-460X/& 2016 Elsevier Ltd. All rights reserved.
256
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
damping, which is relevant for beam resonators operating in a fluid and with a sufficiently small gap between the underlying beam surface and the substrate [1]. When the amplitude of vibration is comparable to the beam–substrate gap d the squeeze film _ damping force per unit length acting along the beam is proportional to the nonlinear damping term wðx; tÞ½d wðx; tÞ 3 , where wðx; tÞ denotes the transversal displacement at position x along the beam length. Two other mechanisms involving noncontact friction that, so far, have been investigated only in the linear regime, can also lead to nonlinear damping of the beam vibration [8,9]. In the noncontact friction mechanism discussed in Ref. [8] the damping originates from the energy loss mechanism observed when an atomic force microscopy tip was set to vibrate vertically on top of a substrate. The origin of this energy loss remains without a satisfactory theoretical explanation, however, the consequences of this energy loss mechanism can be approached phenomenologically, as done in Ref. [8]. The phenomenological models predict that the force per unit of length is proportional to _ wðx; tÞ½d wðx; tÞ α , where α varies between 3 and 4, depending on the specific phenomenological model being considered. For the noncontact acoustic energy loss mechanisms discussed in Ref. [9] we cannot define immediately a nonlinear force per unit length, however, a nonlinear damping is expected because the force acting on the surface of the substrate, and that leads to the acoustic energy loss, varies nonlinearly with the beam substrate separation. Nevertheless, these noncontact friction mechanisms cannot contribute to the observed nonlinear damping because they are relevant only when the gap between the resonator and the substrate is smaller than approximately 100 nm and the experiments involved gaps larger than a few hundreds of nanometers. So far, only two energy loss mechanisms leading to nonlinear damping in micro- and nanoresonators, and that could contribute to the observed large nonlinear damping, have been proposed. Zaitsev et al. [4], in an attempt to explain their experimental findings, considered to which extent the geometric nonlinearities can induce nonlinear damping due to internal viscoelastic damping. From their theoretical results they conclude that this mechanism could account for at most only about 1 percent of the measured nonlinear dissipation. The other mechanism emerged from the investigations of Croy et al. [10] of the vibrations of nanoresonators made of suspended graphene sheets and actuated electrostatically, like those investigated experimentally by Eichler et al. [3]. They modeled the graphene sheet using an elastic continuum model and derived a pair of coupled nonlinear partial differential equations (PDEs) as the equations of motion for the out-of-plane and in-plane vibrations. The dissipation of mechanical energy occurred because the out-of-plane motion was converted into in-plane phonons in the suspended region resulting in phonon emission (acoustic energy loss) through the supporting region. In this resonator this region was located at the surface of the substrate and extended over the area along which the graphene sheet was glued to this surface due to the van der Waals force [10]. In order to calculate the acoustic emissions through the support region it was necessary to discretize the Fourier transformed response of the substrate to the applied stress and obtain a linear system of equations for the in-plane response function that had to be solved numerically. Unfortunately, the parameters characterizing the linear and nonlinear dissipation cannot be easily obtained using this approach, and no analytical expressions for such parameters as a function of geometrical and material parameters have been derived. The mechanism of acoustic energy loss through the attachment region is, in fact, a well-known damping mechanism [11,12] in more conventional micro- and nanoresonators like those investigated by Zaitsev [4] and Imboden [5]. Such resonators consist of a suspended beam clamped to the substrate or other anchoring structures, and are one of the most commonly investigated micro- and nanoresonators for practical applications. However, only the acoustic emissions resulting in linear damping have been investigated theoretically for such resonators resorting either to purely analytical [11,12] or numerical methods, such as the finite element method [13], or a combination of both [14]. In this work we extend the conventional purely analytical treatment of clamping loss in doubly clamped beam resonators and include a source of acoustic energy loss that leads to nonlinear damping. Using this approach we can derive analytical expressions for the amplitude dependent quality factor and the nonlinear damping coefficient of doubly clamped beam resonators. The paper is organized as follows: in Section 2 we present the origin of the acoustic energy loss that results in the nonlinear damping and discuss the assumptions and approximations that have to be made in order to calculate the quality factor and the nonlinear damping parameter. In Section 3 we show how to obtain the approximate solutions of the nonlinear PDE required to calculate the relevant physical quantities from a corresponding reduced order model of the beam. In Sections 4 and 5 we derive analytical expressions for the amplitude dependent quality factor and the nonlinear damping coefficient, respectively. In Section 6 we analyze the relevance of the proposed nonlinear damping mechanism for realistic MEMS/NEMS resonators. In Section 7 we compare the predictions based on our analytical expressions with the experimental results obtained by Zaitsev et al. and conclude in Section 8.
2. Nonlinear damping mechanism When operating in the linear regime the dynamics of thin beams is well described by the Euler–Bernoulli beam equation. However, as the vibration amplitude increases the equation must be modified to include the effects of the elongation that become significant. If the amplitude of vibration is large, but sufficiently small that linear stress–strain law can be assumed, and that the effects of rotary inertia and transverse shear, as well as, nonlinearities due to inertia and curvature can be neglected, the resulting equation of motion is the following nonlinear PDE [1], ! Z EA l 0 2 00 ⁗ € þcw _ ¼ F ðx; t Þ; w dx w þ ρAw (1) EIw N 0 þ 2l 0
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
257
where w ¼ wðx; tÞ denotes the transversal displacement of the beam located between x¼0 and x ¼l. The overdots and primes represent derivatives with respect to time, t, and space, x, respectively, E denotes the Young modulus, I the geometric moment of inertia, A the cross-sectional area, and ρ the beam density. We have also considered that the beam can be subject to an axial load N0 due to residual stress. Included in the equation is a viscous damping term, which accounts for any existing linear damping mechanisms, and the forcing along the beam F ðx; tÞ. The nonlinear damping mechanism that we investigate originates from the geometric nonlinearity term Z EA l 0 2 00 00 w dx w ¼ Nðt Þw : (2) 2l 0 In this equation we have evidenced the time-varying axial force N(t) induced along the beam as it vibrates transversally. This time-varying force is transmitted to the point of attachment of the beam and is directly or indirectly responsible for the production of acoustic waves that carry away a fraction of the beam mechanical energy. In order to understand how this normal force can result into acoustic emissions we have to consider that there are typically two classes of anchoring configurations for suspended beam micro- and nanoresonators: we can have the beams anchored directly to the substrate (within the substrate) or to anchoring structures built on top of it (see Ref. [12] for more details). If the beam is attached to an anchoring structure built on top of the substrate the normal force N(t) is applied at some point of the structure pulling it and inducing, as reaction forces, time-varying shear stresses and bending moments along the region of the surface of the substrate where the anchoring structure is founded. These shear and bending moments then produce the emission of acoustic waves into the substrate. If the beam is clamped directly to the substrate the force N(t) results in a time-varying normal stress at the attachment point which induces the emission of acoustic waves directly into the substrate. Since in the case of anchoring structures built on top of the substrate the determination of the reaction forces and moments induced at the substrate can be rather complex, we restrict our analysis to beams attached directly to the substrate, a configuration that is commonly found in practical devices (examples can be found in Refs. [4,14]) and is much more accessible to a fully analytical investigation. For the configuration of a beam attached directly to the substrate it has been customary to model the beam as attached to either a semi-space or a semi-infinite thin plate [11,12,15] depending upon the details of the attachment region being investigated. Here we are going to investigate the resulting nonlinear damping considering both models for the attachment region, which are depicted in Fig. 1. Before we proceed we have to note that besides the normal force N(t) induced by the geometric nonlinearity, the in-plane and out-of-plane transverse vibrations of the beam induce shear forces and bending moments at the attachment point [11] that result in the well-known linear damping mechanism of clamping loss. While the combined effects of shear and bending may have to be considered in order to determine the total acoustic energy loss resulting in the linear damping [11], as in the case of beams supported by thin plates, the acoustic emissions due to N(t) can be analyzed independently of the emissions due to the other forces [16]. Assuming small nonlinearity, we can still characterize the damping of the beam vibration due to the acoustic emissions induced by N(t) by the corresponding quality factor defined as Q ¼ 2π
U U ¼ω ; ΔU Π
(3)
where U denotes the total mechanical energy of the beam vibrating at frequency ω, ΔU the energy dissipated per oscillation cycle, and Π the average emitted acoustic power. To calculate analytically these physical quantities we adopt the same assumptions made previously in Refs. [11,12,15,17]: (i) the dimensions of the beam are such that the wavelength of the emitted acoustic waves are much larger than the dimension of the attachment point, alloying the use of the point source approximation. (ii) In spite of the elastic vibrations at the attachment point, the vibrations of the beam are determined assuming perfectly rigid attachment, that means the beam satisfies the boundary conditions wð0; tÞ ¼ wðl; tÞ ¼ w0 ð0; tÞ ¼ w0 ðl; tÞ ¼ 0. (iii) The stresses required to satisfy the corresponding boundary conditions are those applied to the
Fig. 1. Schematic diagram showing one end of a beam attached to a semi-space (left) and a thin plate (right).
258
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
surface of the substrate in the region corresponding to the attachment point, and act as radiation sources of acoustic waves that are completely lost into the substrate or, more generally, into large surrounding structures. The calculation of U, ΔU and Π requires the knowledge of wðx; tÞ. Since we cannot obtain an exact solution of Eq. (1) we have to resort to an approximate method of calculation of wðx; tÞ as detailed in the next section.
3. Approximate solution to the PDE In order to choose an adequate method we note first that we are interested in calculating the dissipation and the quality factor for frequencies close to or at the resonant frequencies of the beam. Furthermore, for small nonlinearity, the expected resonant frequencies and modeshapes must be close to those obtained for the corresponding linear PDE, except under special circumstances, such as the occurrence of internal resonances [6], which we assume that are being circumvented with an adequate choice of the system parameters. In this case we can resort to the Garlekin method [1], which allow us, in a mathematically consistent manner, to reduce Eq. (1) to a single nonlinear ordinary differential equation (ODE) which contains adequate approximate expressions for the physical parameters required to calculate U and provides a simple procedure to calculate N(t), required to calculate Π. The Garlekin method is widely used in the modeling of micro and nanoresonators, and has been used to reduce the original nonlinear PDE to a single nonlinear ODE which was then solved numerically [18] or by approximate analytical methods [1,2]. Applying this method, we have that close to the resonance corresponding to the nth eigenmode of the corresponding linear PDE we approximate w by wðx; tÞ ¼ ϕn ðxÞun ðtÞ, where ϕn ðxÞ is nth eigenmode of the linear and homogeneous PDE obtained from Eq. (1) by suppression of the geometric nonlinearity, viscous damping and forcing. More specifically, the ϕn's are solutions of the following time-independent equation: EI ⁗ N 0 00 ϕ ϕ ω2n ϕn ¼ 0; ρA n ρA n where ωn denotes the corresponding eigenfrequencies. The solutions to the above equation have the general form ϕn ðxÞ ¼ coshðκn x=lÞ cos ðκn x=lÞ þ χ n sinhðκn x=lÞ sin ðκn x=lÞ :
(4)
(5)
The parameters κn and χn are functions of the axial stress N0 and, in general, can only be obtained numerically [1]. However, in the limiting case N0 ¼ 0 the ϕn ðxÞ are the well-known modeshapes of the Euler–Bernoulli equation, for which we have κ n ¼ 4:7300; 7:8532; 10:996… and χ n ¼ 0:9825; 1:0008; 0:99996…. The other relevant limiting case is obtained when N 0 -1, that means, when the axial force is much larger than pffiffiffithe elastic restoring force. In this limit the beam behaves like a string and the modeshapes have the simple form ϕn ðxÞ ¼ 2 sin ðnπx=lÞ. At this point we have defined the approximate spatial dependence of wðx; tÞ and we are left with the task of determining the approximate dynamical behavior of the system. The dynamics is contained in the function un(t), whose corresponding equation of motion is the final product of the Galerkin method. To obtain this equation applying this method we first substitute wðx; tÞ ¼ ϕn ðxÞun ðtÞ in Eq. (1). Eq. (4) is used to simplify the resulting equation which is then multiplied by ϕn ðxÞ and integrated along the beam length, and further simplified resorting to the orthogonality of the modeshapes Z l Z 1 ½ϕn ðxÞ2 dx ¼ l ½ϕn ðyÞ2 dy ¼ l; (6) 0
0
where we have evidenced the orthonormality of the modeshapes in the renormalized variable y ¼ x=l. Considering that the resonator is being driven by an harmonic force of the form F ðx; tÞ ¼ F 0 cos ðωn tÞ, which is the case, for instance, in the experiments of Refs. [3–5] we thus obtain the following nonlinear ODE: mu€ n þ2b11 u_ n þ k1 un þk3 u3n ¼ F cos ðωn tÞ;
(7)
where the notation is the same adopted by Zaitsev et al. [4] in order to facilitate the later comparison with their experimental results. In this equation we have that m ¼ ρAl, b11 ¼ cl=2, the effective spring constant is given by N0 I n þ I n1 ; 3 3 l
(8)
ðI n1 Þ2 :
(9)
ϕn ðyÞ dy;
(10)
½ϕ0n ðyÞ2 dy;
(11)
k1 ¼ mω2n ¼
EI l
and the nonlinear spring constant reads k3 ¼
EA 2l
3
n
The force now is given by F ¼ F 0 lI 0 and we have defined I n0 ¼ I n1 ¼
Z
1 0
Z
1 0
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
Z I n3 ¼
1 0
ϕ⁗ n ðyÞϕn ðyÞ dy:
259
(12)
Eq. (7) has the form of the well-known Duffing equation [6] for a nonlinear oscillator, and has been extensively studied. For the purpose of calculating U, required to derive Q predicted by the proposed damping mechanism, as well as others quantities like N(t) and ΔU, we need to know the un(t) which is obtained as the solution of Eq. (7). In order to define an adequate approximation to the function un(t), we first note that in common situations of practical interest the amplitude of motion of the beam is sufficiently small that the resonator performs quasi-harmonic oscillations, as a consequence of the harmonic forcing. That means the amplitude of the terms n 41 in the general steady state solution un ðtÞ ¼
1 X
am cos ðmωn t þ ϕm Þ;
(13)
m¼0
are small compared to a1, the amplitude of the harmonic term (note that a0 is not relevant for our problem). In fact, using the harmonic balance method, the authors of Ref. [19] have evidenced that even in the strong nonlinear regime, where the resonant curve is highly distorted, the amplitude of the next relevant higher harmonic does not exceed approximately 10 percent of the leading harmonic term. He have further checked the amplitude of the next higher harmonic solving numerically Eq. (7) for a resonator operating in the weak nonlinear regime, when the forcing results in a slightly tilted frequency–response curve. We have found that the amplitude of the next higher harmonic does not exceed 1 percent of a1. Therefore, considering that the resonators are operating in the weak nonlinear regime it suffices to take the harmonic approximation un ðtÞ ¼ u0 cos ðωn tÞ in order to calculate U, ΔU, and N(t) because such quantities are proportional to, at least, the square of the amplitude. We note that we can ignore any phase in the definition of un(t) because we are interested only in time averaged quantities.
4. Quality factor 4.1. Q for the proposed damping mechanism The acoustic energy loss in the proposed nonlinear damping mechanism results from the time-varying normal force N(t) applied to the clamping region. Substituting the harmonic approximation wðx; tÞ ¼ u0 ϕn ðxÞ cos ðωn tÞ into the expression for N(t) [see Eq. (2)] and neglecting the contribution of a constant force term we obtain N ðt Þ ¼
EAI n1 4l
2
u20 cos ð2ωn t Þ:
(14)
It is important to note that the forcing is harmonic and has a frequency that is twice that of the vibrating beam. The resulting emitted acoustic power is going to depend upon details of the supporting region. As mentioned previously, here we consider two cases: the beam attached to semispaces and to semi-infinite thin plates having the same thickness of the beam. The first case is an approximation to large anchoring structures. The second case represents approximately the conditions found in many actual micro and nanoresonators, including those investigated in Refs. [4,5], where the fabrication process results in cantilevered supporting structures. For the case of anchoring in the semispaces we use the expression, derived by Miller and Pursey [20], for the acoustic power emitted by a point source vibrating normal to the surface of a semispace. For an harmonically varying force of the form FðtÞ ¼ f cos ðωtÞ they employed the admittance method obtaining pffiffiffiffiffiffiffiffiffiffi ϕ ðγÞ ρC 11 2 2 ω f (15) Πs ¼ 1 4π C 244
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where ρ is the density, γ ¼ C 11 =C 44 , with C 11 ¼ Eð1 νÞ=ð1 ν 2ν2 Þ and C 44 ¼ E=2ð1 þνÞ, where E corresponds to the Young modulus and ν the Poisson constant of the isotropic material of the beam and semispace. The function ϕq ðγÞ is defined as "Z pffiffiffiffiffiffiffiffiffiffiffiffiffi # 1 q p p2 1 dp ; (16) ϕq ðγ Þ ¼ Im F 0 ðp; γÞ 0 with F 0 ðp; γÞ ¼ ð2p2 γ 2 Þ2 4p2
qffiffiffiffiffiffiffiffiffiffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffi p2 1 p2 γ 2 :
(17)
Before we proceed, we note that ϕq ðγÞ has been evaluated differently in the literature as we have already pointed previously [9] for the cases q ¼1 and 3. In Ref. [9] explicit expressions for ϕ1 ðγÞ and ϕ3 ðγÞ have been obtained following the prescriptions of Miller and Pursey [20], and can be very easily generalized for arbitrary positive integer q.
260
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
The emitted acoustic power for a point source acting on the edge of a semi-infinite thin plate of thickness h was derived by Hao and Xu [12] and can be cast in the form 2
Πp ¼
2ϕ0 ðγÞ 1 þν f ω : π Eð1 νÞ h
(18)
The quality factor defined in Eq. (3) must be rewritten in terms of the average emitted powers Π s;p as Q ¼ ωU=ð2Π s;p Þ ¼ ωn U=Π s;p , where the factor of 2 in the denominator, which is a consequence of the energy being lost by both ends of the beam, cancels with the factor of two coming from the fact that ω ¼ 2ωn according to Eq. (14). We can obtain U from the maximum potential energy in terms of the parameters in Eq. (7) as u2 k3 (19) U ¼ 0 k1 þ u20 : 2 2 Using k1 ¼ mω2n and k3 defined in Eq. (9) in this expression for U, the resulting quality factor when the beam is attached to semi-spaces is pffiffiffi 5 4πC 244 ρl ωn πC 244 l pffiffiffiffiffiffiffiffiffiffiffi p ffiffiffiffiffiffiffiffi 2 þ ¼ Q sn;1 þ Q sn;2 ; (20) Q sn ¼ n 2 2 C 11 E A u0 ϕ1 ðγÞ ρ C 11 EA ωn ϕ1 ðγÞðI 1 Þ where we have left the final expression in terms of ωn for the sake of clarity, however, using Eq. (8) this frequency can be s rewritten in terms of the parameters of the beam. We can see that Qn can be separated in two terms. Q sn;1 is the amplitude s 2 dependent contribution, and goes as u0 . The other term contributing to Qn is amplitude independent. It results from the ratio between the energy stored in the nonlinear component of the spring constant to the energy lost through the nonlinear s damping mechanism, and determines the smallest value of Qn that is attained as the amplitude of vibration increases. We note, however, that as the amplitude increases higher order contributions to the effective elastic constant and the dissipated energy can become relevant and change the minimum quality factor before the beam suffers permanent deformation or rupture. Using the expression for Πp, Eq. (18), the quality factor Q ¼ ωn U=Π p for a beam clamped to thin plates of thickness h results to be 5
Q pn ¼
πρð1 νÞhl ω2n ϕ0 ðγÞðI n1 Þ2 Eð1 þ νÞA
u20
þ
πð1 νÞhl ¼ Q pn;1 þQ pn;2 ; 4ϕ0 ðγÞð1 þ νÞA
(21)
Again we have the amplitude dependent and the amplitude independent contributions to the total quality factor. The more s p noticeable difference between Qn and Qn is the dependence on the angular frequency ωn. Partially as a consequence of this s p distinct dependence on ωn, Qn and Qn have very distinct dependences on the geometrical parameters of the beam. The dependence on u0 2 observed in our expressions for Q sn;1 and Q pn;1 is the same one obtained in previous works that considered the addition of nonlinear damping terms of the form b31 u2n u_ n or b32 u_ 3n to the Duffing-like equations describing the nonlinear dynamics of resonant beams [2–5,10]. We note, however, that as no nonlinear damping was present in the original PDE, Eq. (1), such a term is also missing in the reduced order model, Eq. (7). Unfortunately, our approach does not allow us to determine which one of these two terms should be present in the reduced order model. Nevertheless, many aspects of the dynamics of the resonators do not depend on which of these two terms is present in the equation of motion, as already noticed in Refs. [2,4], because their effects are combined linearly in the same effective nonlinear parameters. In order to have a better understanding of the results we have obtained, in the next subsection we are going to add to Eq. (7) a term of the form b31 u2n u_ n þ b32 u_ 3n , and consider its effects on the quality factor of a resonator described by a modified equation of motion for the reduced order model, namely mu€ n þ ð2b11 þ b31 u2n þ b32 u_ 2n Þu_ n þk1 un þ k3 u3n ¼ F cos ðωn tÞ:
(22)
4.2. Q in the reduced order model We can easily calculate Q as a function of the parameters in Eq. (22) in the harmonic approximation. We already have U, given in Eq. (19), and the missing expression for ΔU is obtained from the work done by the dissipative forces along one oscillation cycle as b31 2 3 u0 þ b32 ω2n u20 : (23) ΔU ¼ πωn u20 2b11 þ 4 4 Therefore, the quality factor in terms of the parameters of the reduced order model is Q ¼ 2π
U ¼ ΔU
k1 þ ðk3 =2Þu20
: u2 ωn 2b11 þ b31 þ 3b32 ω2n 0 4
(24)
This result coincides with an analogous expression obtained in Ref. [10] using the slow envelope approximation. It is interesting to analyze some general features of this expression for Q. In the absence of nonlinear terms, we recover the usual expression for the
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
261
1.0
max
0.8 0.6 0.4 0.2
0.97
0.98
0.99
1
1.01
1.02
1.03
Fig. 2. Normalized frequency–response curves obtained by solving Eq. (26), with Q 0 ¼ 100, and Γ ¼ 10 8 for different forcing A2 ¼ 1; 2; 5, and 10 (Γ and A in arbitrary units). The broader curves, that correspond to a smaller quality factor, are obtained for larger forcing A that results in larger amplitudes.
amplitude independent Q. If, in addition, we had only the nonlinear spring term, Q would increase with the amplitude, because more mechanical energy is stored in the resonator for a given amplitude. When only b11 ¼ 0, the quality factor can be separated in a constant term and one that depends on the amplitude as u0 2 . Since this result is important for the comparison with Q predicted by the nonlinear acoustic energy loss mechanism, we present it explicitly Q ¼ Q 1 þQ 2 ¼
4k1 2
ωn bu0
þ
2k3 ; ωn b
(25)
where b ¼ b31 þ 3b32 ω2n , and we have separated Q in the amplitude dependent and amplitude independent terms, Q1 and Q2, respectively. The concept of an amplitude dependent quality factor was first introduced by Zhong-Kun and Luo [21] during the theoretical and experimental investigation of the quality factor of a torsion pendulum, where the effects of nonlinear damping were considered. More recently, Eichler et al. [3] observed experimentally a dependence with the driving amplitude of the quality factor of electrostatically driven nanotubes and graphene sheets. They developed a theoretical model which incorporated nonlinear damping and led to an expression for Q that coincides with Q1 presented in Eq. (25). The concept of an amplitude dependent Q is clarified in Fig. 2 where we have plotted the theoretical frequency–response curves for a system described by Eq. (22), obtained for different values of the forcing amplitude F and a fixed nonlinear damping coefficient b, and assuming k3 ¼ 0. The curves have been obtained by solving for u0 the equation u20 ¼
A2 2 2 ω ωn 3 k3 2 1 Q 0 1 þΓu20 u0 þ 8 k1 2 ωn
(26)
where A ¼ F=ð2k1 Þ, Q 0 1 ¼ 2b11 ωn =k1 , and Γ ¼ bωn =ð4k1 Þ. This equation for u0 was adapted to our notation from the result derived in Ref. [2] using secular perturbation theory for a nonlinear Duffing oscillator with nonlinear damping. It is interesting to note that the authors of Ref. [2] have derived an expression which evidenced the existence of an effective 1 amplitude dependent quality factor Q eff ¼ Q 0 1 þΓu20 , however, this interpretation was not given in that work.
5. Nonlinear damping coefficient s
p
From the above results for Qn and Qn we can derive the more fundamental nonlinear damping parameter b. This can be done very easily comparing the expressions for Q sn;2 and Q pn;2 with that of Q2 defined in Eq. (25) and the results are pffiffiffiffiffiffiffiffiffiffi ϕ ðγÞðI n1 Þ2 ρC 11 E2 A2 s ; (27) bn ¼ 1 4 πC 244 l for clamping on a semispace and p
bn ¼
4ϕ0 ðγÞðI n1 Þ2 ð1 þ νÞEA2 4
πð1 νÞhl ωn s
;
(28)
for clamping on a thin plate. It is interesting to note that bn does not depend explicitly on the frequency and, therefore, has s only a small dependence on the axial force N0 through the term ðI n1 Þ2 , while bn has a stronger dependence on N0 through the term ωn appearing in the denominator. The two b's also have very distinct dependencies on the geometrical parameters of the beam. In order to see that, let us s 2 4 consider a beam with rectangular cross section having width w and thickness h. In this case we have bn p h w2 =l . The
262
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265 p
particular dependence of bn on the geometrical parameters depends upon N0 and the polarity ofp the beam vibration. Let us ffiffiffiffiffiffiffiffiffiffiffiffi consider the two limiting cases for N0 we have considered previously. For N 0 ¼ 0 we have ωn ¼ EI=ρAðκ n =lÞ2 , with A ¼ hw 3 3 p 2 and I ¼ wh =12 in the case of out-of-plane vibrations, and I ¼ hw =12 for in-plane vibrations. Therefore, bp w2 =l np ffiffiffiffiffiffiffiffiffiffiffi ffi for outp 2 of-plane vibrations and bn p hw=l for in-plane vibrations. In the limit N 0 -1 we have that ωn ¼ ðnπ=lÞ N=ρA and, conpffiffiffiffiffiffiffiffiffiffiffi p p 3 2 sequently, bn p h w5 =l , regardless of the polarization. We can see that in any of the cases we have considered that bn and s bn have quite distinct dependencies on h; w and l. This fact evidences the need to model adequately the clamping region, since discrepancies of orders of magnitude can be obtained in the theoretical estimation of the nonlinear damping parameter and, consequently, on the quality factor. It is also evident that this acoustic nonlinear damping mechanism is more relevant for short beams with large cross sectional area, independently of the specific clamping.
6. Analysis of the results In this section we briefly analyze the potential relevance of the proposed nonlinear damping mechanism to realistic MEMS/NEMS resonators. Let us consider a MEMS resonator made from polysilicon (E ¼170 GPa, C 44 ¼ 69:7 GPa, C 11 ¼ 194:1 GPa, and ρ ¼ 2:3 103 kg m 3 ), which is a common material for MEMS/NEMS fabrication, with thickness h ¼ 1 μm and l ¼ 30 h, w ¼ 5 h, and Q 0 ¼ 103 . In Fig. 3(a) we plot the frequency response curve obtained solving Eq. (26) for the first mode (n ¼1), using k1 and k3 defined in Eqs. (8) and (9), respectively and Γ calculated using our analytical expressions for bs1 and bp1. The continuous lines represent the response in the absence of nonlinear dissipation [Γ ¼ 0 in Eq. (26)]. The dashed curves have been obtained assuming the resonator is clamped to a plate. We can see that there is a relevant change in the frequency response curve, particularly at larger amplitudes of vibration, which present a decrease of the amplitudes close to its maximum. The curves obtained considering the clamping to a semi-space were so close to those with no dissipation that they have been omitted. A more general result comes from the comparison of the additional dissipation that results from the nonlinear damping compared to the dissipation due to linear damping. While the linear damping Q0 results in a dissipation δ0 ¼ Q 0 1 that is independent of the amplitude of vibration, the nonlinear damping results in a dissipation that grows with the amplitude. The total dissipation is the sum of the linear and nonlinear contributions and results in the effective quality factor Q eff s defined after Eq. (26). For resonators made from polysilicon we obtain the following expressions for the dissipation from Q1 p and Q1 ! 3 3 1 1 0:131l 0:144l þ δs1 ¼ s ¼ ; (29) 2 Q1 u20 w h w and δp1
1 ¼ p¼ Q1
2
1:14l 1:26h l þ 2 w u0 w
!1 :
(30)
From these equations we can see that the dissipation is independent of the actual dimensions and amplitude of vibration of the resonators, and once the relative dimensions are fixed, the dissipation becomes a function only of the normalized p displacement u0 =h. For this reason, in Fig. 3(b) we plot the nonlinear dissipation δnl ¼ δs1 or δ1 as a function of u0 =h for two
40
u0
30 20 10 0
0.996
0.998
1
1.002
1.004
10
3
10
4
10
5
10
6
0.02
0.04
0.06
u0
0.08
0.10
Fig. 3. (a) Frequency–response curves for a realistic MEMS resonator made from polysilicon with dimensions h ¼ 1 μm, l ¼ 30 μm, and w ¼ 5 μm oscillating in its first mode n¼1 and having an intrinsic quality factor Q 0 ¼ 103 . The lower, middle, and higher curves have been obtained for a forcing F ¼ 2 10 8 ; 3 10 8 , and 4:5 10 8 N, respectively. Continuous lines represent the response in the absence of nonlinear damping (b ¼0) and the p dashed curves have been obtained considering the nonlinear damping assuming clamping to a plate ðb ¼ b1 Þ. (b) Nonlinear dissipation δnl ¼ 1=Q nl as a function of the normalized amplitude of vibration u0 =h in the first mode (n¼ 1) for a resonator clamped to a plate (continuous) and to a semi-space (dashed). The black and gray curves have been obtained for resonators with relative dimensions l ¼ 10 h, and w ¼ 3 h and l ¼ 30 h, and w ¼ 5 h, respectively.
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
263
representative sets of relative dimensions, namely, l ¼ 30 h, and w ¼ 5 h and l ¼ 10 h, and w ¼ 3 h. The aspect ratio l=h ¼ 10 assumed in this last case is approximately the smallest value for which results obtained using the Euler–Bernoulli beam theory can be trusted [22]. The range of u0 =h was chosen as representative of a region where the amplitudes are sufficiently large that nonlinear effects are relevant, but the resonators are not in a strong nonlinear regime. For example, in Fig. 3(a) the maximum displacement for each curve ranges from u0 =h 0:015 up to u0 =h 0:035. The last maximum value corresponds to the amplitude close to the onset of bistability for the resonator considered for Fig. 3(a) with Q 0 ¼ 103 . This amplitude agrees well with the predicted theoretical value of the critical amplitude for the onset of bistability which, ignoring the small effect of nonlinear damping and using our notation, is given by [2] 34 sffiffiffiffiffiffiffiffi 4 k1 c u0 ¼ : (31) 3 Qk3 For the resonator considered for the results presented in Fig. 3(a) it results to be uc0 ¼ 0:029 h. According to this equation, if Q0 were of the order of 102 the critical amplitude would increase to approximately uc0 ¼ 0:092 for the same resonator, a prediction that is slightly smaller than the results we have obtained by solving Eq. (26) numerically, and that was used to set the maximum value of u0 =h in Fig. 3(b). Therefore we can compare the dissipation presented in Fig. 3(b) with the linear dissipation of a given device and estimate the potential relevance of the nonlinear damping mechanism. We can, for example, see that for u0 =h 0:035, at the onset of bistability, a nonlinear dissipation of δp1 10 4 is obtained for a resonator with l ¼ 30 h and w ¼ 5 h, which amounts to about 10 percent of the dissipation δ0 ¼ 10 3 of a resonator with Q 0 ¼ 103 . For the shorter resonator the dissipation increases to δp1 2:5 10 4 , and amounts to 25 percent of the linear dissipation. We can also see that the dissipation for beams clamped to semi-spaces is much smaller, and would amount to no more than a few percent of the linear dissipation considered above. Similar results hold for the nonlinear dissipation predicted at the onset o bistability for resonators with Q 0 ¼ 102 as can be seen from the curves at u0 =h 0:1.
7. Comparison with experiment Since our theoretical investigation has been motivated by recent measurements of nonlinear damping in micro and nanoresonators, it is interesting to evaluate the relevance of this nonlinear energy loss channel to the observed large nonlinear damping. For that purpose, we compare our predictions with the experimental results of Zaitsev et al. [4]. In this work detailed experimental results are presented for a rectangular cross section suspended beam microresonator made of the Pd0.15Au0.85 alloy, and dimensions l¼125 μm, w¼ 1 μm, and h¼0.3 μm. The microresonator was excited in its first mode, whose resonant frequency was measured to be f 1 ¼ 885:5 kHz, and performed in-plane vibrations (parallel to the dimension w). They have extracted experimentally several physical parameters that are related to the parameters defined in Eq. (22). One particularly relevant parameter is p, that is related to the coefficients in Eq. (22) by p¼
bωn 2 1 ¼ ; 3k3 3 Q 2
(32)
where the last equality results from the comparison of the definition of p with the expression for Q given in Eq. (25). Considering that the wet etching thins the support region close to the attachment region [23], we approximate the support structure by a thin plate in order to derive an analytical expression for p. Substituting into the previous equation the expression for Q pn;2 we have pt ¼
8ϕ0 ðγÞð1 þ νÞ w : 3πð1 νÞ l
(33)
In order to calculate pt we are going to consider that ν ¼ 0:43. This value is the weighted average (weighted by the atomic fractions) of the corresponding values for ν of the bulk Pd ðν ¼ 0:39Þ and that of the bulk Au ðν ¼ 0:44Þ. For ν ¼ 0:43 we have that γ ¼ 2:854 and ϕ0 ð2:854Þ ¼ 0:10996, and Eq. (33) predicts the value pt ¼ 0:002. This is only a small fraction of the value found experimentally [4] p¼ 0.109. We have to note that this value for pt was obtained assuming that the elastic properties of the beam are those of an homogeneous beam, having the properties of the bulk PdAu alloy. However, the fabrication process of the microbeams may result in structures with effective properties that differ significantly from those of bulk materials. For instance, the deposition process of thin films may result in PdAu alloys with densities and elastic properties that differ significantly from those of the bulk material [24]. Furthermore, the electron beam lithography used in the microfabrication process resulted in the deposition of an amorphous carbon layer on the beam surface [23]. This material typically has a density much smaller than that of PdAu alloys but a much larger Young modulus. The theoretical modeling of the resulting mechanical system can be rather complex, however, it is usually the case that effective mechanical properties can be defined for the system such as an effective flexural rigidity and effective density [1]. Therefore, a more precise comparison between theory and experiment would be obtained using the relevant effective physical properties of the microbeam. Fortunately, the detailed experimental analysis of Zaitsev et al. [4] provides us with the possibility to determine the experimental value of the effective nonlinear spring constant k3, which is a key ingredient in the derivation of pt.
264
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
Considering the value for the parameter α3 ¼ k3 =m estimated by Zaitsev et al. [4], α3 2 1023 m 2 s 2 , we can extract an approximate experimental value for k3. We estimate the mass of the PdAu alloy as was done by Zaitsev et al. [4] and find m ¼ 6:8 10 13 kg. This value is obtained considering the dimensions of the beam given above, and assuming a density for the PdAu alloy of ρ ¼ 18:2 103 kg m 3 , which is the weighted average of the Au density of 19:3 103 kg m 3 and that of e Pd, 12:0 103 kg m 3 . We thus obtain the experimental value k3 ¼ 1:36 1011 N m 3 . This experimental value can be t compared with the expected theoretical value k3, calculated using the definition presented after Eq. (7). For that purpose we use the Young modulus E ¼95 GPa, estimated from the experimental results [25] for E of PdxAu1 x alloys with Pd atomic fraction x. We have also considered the fact that the beam dynamics is dominated by the axial stress since we have estimated theoretically that the contribution of the elastic deformation term to k1 represents only about 2 percent of its total t e t experimental value k1 ¼ mω2 22 N m 1 . Therefore, calculating k3 in the limit N 0 -1 we obtain k3 ¼ 5 1011 N m 3 , which is considerably larger than the experimental value. e The expression for pe, the theoretical value estimated using the experimental value of k3, is 2
pe ¼
4π 3 ϕ0 ðγÞEð1 þνÞhw e 4 3k3 ð1 νÞl
:
(34)
e
Considering the value for k3 given above, we obtain pe ¼ 0:01. This is significantly larger than the purely theoretical value and is a better estimate of the actual value of this parameter as a result of the proposed nonlinear damping mechanism. We can also evaluate the more fundamental parameter b and compare with its experimental value. The experimental b is obtained from the parameter γ 3 ¼ b=m estimated by Zaitsev et al. [4] to be γ 3 ¼ 1 1016 m 2 s 1 . Therefore, we obtain e b ¼ 6:8 103 kg m 2 s 1 . The theoretical value, obtained from Eq. (28), is independent of k3 and, because we use as input the experimental value of the frequency, can be considered as the more realistic estimate of b. We thus have t b ¼ 0:72 103 kg m 2 s 1 . These theoretical results for p and b indicate that the acoustic nonlinear damping mechanism can be responsible for about 10 percent of the observed nonlinear damping. While this contribution can change if we perform a more detailed analysis of the experimental setup, it seems that, in this particular microresonator, the acoustic nonlinear damping can be relevant but it is not the prevailing nonlinear energy loss channel and other unknown nonlinear damping mechanisms should be present that must still be unveiled.
8. Conclusion While we have derived analytical expressions for Q and b for two specific models of the clamping region our approach can be applied if more specific models of the clamping region are required for a precise calculation of the nonlinear damping. The fundamental idea in our approach is to recognize that a beam vibrating with large amplitude develops a timevarying axial stress that is transmitted to the clamping region. For more general beams, for instance, with a variable crosssection, N(t) can be obtained analytically for a great variety of beam designs, or it could be determined resorting to well established computational methods, such as the finite element method. The calculation of the emitted acoustical power due to the resulting N(t) for more general configurations of the clamping region is certainly more complicated. However, computational approaches have been developed that can use both the analytical or fully computational results for N(t) as input for the calculation of the dissipated acoustical energy [13,16]. We believe that the use of the fundamental idea presented here in conjunction with other computational methods can be helpful in determining the induced nonlinear acoustical energy loss for more complex microscopic or macroscopic structures. It is also important to note that the damping resulting from the proposed nonlinear mechanism must, in general, be combined with other linear or nonlinear damping mechanisms relevant for a given suspended beam resonator. Having the results for the corresponding quality factors for each damping mechanism, the total loss (Q 1) is simply the sum of the loss predicted for each damping mechanism [7]. To conclude we note that the acoustic nonlinear damping mechanism proposed here has a strong dependence on the s p geometrical parameters of the doubly clamped resonators as can be seen, for instance, from the expressions for bn and bn. Furthermore, this dependence changes considerably depending on the geometry of the supporting region. As a consequence, the relevance of this mechanism in a particular setup is going to depend on the geometry and on details of the clamping region. This situation is analogous to that observed for several linear damping mechanisms [9,7,11,12], whose individual contribution to the total linear damping depends on the detailed geometry of the resonator. In the particular case of the proposed mechanism the expressions for Q and b indicate that it can be the prevailing source of nonlinear damping for short resonators with large cross sectional area, and should be considered in the analysis of the nonlinear dynamics of such devices.
Acknowledgments The author is thankful to S. Zaitsev for providing valuable information about the experiment reported in Ref. [4].
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
265
References [1] M.I. Younis, MEMS Linear and Nonlinear Statics and Dynamics, Springer, New York, USA, 2011. [2] R. Lifshitz, M.C. Cross, Nonlinear dynamics of nanomechanical resonators, in: G. Radons, B. Rumpf, H.G. Schuster (Eds.), Nonlinear Dynamics of Nanosystems, Wiley-VCH Verlag, Weinheim, Germany, 2010, pp. 221–266. [3] A. Eichler, J. Moser, J. Chaste, M. Zdrojek, I. Wilson-Rae, A. Bachtold, Nonlinear damping in mechanical resonators made from carbon nanotubes and graphene, Nature Nanotechnology 6 (2011) 339–342. [4] S. Zaitsev, O. Shtempluck, E. Buks, O. Gottlieb, Nonlinear damping in a micromechanical oscillator, Nonlinear Dynamics 67 (2012) 859–883. [5] M. Imboden, O.A. Williams, P. Mohanty, Nonlinear dissipation in diamond nanoelectromechanical resonators, Applied Physics Letters 102 (2013) 103502; M. Imboden, O.A. Williams, P. Mohanty, Observation of nonlinear dissipation in piezoresistive diamond nanomechanical resonators by heterodyne down-mixing, Nano Letters 13 (2013) 4014–4019. [6] A.H. Nayfeh, D.T. Mook, Nonlinear Oscillations, Wiley, New York, USA, 1979. [7] M. Imboden, P. Mohanty, Dissipation in nanoelectromechanical systems, Physics Reports 534 (2014) 89–146. [8] A. Gusso, Phenomenological modeling of long range noncontact friction in micro- and nanoresonators, Journal of Applied Physics 110 (2011) 064512. [9] A. Gusso, Energy loss mechanism for suspended micro- and nanoresonators due to the Casimir force, Physical Review B 81 (2010) 035425; A. Gusso, Acoustic electromechanical energy loss mechanism for suspended micro- and nanoelectromechanical resonators, Applied Physics Letters 96 (2010) 193504. [10] A. Croy, D. Midtvedt, A. Isacsson, J.M. Kinaret, Nonlinear damping in graphene resonators, Physical Review B 86 (2012) 235435. [11] J.A. Judge, D.M. Photiadis, J.F. Vignola, B.H. Houston, J. Jarzynski, Attachment loss of micromechanical and nanomechanical resonators in the limits of thick and thin support structures, Journal of Applied Physics 101 (2007) 013521. [12] Z. Hao, Y. Xu, Vibration displacement on substrate due to time-harmonic stress sources from a micromechanical resonator, Journal of Sound and Vibration 322 (2009) 196–215. [13] D.S. Bindel, S. Govindjee, Elastic PMLs for resonator anchor loss simulation, International Journal for Numerical Methods in Engineering 64 (2005) 789–818. [14] G.D. Cole, I. Wilson-Rae, K. Werbach, M.R. Vanner, M. Aspelmeyer, Phonon-tunneling dissipation in mechanical resonators, Nature Communications 231 (2) (2011) 1–8. [15] M.C. Cross, R. Lifshitz, Elastic wave transmission at an abrupt junction in a thin plate with applications to heat transport and vibrations in mesoscopic systems, Physical Review B 64 (2001) 085324. [16] B. Chouvion, S. Mcwilliam, A.A. Popov, C. Fox, Review and comparison of different support loss models for microelectromechanical systems resonators undergoing in-plane vibration, Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science 226 (2011) 283–295. [17] D.M. Photiadis, J.A. Judge, Attachment losses of high Q oscillators, Applied Physics Letters 85 (3) (2004) 482–484. [18] T.D. Amorim, W.G. Dantas, A. Gusso, Analysis of the chaotic regime of MEMS/NEMS fixed-fixed beam resonators using an improved 1DOF model, Nonlinear Dynamics 79 (2015) 967–981. [19] Z.K. Peng, Z.Q. Lang, S.A. Billings, G.R. Tomlinson, Comparisons between harmonic balance and nonlinear output frequency response function in nonlinear system analysis, Journal of Sound and Vibration 311 (2008) 56–73. [20] G.F. Miller, H. Pursey, The field and radiation impedance of mechanical radiators on the free surface of a semi-infinite isotropic solid, Proceedings of the Royal Society of London Series A 223 (1954) 521–541; G.F. Miller, H. Pursey, On the partition of energy between elastic waves in a semi-infinite solid, Proceedings of the Royal Society of London Series A 233 (1955) 55–69. [21] Z.-K. Hu, J. Luo, Amplitude dependence of quality factor of the torsion pendulum, Physics Letters A 268 (2000) 255–259. [22] K.F. Graff, Wave Motion in Elastic Solids, Dover, New York, USA, 1991. [23] S. Zaitsev, Forced and self-excited oscillations in nonlinear micromechanical resonators, PhD Thesis, Technion, Israel, 2011. [24] S. Yalcin, R. Avci, Characterization of PdAu thin films on oxidized silicon wafers: interdiffusion and reaction, Applied Surface Science 214 (2003) 319–337. [25] J. Lohmiller, Investigation of deformation mechanisms in nanocrystalline metals and alloys by in situ synchrotron X-ray diffraction, Scientific Publishing, Karlsruher, Germany, 2012.
Contents lists available at ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
Nonlinear damping in doubly clamped beam resonators due to the attachment loss induced by the geometric nonlinearity André Gusso Departamento de Ciências Exatas-EEIMVR, Universidade Federal Fluminense, Volta Redonda, RJ 27255-125, Brazil
a r t i c l e i n f o
abstract
Article history: Received 21 September 2015 Received in revised form 2 February 2016 Accepted 25 February 2016 Handling Editor: L.N. Virgin Available online 11 March 2016
A nonlinear damping mechanism relevant for doubly clamped beam resonators vibrating transversally is proposed and investigated theoretically. The energy loss is a consequence of the axial stress induced along the beam due to the geometric nonlinearity. As the beam vibrates a time varying normal stress is induced at the attachment point which results in acoustic energy loss. Analytical expressions for the resulting amplitude dependent quality factor and the nonlinear damping parameter in a reduced order model are derived considering supports modeled as semi-spaces and as semi-infinite thin plates. The results are expected to be particularly relevant in the analysis of the nonlinear dynamics of suspended beam micro- and nanoresonators, but are not restricted to these particular devices, being valid for similar macroscopic systems. & 2016 Elsevier Ltd. All rights reserved.
Keywords: Nonlinear damping Clamping loss Geometric nonlinearity Microresonator Nanoresonator
1. Introduction The investigation of the nonlinear dynamics of micro- and nanodevices (MEMS and NEMS) has been receiving great attention because of the need to establish the dynamical range of devices required to operate in the linear regime, as well as, because of potential applications of the nonlinear phenomena in new devices [1,2]. The nonlinearity in MEMS and NEMS can originate from many sources, the most common being the electrostatic force and the geometric nonlinearity. More recently, nonlinear damping has also been considered as a relevant source of nonlinear behavior [2] and attracted special attention due to the experimental observation of strong nonlinear damping in micro and nanoresonators [3–5]. One interesting aspect of these observations is that well-known nonlinear damping mechanisms that are present in usual macroscopic engineering problems [6] are absent in these experiments. For instance, the drag force of fluids, that can be proportional to the velocity squared, did not contribute because the experiments were carried out in high vacuum. Hence, new nonlinear damping mechanisms must be responsible for the observed dissipation. We can expect that the large nonlinear dissipation observed in micro- and nanoresonators results from the simultaneous contribution of more than one nonlinear damping mechanism, in analogy to what has been observed with the dissipation in the linear regime [7]. In this last case, the linear damping can result from the contribution of several known mechanisms such as clamping loss, two-level systems (quantum dissipation), mechanical defects, thermoelastic damping [7] and noncontact friction [8,9]. It can also be expected that the observed large nonlinear damping could result from novel underlying physical process or from the same physical processes that lead to the linear damping mechanism mentioned above. There are also some recently proposed damping mechanisms that are specific to micro- and nanoresonators and that can result in nonlinear damping, but cannot explain the observed large nonlinear damping. A well-known case is that of squeeze film E-mail address: [email protected] http://dx.doi.org/10.1016/j.jsv.2016.02.044 0022-460X/& 2016 Elsevier Ltd. All rights reserved.
256
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
damping, which is relevant for beam resonators operating in a fluid and with a sufficiently small gap between the underlying beam surface and the substrate [1]. When the amplitude of vibration is comparable to the beam–substrate gap d the squeeze film _ damping force per unit length acting along the beam is proportional to the nonlinear damping term wðx; tÞ½d wðx; tÞ 3 , where wðx; tÞ denotes the transversal displacement at position x along the beam length. Two other mechanisms involving noncontact friction that, so far, have been investigated only in the linear regime, can also lead to nonlinear damping of the beam vibration [8,9]. In the noncontact friction mechanism discussed in Ref. [8] the damping originates from the energy loss mechanism observed when an atomic force microscopy tip was set to vibrate vertically on top of a substrate. The origin of this energy loss remains without a satisfactory theoretical explanation, however, the consequences of this energy loss mechanism can be approached phenomenologically, as done in Ref. [8]. The phenomenological models predict that the force per unit of length is proportional to _ wðx; tÞ½d wðx; tÞ α , where α varies between 3 and 4, depending on the specific phenomenological model being considered. For the noncontact acoustic energy loss mechanisms discussed in Ref. [9] we cannot define immediately a nonlinear force per unit length, however, a nonlinear damping is expected because the force acting on the surface of the substrate, and that leads to the acoustic energy loss, varies nonlinearly with the beam substrate separation. Nevertheless, these noncontact friction mechanisms cannot contribute to the observed nonlinear damping because they are relevant only when the gap between the resonator and the substrate is smaller than approximately 100 nm and the experiments involved gaps larger than a few hundreds of nanometers. So far, only two energy loss mechanisms leading to nonlinear damping in micro- and nanoresonators, and that could contribute to the observed large nonlinear damping, have been proposed. Zaitsev et al. [4], in an attempt to explain their experimental findings, considered to which extent the geometric nonlinearities can induce nonlinear damping due to internal viscoelastic damping. From their theoretical results they conclude that this mechanism could account for at most only about 1 percent of the measured nonlinear dissipation. The other mechanism emerged from the investigations of Croy et al. [10] of the vibrations of nanoresonators made of suspended graphene sheets and actuated electrostatically, like those investigated experimentally by Eichler et al. [3]. They modeled the graphene sheet using an elastic continuum model and derived a pair of coupled nonlinear partial differential equations (PDEs) as the equations of motion for the out-of-plane and in-plane vibrations. The dissipation of mechanical energy occurred because the out-of-plane motion was converted into in-plane phonons in the suspended region resulting in phonon emission (acoustic energy loss) through the supporting region. In this resonator this region was located at the surface of the substrate and extended over the area along which the graphene sheet was glued to this surface due to the van der Waals force [10]. In order to calculate the acoustic emissions through the support region it was necessary to discretize the Fourier transformed response of the substrate to the applied stress and obtain a linear system of equations for the in-plane response function that had to be solved numerically. Unfortunately, the parameters characterizing the linear and nonlinear dissipation cannot be easily obtained using this approach, and no analytical expressions for such parameters as a function of geometrical and material parameters have been derived. The mechanism of acoustic energy loss through the attachment region is, in fact, a well-known damping mechanism [11,12] in more conventional micro- and nanoresonators like those investigated by Zaitsev [4] and Imboden [5]. Such resonators consist of a suspended beam clamped to the substrate or other anchoring structures, and are one of the most commonly investigated micro- and nanoresonators for practical applications. However, only the acoustic emissions resulting in linear damping have been investigated theoretically for such resonators resorting either to purely analytical [11,12] or numerical methods, such as the finite element method [13], or a combination of both [14]. In this work we extend the conventional purely analytical treatment of clamping loss in doubly clamped beam resonators and include a source of acoustic energy loss that leads to nonlinear damping. Using this approach we can derive analytical expressions for the amplitude dependent quality factor and the nonlinear damping coefficient of doubly clamped beam resonators. The paper is organized as follows: in Section 2 we present the origin of the acoustic energy loss that results in the nonlinear damping and discuss the assumptions and approximations that have to be made in order to calculate the quality factor and the nonlinear damping parameter. In Section 3 we show how to obtain the approximate solutions of the nonlinear PDE required to calculate the relevant physical quantities from a corresponding reduced order model of the beam. In Sections 4 and 5 we derive analytical expressions for the amplitude dependent quality factor and the nonlinear damping coefficient, respectively. In Section 6 we analyze the relevance of the proposed nonlinear damping mechanism for realistic MEMS/NEMS resonators. In Section 7 we compare the predictions based on our analytical expressions with the experimental results obtained by Zaitsev et al. and conclude in Section 8.
2. Nonlinear damping mechanism When operating in the linear regime the dynamics of thin beams is well described by the Euler–Bernoulli beam equation. However, as the vibration amplitude increases the equation must be modified to include the effects of the elongation that become significant. If the amplitude of vibration is large, but sufficiently small that linear stress–strain law can be assumed, and that the effects of rotary inertia and transverse shear, as well as, nonlinearities due to inertia and curvature can be neglected, the resulting equation of motion is the following nonlinear PDE [1], ! Z EA l 0 2 00 ⁗ € þcw _ ¼ F ðx; t Þ; w dx w þ ρAw (1) EIw N 0 þ 2l 0
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
257
where w ¼ wðx; tÞ denotes the transversal displacement of the beam located between x¼0 and x ¼l. The overdots and primes represent derivatives with respect to time, t, and space, x, respectively, E denotes the Young modulus, I the geometric moment of inertia, A the cross-sectional area, and ρ the beam density. We have also considered that the beam can be subject to an axial load N0 due to residual stress. Included in the equation is a viscous damping term, which accounts for any existing linear damping mechanisms, and the forcing along the beam F ðx; tÞ. The nonlinear damping mechanism that we investigate originates from the geometric nonlinearity term Z EA l 0 2 00 00 w dx w ¼ Nðt Þw : (2) 2l 0 In this equation we have evidenced the time-varying axial force N(t) induced along the beam as it vibrates transversally. This time-varying force is transmitted to the point of attachment of the beam and is directly or indirectly responsible for the production of acoustic waves that carry away a fraction of the beam mechanical energy. In order to understand how this normal force can result into acoustic emissions we have to consider that there are typically two classes of anchoring configurations for suspended beam micro- and nanoresonators: we can have the beams anchored directly to the substrate (within the substrate) or to anchoring structures built on top of it (see Ref. [12] for more details). If the beam is attached to an anchoring structure built on top of the substrate the normal force N(t) is applied at some point of the structure pulling it and inducing, as reaction forces, time-varying shear stresses and bending moments along the region of the surface of the substrate where the anchoring structure is founded. These shear and bending moments then produce the emission of acoustic waves into the substrate. If the beam is clamped directly to the substrate the force N(t) results in a time-varying normal stress at the attachment point which induces the emission of acoustic waves directly into the substrate. Since in the case of anchoring structures built on top of the substrate the determination of the reaction forces and moments induced at the substrate can be rather complex, we restrict our analysis to beams attached directly to the substrate, a configuration that is commonly found in practical devices (examples can be found in Refs. [4,14]) and is much more accessible to a fully analytical investigation. For the configuration of a beam attached directly to the substrate it has been customary to model the beam as attached to either a semi-space or a semi-infinite thin plate [11,12,15] depending upon the details of the attachment region being investigated. Here we are going to investigate the resulting nonlinear damping considering both models for the attachment region, which are depicted in Fig. 1. Before we proceed we have to note that besides the normal force N(t) induced by the geometric nonlinearity, the in-plane and out-of-plane transverse vibrations of the beam induce shear forces and bending moments at the attachment point [11] that result in the well-known linear damping mechanism of clamping loss. While the combined effects of shear and bending may have to be considered in order to determine the total acoustic energy loss resulting in the linear damping [11], as in the case of beams supported by thin plates, the acoustic emissions due to N(t) can be analyzed independently of the emissions due to the other forces [16]. Assuming small nonlinearity, we can still characterize the damping of the beam vibration due to the acoustic emissions induced by N(t) by the corresponding quality factor defined as Q ¼ 2π
U U ¼ω ; ΔU Π
(3)
where U denotes the total mechanical energy of the beam vibrating at frequency ω, ΔU the energy dissipated per oscillation cycle, and Π the average emitted acoustic power. To calculate analytically these physical quantities we adopt the same assumptions made previously in Refs. [11,12,15,17]: (i) the dimensions of the beam are such that the wavelength of the emitted acoustic waves are much larger than the dimension of the attachment point, alloying the use of the point source approximation. (ii) In spite of the elastic vibrations at the attachment point, the vibrations of the beam are determined assuming perfectly rigid attachment, that means the beam satisfies the boundary conditions wð0; tÞ ¼ wðl; tÞ ¼ w0 ð0; tÞ ¼ w0 ðl; tÞ ¼ 0. (iii) The stresses required to satisfy the corresponding boundary conditions are those applied to the
Fig. 1. Schematic diagram showing one end of a beam attached to a semi-space (left) and a thin plate (right).
258
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
surface of the substrate in the region corresponding to the attachment point, and act as radiation sources of acoustic waves that are completely lost into the substrate or, more generally, into large surrounding structures. The calculation of U, ΔU and Π requires the knowledge of wðx; tÞ. Since we cannot obtain an exact solution of Eq. (1) we have to resort to an approximate method of calculation of wðx; tÞ as detailed in the next section.
3. Approximate solution to the PDE In order to choose an adequate method we note first that we are interested in calculating the dissipation and the quality factor for frequencies close to or at the resonant frequencies of the beam. Furthermore, for small nonlinearity, the expected resonant frequencies and modeshapes must be close to those obtained for the corresponding linear PDE, except under special circumstances, such as the occurrence of internal resonances [6], which we assume that are being circumvented with an adequate choice of the system parameters. In this case we can resort to the Garlekin method [1], which allow us, in a mathematically consistent manner, to reduce Eq. (1) to a single nonlinear ordinary differential equation (ODE) which contains adequate approximate expressions for the physical parameters required to calculate U and provides a simple procedure to calculate N(t), required to calculate Π. The Garlekin method is widely used in the modeling of micro and nanoresonators, and has been used to reduce the original nonlinear PDE to a single nonlinear ODE which was then solved numerically [18] or by approximate analytical methods [1,2]. Applying this method, we have that close to the resonance corresponding to the nth eigenmode of the corresponding linear PDE we approximate w by wðx; tÞ ¼ ϕn ðxÞun ðtÞ, where ϕn ðxÞ is nth eigenmode of the linear and homogeneous PDE obtained from Eq. (1) by suppression of the geometric nonlinearity, viscous damping and forcing. More specifically, the ϕn's are solutions of the following time-independent equation: EI ⁗ N 0 00 ϕ ϕ ω2n ϕn ¼ 0; ρA n ρA n where ωn denotes the corresponding eigenfrequencies. The solutions to the above equation have the general form ϕn ðxÞ ¼ coshðκn x=lÞ cos ðκn x=lÞ þ χ n sinhðκn x=lÞ sin ðκn x=lÞ :
(4)
(5)
The parameters κn and χn are functions of the axial stress N0 and, in general, can only be obtained numerically [1]. However, in the limiting case N0 ¼ 0 the ϕn ðxÞ are the well-known modeshapes of the Euler–Bernoulli equation, for which we have κ n ¼ 4:7300; 7:8532; 10:996… and χ n ¼ 0:9825; 1:0008; 0:99996…. The other relevant limiting case is obtained when N 0 -1, that means, when the axial force is much larger than pffiffiffithe elastic restoring force. In this limit the beam behaves like a string and the modeshapes have the simple form ϕn ðxÞ ¼ 2 sin ðnπx=lÞ. At this point we have defined the approximate spatial dependence of wðx; tÞ and we are left with the task of determining the approximate dynamical behavior of the system. The dynamics is contained in the function un(t), whose corresponding equation of motion is the final product of the Galerkin method. To obtain this equation applying this method we first substitute wðx; tÞ ¼ ϕn ðxÞun ðtÞ in Eq. (1). Eq. (4) is used to simplify the resulting equation which is then multiplied by ϕn ðxÞ and integrated along the beam length, and further simplified resorting to the orthogonality of the modeshapes Z l Z 1 ½ϕn ðxÞ2 dx ¼ l ½ϕn ðyÞ2 dy ¼ l; (6) 0
0
where we have evidenced the orthonormality of the modeshapes in the renormalized variable y ¼ x=l. Considering that the resonator is being driven by an harmonic force of the form F ðx; tÞ ¼ F 0 cos ðωn tÞ, which is the case, for instance, in the experiments of Refs. [3–5] we thus obtain the following nonlinear ODE: mu€ n þ2b11 u_ n þ k1 un þk3 u3n ¼ F cos ðωn tÞ;
(7)
where the notation is the same adopted by Zaitsev et al. [4] in order to facilitate the later comparison with their experimental results. In this equation we have that m ¼ ρAl, b11 ¼ cl=2, the effective spring constant is given by N0 I n þ I n1 ; 3 3 l
(8)
ðI n1 Þ2 :
(9)
ϕn ðyÞ dy;
(10)
½ϕ0n ðyÞ2 dy;
(11)
k1 ¼ mω2n ¼
EI l
and the nonlinear spring constant reads k3 ¼
EA 2l
3
n
The force now is given by F ¼ F 0 lI 0 and we have defined I n0 ¼ I n1 ¼
Z
1 0
Z
1 0
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
Z I n3 ¼
1 0
ϕ⁗ n ðyÞϕn ðyÞ dy:
259
(12)
Eq. (7) has the form of the well-known Duffing equation [6] for a nonlinear oscillator, and has been extensively studied. For the purpose of calculating U, required to derive Q predicted by the proposed damping mechanism, as well as others quantities like N(t) and ΔU, we need to know the un(t) which is obtained as the solution of Eq. (7). In order to define an adequate approximation to the function un(t), we first note that in common situations of practical interest the amplitude of motion of the beam is sufficiently small that the resonator performs quasi-harmonic oscillations, as a consequence of the harmonic forcing. That means the amplitude of the terms n 41 in the general steady state solution un ðtÞ ¼
1 X
am cos ðmωn t þ ϕm Þ;
(13)
m¼0
are small compared to a1, the amplitude of the harmonic term (note that a0 is not relevant for our problem). In fact, using the harmonic balance method, the authors of Ref. [19] have evidenced that even in the strong nonlinear regime, where the resonant curve is highly distorted, the amplitude of the next relevant higher harmonic does not exceed approximately 10 percent of the leading harmonic term. He have further checked the amplitude of the next higher harmonic solving numerically Eq. (7) for a resonator operating in the weak nonlinear regime, when the forcing results in a slightly tilted frequency–response curve. We have found that the amplitude of the next higher harmonic does not exceed 1 percent of a1. Therefore, considering that the resonators are operating in the weak nonlinear regime it suffices to take the harmonic approximation un ðtÞ ¼ u0 cos ðωn tÞ in order to calculate U, ΔU, and N(t) because such quantities are proportional to, at least, the square of the amplitude. We note that we can ignore any phase in the definition of un(t) because we are interested only in time averaged quantities.
4. Quality factor 4.1. Q for the proposed damping mechanism The acoustic energy loss in the proposed nonlinear damping mechanism results from the time-varying normal force N(t) applied to the clamping region. Substituting the harmonic approximation wðx; tÞ ¼ u0 ϕn ðxÞ cos ðωn tÞ into the expression for N(t) [see Eq. (2)] and neglecting the contribution of a constant force term we obtain N ðt Þ ¼
EAI n1 4l
2
u20 cos ð2ωn t Þ:
(14)
It is important to note that the forcing is harmonic and has a frequency that is twice that of the vibrating beam. The resulting emitted acoustic power is going to depend upon details of the supporting region. As mentioned previously, here we consider two cases: the beam attached to semispaces and to semi-infinite thin plates having the same thickness of the beam. The first case is an approximation to large anchoring structures. The second case represents approximately the conditions found in many actual micro and nanoresonators, including those investigated in Refs. [4,5], where the fabrication process results in cantilevered supporting structures. For the case of anchoring in the semispaces we use the expression, derived by Miller and Pursey [20], for the acoustic power emitted by a point source vibrating normal to the surface of a semispace. For an harmonically varying force of the form FðtÞ ¼ f cos ðωtÞ they employed the admittance method obtaining pffiffiffiffiffiffiffiffiffiffi ϕ ðγÞ ρC 11 2 2 ω f (15) Πs ¼ 1 4π C 244
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where ρ is the density, γ ¼ C 11 =C 44 , with C 11 ¼ Eð1 νÞ=ð1 ν 2ν2 Þ and C 44 ¼ E=2ð1 þνÞ, where E corresponds to the Young modulus and ν the Poisson constant of the isotropic material of the beam and semispace. The function ϕq ðγÞ is defined as "Z pffiffiffiffiffiffiffiffiffiffiffiffiffi # 1 q p p2 1 dp ; (16) ϕq ðγ Þ ¼ Im F 0 ðp; γÞ 0 with F 0 ðp; γÞ ¼ ð2p2 γ 2 Þ2 4p2
qffiffiffiffiffiffiffiffiffiffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffi p2 1 p2 γ 2 :
(17)
Before we proceed, we note that ϕq ðγÞ has been evaluated differently in the literature as we have already pointed previously [9] for the cases q ¼1 and 3. In Ref. [9] explicit expressions for ϕ1 ðγÞ and ϕ3 ðγÞ have been obtained following the prescriptions of Miller and Pursey [20], and can be very easily generalized for arbitrary positive integer q.
260
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
The emitted acoustic power for a point source acting on the edge of a semi-infinite thin plate of thickness h was derived by Hao and Xu [12] and can be cast in the form 2
Πp ¼
2ϕ0 ðγÞ 1 þν f ω : π Eð1 νÞ h
(18)
The quality factor defined in Eq. (3) must be rewritten in terms of the average emitted powers Π s;p as Q ¼ ωU=ð2Π s;p Þ ¼ ωn U=Π s;p , where the factor of 2 in the denominator, which is a consequence of the energy being lost by both ends of the beam, cancels with the factor of two coming from the fact that ω ¼ 2ωn according to Eq. (14). We can obtain U from the maximum potential energy in terms of the parameters in Eq. (7) as u2 k3 (19) U ¼ 0 k1 þ u20 : 2 2 Using k1 ¼ mω2n and k3 defined in Eq. (9) in this expression for U, the resulting quality factor when the beam is attached to semi-spaces is pffiffiffi 5 4πC 244 ρl ωn πC 244 l pffiffiffiffiffiffiffiffiffiffiffi p ffiffiffiffiffiffiffiffi 2 þ ¼ Q sn;1 þ Q sn;2 ; (20) Q sn ¼ n 2 2 C 11 E A u0 ϕ1 ðγÞ ρ C 11 EA ωn ϕ1 ðγÞðI 1 Þ where we have left the final expression in terms of ωn for the sake of clarity, however, using Eq. (8) this frequency can be s rewritten in terms of the parameters of the beam. We can see that Qn can be separated in two terms. Q sn;1 is the amplitude s 2 dependent contribution, and goes as u0 . The other term contributing to Qn is amplitude independent. It results from the ratio between the energy stored in the nonlinear component of the spring constant to the energy lost through the nonlinear s damping mechanism, and determines the smallest value of Qn that is attained as the amplitude of vibration increases. We note, however, that as the amplitude increases higher order contributions to the effective elastic constant and the dissipated energy can become relevant and change the minimum quality factor before the beam suffers permanent deformation or rupture. Using the expression for Πp, Eq. (18), the quality factor Q ¼ ωn U=Π p for a beam clamped to thin plates of thickness h results to be 5
Q pn ¼
πρð1 νÞhl ω2n ϕ0 ðγÞðI n1 Þ2 Eð1 þ νÞA
u20
þ
πð1 νÞhl ¼ Q pn;1 þQ pn;2 ; 4ϕ0 ðγÞð1 þ νÞA
(21)
Again we have the amplitude dependent and the amplitude independent contributions to the total quality factor. The more s p noticeable difference between Qn and Qn is the dependence on the angular frequency ωn. Partially as a consequence of this s p distinct dependence on ωn, Qn and Qn have very distinct dependences on the geometrical parameters of the beam. The dependence on u0 2 observed in our expressions for Q sn;1 and Q pn;1 is the same one obtained in previous works that considered the addition of nonlinear damping terms of the form b31 u2n u_ n or b32 u_ 3n to the Duffing-like equations describing the nonlinear dynamics of resonant beams [2–5,10]. We note, however, that as no nonlinear damping was present in the original PDE, Eq. (1), such a term is also missing in the reduced order model, Eq. (7). Unfortunately, our approach does not allow us to determine which one of these two terms should be present in the reduced order model. Nevertheless, many aspects of the dynamics of the resonators do not depend on which of these two terms is present in the equation of motion, as already noticed in Refs. [2,4], because their effects are combined linearly in the same effective nonlinear parameters. In order to have a better understanding of the results we have obtained, in the next subsection we are going to add to Eq. (7) a term of the form b31 u2n u_ n þ b32 u_ 3n , and consider its effects on the quality factor of a resonator described by a modified equation of motion for the reduced order model, namely mu€ n þ ð2b11 þ b31 u2n þ b32 u_ 2n Þu_ n þk1 un þ k3 u3n ¼ F cos ðωn tÞ:
(22)
4.2. Q in the reduced order model We can easily calculate Q as a function of the parameters in Eq. (22) in the harmonic approximation. We already have U, given in Eq. (19), and the missing expression for ΔU is obtained from the work done by the dissipative forces along one oscillation cycle as b31 2 3 u0 þ b32 ω2n u20 : (23) ΔU ¼ πωn u20 2b11 þ 4 4 Therefore, the quality factor in terms of the parameters of the reduced order model is Q ¼ 2π
U ¼ ΔU
k1 þ ðk3 =2Þu20
: u2 ωn 2b11 þ b31 þ 3b32 ω2n 0 4
(24)
This result coincides with an analogous expression obtained in Ref. [10] using the slow envelope approximation. It is interesting to analyze some general features of this expression for Q. In the absence of nonlinear terms, we recover the usual expression for the
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
261
1.0
max
0.8 0.6 0.4 0.2
0.97
0.98
0.99
1
1.01
1.02
1.03
Fig. 2. Normalized frequency–response curves obtained by solving Eq. (26), with Q 0 ¼ 100, and Γ ¼ 10 8 for different forcing A2 ¼ 1; 2; 5, and 10 (Γ and A in arbitrary units). The broader curves, that correspond to a smaller quality factor, are obtained for larger forcing A that results in larger amplitudes.
amplitude independent Q. If, in addition, we had only the nonlinear spring term, Q would increase with the amplitude, because more mechanical energy is stored in the resonator for a given amplitude. When only b11 ¼ 0, the quality factor can be separated in a constant term and one that depends on the amplitude as u0 2 . Since this result is important for the comparison with Q predicted by the nonlinear acoustic energy loss mechanism, we present it explicitly Q ¼ Q 1 þQ 2 ¼
4k1 2
ωn bu0
þ
2k3 ; ωn b
(25)
where b ¼ b31 þ 3b32 ω2n , and we have separated Q in the amplitude dependent and amplitude independent terms, Q1 and Q2, respectively. The concept of an amplitude dependent quality factor was first introduced by Zhong-Kun and Luo [21] during the theoretical and experimental investigation of the quality factor of a torsion pendulum, where the effects of nonlinear damping were considered. More recently, Eichler et al. [3] observed experimentally a dependence with the driving amplitude of the quality factor of electrostatically driven nanotubes and graphene sheets. They developed a theoretical model which incorporated nonlinear damping and led to an expression for Q that coincides with Q1 presented in Eq. (25). The concept of an amplitude dependent Q is clarified in Fig. 2 where we have plotted the theoretical frequency–response curves for a system described by Eq. (22), obtained for different values of the forcing amplitude F and a fixed nonlinear damping coefficient b, and assuming k3 ¼ 0. The curves have been obtained by solving for u0 the equation u20 ¼
A2 2 2 ω ωn 3 k3 2 1 Q 0 1 þΓu20 u0 þ 8 k1 2 ωn
(26)
where A ¼ F=ð2k1 Þ, Q 0 1 ¼ 2b11 ωn =k1 , and Γ ¼ bωn =ð4k1 Þ. This equation for u0 was adapted to our notation from the result derived in Ref. [2] using secular perturbation theory for a nonlinear Duffing oscillator with nonlinear damping. It is interesting to note that the authors of Ref. [2] have derived an expression which evidenced the existence of an effective 1 amplitude dependent quality factor Q eff ¼ Q 0 1 þΓu20 , however, this interpretation was not given in that work.
5. Nonlinear damping coefficient s
p
From the above results for Qn and Qn we can derive the more fundamental nonlinear damping parameter b. This can be done very easily comparing the expressions for Q sn;2 and Q pn;2 with that of Q2 defined in Eq. (25) and the results are pffiffiffiffiffiffiffiffiffiffi ϕ ðγÞðI n1 Þ2 ρC 11 E2 A2 s ; (27) bn ¼ 1 4 πC 244 l for clamping on a semispace and p
bn ¼
4ϕ0 ðγÞðI n1 Þ2 ð1 þ νÞEA2 4
πð1 νÞhl ωn s
;
(28)
for clamping on a thin plate. It is interesting to note that bn does not depend explicitly on the frequency and, therefore, has s only a small dependence on the axial force N0 through the term ðI n1 Þ2 , while bn has a stronger dependence on N0 through the term ωn appearing in the denominator. The two b's also have very distinct dependencies on the geometrical parameters of the beam. In order to see that, let us s 2 4 consider a beam with rectangular cross section having width w and thickness h. In this case we have bn p h w2 =l . The
262
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265 p
particular dependence of bn on the geometrical parameters depends upon N0 and the polarity ofp the beam vibration. Let us ffiffiffiffiffiffiffiffiffiffiffiffi consider the two limiting cases for N0 we have considered previously. For N 0 ¼ 0 we have ωn ¼ EI=ρAðκ n =lÞ2 , with A ¼ hw 3 3 p 2 and I ¼ wh =12 in the case of out-of-plane vibrations, and I ¼ hw =12 for in-plane vibrations. Therefore, bp w2 =l np ffiffiffiffiffiffiffiffiffiffiffi ffi for outp 2 of-plane vibrations and bn p hw=l for in-plane vibrations. In the limit N 0 -1 we have that ωn ¼ ðnπ=lÞ N=ρA and, conpffiffiffiffiffiffiffiffiffiffiffi p p 3 2 sequently, bn p h w5 =l , regardless of the polarization. We can see that in any of the cases we have considered that bn and s bn have quite distinct dependencies on h; w and l. This fact evidences the need to model adequately the clamping region, since discrepancies of orders of magnitude can be obtained in the theoretical estimation of the nonlinear damping parameter and, consequently, on the quality factor. It is also evident that this acoustic nonlinear damping mechanism is more relevant for short beams with large cross sectional area, independently of the specific clamping.
6. Analysis of the results In this section we briefly analyze the potential relevance of the proposed nonlinear damping mechanism to realistic MEMS/NEMS resonators. Let us consider a MEMS resonator made from polysilicon (E ¼170 GPa, C 44 ¼ 69:7 GPa, C 11 ¼ 194:1 GPa, and ρ ¼ 2:3 103 kg m 3 ), which is a common material for MEMS/NEMS fabrication, with thickness h ¼ 1 μm and l ¼ 30 h, w ¼ 5 h, and Q 0 ¼ 103 . In Fig. 3(a) we plot the frequency response curve obtained solving Eq. (26) for the first mode (n ¼1), using k1 and k3 defined in Eqs. (8) and (9), respectively and Γ calculated using our analytical expressions for bs1 and bp1. The continuous lines represent the response in the absence of nonlinear dissipation [Γ ¼ 0 in Eq. (26)]. The dashed curves have been obtained assuming the resonator is clamped to a plate. We can see that there is a relevant change in the frequency response curve, particularly at larger amplitudes of vibration, which present a decrease of the amplitudes close to its maximum. The curves obtained considering the clamping to a semi-space were so close to those with no dissipation that they have been omitted. A more general result comes from the comparison of the additional dissipation that results from the nonlinear damping compared to the dissipation due to linear damping. While the linear damping Q0 results in a dissipation δ0 ¼ Q 0 1 that is independent of the amplitude of vibration, the nonlinear damping results in a dissipation that grows with the amplitude. The total dissipation is the sum of the linear and nonlinear contributions and results in the effective quality factor Q eff s defined after Eq. (26). For resonators made from polysilicon we obtain the following expressions for the dissipation from Q1 p and Q1 ! 3 3 1 1 0:131l 0:144l þ δs1 ¼ s ¼ ; (29) 2 Q1 u20 w h w and δp1
1 ¼ p¼ Q1
2
1:14l 1:26h l þ 2 w u0 w
!1 :
(30)
From these equations we can see that the dissipation is independent of the actual dimensions and amplitude of vibration of the resonators, and once the relative dimensions are fixed, the dissipation becomes a function only of the normalized p displacement u0 =h. For this reason, in Fig. 3(b) we plot the nonlinear dissipation δnl ¼ δs1 or δ1 as a function of u0 =h for two
40
u0
30 20 10 0
0.996
0.998
1
1.002
1.004
10
3
10
4
10
5
10
6
0.02
0.04
0.06
u0
0.08
0.10
Fig. 3. (a) Frequency–response curves for a realistic MEMS resonator made from polysilicon with dimensions h ¼ 1 μm, l ¼ 30 μm, and w ¼ 5 μm oscillating in its first mode n¼1 and having an intrinsic quality factor Q 0 ¼ 103 . The lower, middle, and higher curves have been obtained for a forcing F ¼ 2 10 8 ; 3 10 8 , and 4:5 10 8 N, respectively. Continuous lines represent the response in the absence of nonlinear damping (b ¼0) and the p dashed curves have been obtained considering the nonlinear damping assuming clamping to a plate ðb ¼ b1 Þ. (b) Nonlinear dissipation δnl ¼ 1=Q nl as a function of the normalized amplitude of vibration u0 =h in the first mode (n¼ 1) for a resonator clamped to a plate (continuous) and to a semi-space (dashed). The black and gray curves have been obtained for resonators with relative dimensions l ¼ 10 h, and w ¼ 3 h and l ¼ 30 h, and w ¼ 5 h, respectively.
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
263
representative sets of relative dimensions, namely, l ¼ 30 h, and w ¼ 5 h and l ¼ 10 h, and w ¼ 3 h. The aspect ratio l=h ¼ 10 assumed in this last case is approximately the smallest value for which results obtained using the Euler–Bernoulli beam theory can be trusted [22]. The range of u0 =h was chosen as representative of a region where the amplitudes are sufficiently large that nonlinear effects are relevant, but the resonators are not in a strong nonlinear regime. For example, in Fig. 3(a) the maximum displacement for each curve ranges from u0 =h 0:015 up to u0 =h 0:035. The last maximum value corresponds to the amplitude close to the onset of bistability for the resonator considered for Fig. 3(a) with Q 0 ¼ 103 . This amplitude agrees well with the predicted theoretical value of the critical amplitude for the onset of bistability which, ignoring the small effect of nonlinear damping and using our notation, is given by [2] 34 sffiffiffiffiffiffiffiffi 4 k1 c u0 ¼ : (31) 3 Qk3 For the resonator considered for the results presented in Fig. 3(a) it results to be uc0 ¼ 0:029 h. According to this equation, if Q0 were of the order of 102 the critical amplitude would increase to approximately uc0 ¼ 0:092 for the same resonator, a prediction that is slightly smaller than the results we have obtained by solving Eq. (26) numerically, and that was used to set the maximum value of u0 =h in Fig. 3(b). Therefore we can compare the dissipation presented in Fig. 3(b) with the linear dissipation of a given device and estimate the potential relevance of the nonlinear damping mechanism. We can, for example, see that for u0 =h 0:035, at the onset of bistability, a nonlinear dissipation of δp1 10 4 is obtained for a resonator with l ¼ 30 h and w ¼ 5 h, which amounts to about 10 percent of the dissipation δ0 ¼ 10 3 of a resonator with Q 0 ¼ 103 . For the shorter resonator the dissipation increases to δp1 2:5 10 4 , and amounts to 25 percent of the linear dissipation. We can also see that the dissipation for beams clamped to semi-spaces is much smaller, and would amount to no more than a few percent of the linear dissipation considered above. Similar results hold for the nonlinear dissipation predicted at the onset o bistability for resonators with Q 0 ¼ 102 as can be seen from the curves at u0 =h 0:1.
7. Comparison with experiment Since our theoretical investigation has been motivated by recent measurements of nonlinear damping in micro and nanoresonators, it is interesting to evaluate the relevance of this nonlinear energy loss channel to the observed large nonlinear damping. For that purpose, we compare our predictions with the experimental results of Zaitsev et al. [4]. In this work detailed experimental results are presented for a rectangular cross section suspended beam microresonator made of the Pd0.15Au0.85 alloy, and dimensions l¼125 μm, w¼ 1 μm, and h¼0.3 μm. The microresonator was excited in its first mode, whose resonant frequency was measured to be f 1 ¼ 885:5 kHz, and performed in-plane vibrations (parallel to the dimension w). They have extracted experimentally several physical parameters that are related to the parameters defined in Eq. (22). One particularly relevant parameter is p, that is related to the coefficients in Eq. (22) by p¼
bωn 2 1 ¼ ; 3k3 3 Q 2
(32)
where the last equality results from the comparison of the definition of p with the expression for Q given in Eq. (25). Considering that the wet etching thins the support region close to the attachment region [23], we approximate the support structure by a thin plate in order to derive an analytical expression for p. Substituting into the previous equation the expression for Q pn;2 we have pt ¼
8ϕ0 ðγÞð1 þ νÞ w : 3πð1 νÞ l
(33)
In order to calculate pt we are going to consider that ν ¼ 0:43. This value is the weighted average (weighted by the atomic fractions) of the corresponding values for ν of the bulk Pd ðν ¼ 0:39Þ and that of the bulk Au ðν ¼ 0:44Þ. For ν ¼ 0:43 we have that γ ¼ 2:854 and ϕ0 ð2:854Þ ¼ 0:10996, and Eq. (33) predicts the value pt ¼ 0:002. This is only a small fraction of the value found experimentally [4] p¼ 0.109. We have to note that this value for pt was obtained assuming that the elastic properties of the beam are those of an homogeneous beam, having the properties of the bulk PdAu alloy. However, the fabrication process of the microbeams may result in structures with effective properties that differ significantly from those of bulk materials. For instance, the deposition process of thin films may result in PdAu alloys with densities and elastic properties that differ significantly from those of the bulk material [24]. Furthermore, the electron beam lithography used in the microfabrication process resulted in the deposition of an amorphous carbon layer on the beam surface [23]. This material typically has a density much smaller than that of PdAu alloys but a much larger Young modulus. The theoretical modeling of the resulting mechanical system can be rather complex, however, it is usually the case that effective mechanical properties can be defined for the system such as an effective flexural rigidity and effective density [1]. Therefore, a more precise comparison between theory and experiment would be obtained using the relevant effective physical properties of the microbeam. Fortunately, the detailed experimental analysis of Zaitsev et al. [4] provides us with the possibility to determine the experimental value of the effective nonlinear spring constant k3, which is a key ingredient in the derivation of pt.
264
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
Considering the value for the parameter α3 ¼ k3 =m estimated by Zaitsev et al. [4], α3 2 1023 m 2 s 2 , we can extract an approximate experimental value for k3. We estimate the mass of the PdAu alloy as was done by Zaitsev et al. [4] and find m ¼ 6:8 10 13 kg. This value is obtained considering the dimensions of the beam given above, and assuming a density for the PdAu alloy of ρ ¼ 18:2 103 kg m 3 , which is the weighted average of the Au density of 19:3 103 kg m 3 and that of e Pd, 12:0 103 kg m 3 . We thus obtain the experimental value k3 ¼ 1:36 1011 N m 3 . This experimental value can be t compared with the expected theoretical value k3, calculated using the definition presented after Eq. (7). For that purpose we use the Young modulus E ¼95 GPa, estimated from the experimental results [25] for E of PdxAu1 x alloys with Pd atomic fraction x. We have also considered the fact that the beam dynamics is dominated by the axial stress since we have estimated theoretically that the contribution of the elastic deformation term to k1 represents only about 2 percent of its total t e t experimental value k1 ¼ mω2 22 N m 1 . Therefore, calculating k3 in the limit N 0 -1 we obtain k3 ¼ 5 1011 N m 3 , which is considerably larger than the experimental value. e The expression for pe, the theoretical value estimated using the experimental value of k3, is 2
pe ¼
4π 3 ϕ0 ðγÞEð1 þνÞhw e 4 3k3 ð1 νÞl
:
(34)
e
Considering the value for k3 given above, we obtain pe ¼ 0:01. This is significantly larger than the purely theoretical value and is a better estimate of the actual value of this parameter as a result of the proposed nonlinear damping mechanism. We can also evaluate the more fundamental parameter b and compare with its experimental value. The experimental b is obtained from the parameter γ 3 ¼ b=m estimated by Zaitsev et al. [4] to be γ 3 ¼ 1 1016 m 2 s 1 . Therefore, we obtain e b ¼ 6:8 103 kg m 2 s 1 . The theoretical value, obtained from Eq. (28), is independent of k3 and, because we use as input the experimental value of the frequency, can be considered as the more realistic estimate of b. We thus have t b ¼ 0:72 103 kg m 2 s 1 . These theoretical results for p and b indicate that the acoustic nonlinear damping mechanism can be responsible for about 10 percent of the observed nonlinear damping. While this contribution can change if we perform a more detailed analysis of the experimental setup, it seems that, in this particular microresonator, the acoustic nonlinear damping can be relevant but it is not the prevailing nonlinear energy loss channel and other unknown nonlinear damping mechanisms should be present that must still be unveiled.
8. Conclusion While we have derived analytical expressions for Q and b for two specific models of the clamping region our approach can be applied if more specific models of the clamping region are required for a precise calculation of the nonlinear damping. The fundamental idea in our approach is to recognize that a beam vibrating with large amplitude develops a timevarying axial stress that is transmitted to the clamping region. For more general beams, for instance, with a variable crosssection, N(t) can be obtained analytically for a great variety of beam designs, or it could be determined resorting to well established computational methods, such as the finite element method. The calculation of the emitted acoustical power due to the resulting N(t) for more general configurations of the clamping region is certainly more complicated. However, computational approaches have been developed that can use both the analytical or fully computational results for N(t) as input for the calculation of the dissipated acoustical energy [13,16]. We believe that the use of the fundamental idea presented here in conjunction with other computational methods can be helpful in determining the induced nonlinear acoustical energy loss for more complex microscopic or macroscopic structures. It is also important to note that the damping resulting from the proposed nonlinear mechanism must, in general, be combined with other linear or nonlinear damping mechanisms relevant for a given suspended beam resonator. Having the results for the corresponding quality factors for each damping mechanism, the total loss (Q 1) is simply the sum of the loss predicted for each damping mechanism [7]. To conclude we note that the acoustic nonlinear damping mechanism proposed here has a strong dependence on the s p geometrical parameters of the doubly clamped resonators as can be seen, for instance, from the expressions for bn and bn. Furthermore, this dependence changes considerably depending on the geometry of the supporting region. As a consequence, the relevance of this mechanism in a particular setup is going to depend on the geometry and on details of the clamping region. This situation is analogous to that observed for several linear damping mechanisms [9,7,11,12], whose individual contribution to the total linear damping depends on the detailed geometry of the resonator. In the particular case of the proposed mechanism the expressions for Q and b indicate that it can be the prevailing source of nonlinear damping for short resonators with large cross sectional area, and should be considered in the analysis of the nonlinear dynamics of such devices.
Acknowledgments The author is thankful to S. Zaitsev for providing valuable information about the experiment reported in Ref. [4].
A. Gusso / Journal of Sound and Vibration 372 (2016) 255–265
265
References [1] M.I. Younis, MEMS Linear and Nonlinear Statics and Dynamics, Springer, New York, USA, 2011. [2] R. Lifshitz, M.C. Cross, Nonlinear dynamics of nanomechanical resonators, in: G. Radons, B. Rumpf, H.G. Schuster (Eds.), Nonlinear Dynamics of Nanosystems, Wiley-VCH Verlag, Weinheim, Germany, 2010, pp. 221–266. [3] A. Eichler, J. Moser, J. Chaste, M. Zdrojek, I. Wilson-Rae, A. Bachtold, Nonlinear damping in mechanical resonators made from carbon nanotubes and graphene, Nature Nanotechnology 6 (2011) 339–342. [4] S. Zaitsev, O. Shtempluck, E. Buks, O. Gottlieb, Nonlinear damping in a micromechanical oscillator, Nonlinear Dynamics 67 (2012) 859–883. [5] M. Imboden, O.A. Williams, P. Mohanty, Nonlinear dissipation in diamond nanoelectromechanical resonators, Applied Physics Letters 102 (2013) 103502; M. Imboden, O.A. Williams, P. Mohanty, Observation of nonlinear dissipation in piezoresistive diamond nanomechanical resonators by heterodyne down-mixing, Nano Letters 13 (2013) 4014–4019. [6] A.H. Nayfeh, D.T. Mook, Nonlinear Oscillations, Wiley, New York, USA, 1979. [7] M. Imboden, P. Mohanty, Dissipation in nanoelectromechanical systems, Physics Reports 534 (2014) 89–146. [8] A. Gusso, Phenomenological modeling of long range noncontact friction in micro- and nanoresonators, Journal of Applied Physics 110 (2011) 064512. [9] A. Gusso, Energy loss mechanism for suspended micro- and nanoresonators due to the Casimir force, Physical Review B 81 (2010) 035425; A. Gusso, Acoustic electromechanical energy loss mechanism for suspended micro- and nanoelectromechanical resonators, Applied Physics Letters 96 (2010) 193504. [10] A. Croy, D. Midtvedt, A. Isacsson, J.M. Kinaret, Nonlinear damping in graphene resonators, Physical Review B 86 (2012) 235435. [11] J.A. Judge, D.M. Photiadis, J.F. Vignola, B.H. Houston, J. Jarzynski, Attachment loss of micromechanical and nanomechanical resonators in the limits of thick and thin support structures, Journal of Applied Physics 101 (2007) 013521. [12] Z. Hao, Y. Xu, Vibration displacement on substrate due to time-harmonic stress sources from a micromechanical resonator, Journal of Sound and Vibration 322 (2009) 196–215. [13] D.S. Bindel, S. Govindjee, Elastic PMLs for resonator anchor loss simulation, International Journal for Numerical Methods in Engineering 64 (2005) 789–818. [14] G.D. Cole, I. Wilson-Rae, K. Werbach, M.R. Vanner, M. Aspelmeyer, Phonon-tunneling dissipation in mechanical resonators, Nature Communications 231 (2) (2011) 1–8. [15] M.C. Cross, R. Lifshitz, Elastic wave transmission at an abrupt junction in a thin plate with applications to heat transport and vibrations in mesoscopic systems, Physical Review B 64 (2001) 085324. [16] B. Chouvion, S. Mcwilliam, A.A. Popov, C. Fox, Review and comparison of different support loss models for microelectromechanical systems resonators undergoing in-plane vibration, Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science 226 (2011) 283–295. [17] D.M. Photiadis, J.A. Judge, Attachment losses of high Q oscillators, Applied Physics Letters 85 (3) (2004) 482–484. [18] T.D. Amorim, W.G. Dantas, A. Gusso, Analysis of the chaotic regime of MEMS/NEMS fixed-fixed beam resonators using an improved 1DOF model, Nonlinear Dynamics 79 (2015) 967–981. [19] Z.K. Peng, Z.Q. Lang, S.A. Billings, G.R. Tomlinson, Comparisons between harmonic balance and nonlinear output frequency response function in nonlinear system analysis, Journal of Sound and Vibration 311 (2008) 56–73. [20] G.F. Miller, H. Pursey, The field and radiation impedance of mechanical radiators on the free surface of a semi-infinite isotropic solid, Proceedings of the Royal Society of London Series A 223 (1954) 521–541; G.F. Miller, H. Pursey, On the partition of energy between elastic waves in a semi-infinite solid, Proceedings of the Royal Society of London Series A 233 (1955) 55–69. [21] Z.-K. Hu, J. Luo, Amplitude dependence of quality factor of the torsion pendulum, Physics Letters A 268 (2000) 255–259. [22] K.F. Graff, Wave Motion in Elastic Solids, Dover, New York, USA, 1991. [23] S. Zaitsev, Forced and self-excited oscillations in nonlinear micromechanical resonators, PhD Thesis, Technion, Israel, 2011. [24] S. Yalcin, R. Avci, Characterization of PdAu thin films on oxidized silicon wafers: interdiffusion and reaction, Applied Surface Science 214 (2003) 319–337. [25] J. Lohmiller, Investigation of deformation mechanisms in nanocrystalline metals and alloys by in situ synchrotron X-ray diffraction, Scientific Publishing, Karlsruher, Germany, 2012.