Tailoring the nonlinear response of MEMS resonators using shape optimization

Resonant micro-electromechanical-systems (MEMS) have grown in popularity over the past several decades due to their numerous applications, which include stable frequency generation and timing,^1,2 robust mass and angular rate sensing,^3,4 precise signal filtering,^5–7 and energy harvesting.⁸ While the majority of MEMS resonators are designed to operate in their linear dynamic range, research has shown that utilizing nonlinearity in resonant MEMS sensors can significantly improve the performance of these systems in some of these applications.^1,9 Nonlinear stiffness effects in MEMS resonators commonly arise from multiple sources, such as finite deformations that lead to nonlinear strain-displacement relationships,¹⁰ and the nonlinear nature of the electrostatic forces in capacitive MEMS.^4–7,11,12

The ability to manipulate the nonlinearity in resonant MEMS and NEMS (nano-EMS) allows one to design devices with greater control over their dynamical responses.¹³ For example, it has been shown that one can relax the constraints in mode mismatch in MEMS gyroscopes due to fabrication errors by independently tuning the linear and/or cubic stiffness coefficients of the capacitive drive,^4,14,15 and one can achieve optimal drive conditions for micro-resonators in order to enhance their dynamic range.^16,17 Although a large number of research efforts are focused on the understanding of nonlinear dynamics in MEMS and their applications to sensors and actuators, only a few have focused on the systematic optimization of nonlinearities in MEMS to achieve desirable outputs. Ye et al. demonstrated optimization of interdigitated comb finger actuators to achieve linear, quadratic, and cubic driving force profiles,¹⁸ and, more recently, Guo designed quadratic shaped comb fingers to achieve large displacement parametric resonance, overcoming the limited movement in non-interdigitated comb drive.^6,12 These studies used electrostatic effects to achieve desired nonlinear forces. Dou et al. investigated shape optimization in MEMS resonators using a gradient-based method in order to achieve targeted mechanical nonlinear coefficients.^19,20 The numerical results in Refs. 19 and 21, when applied to a clamped-clamped (bridge) beam, showed that the removal/addition of material (by changing the cross-sectional area) from the areas where the slope of the resonator mode shape is maximal results in a decrease/increase of the geometric contribution to the resonator cubic nonlinearity. This effect becomes more pronounced as the difference in the maximum and minimum beam thickness increases. The approach allows one to achieve combinations of frequencies and nonlinearities not obtainable with uniform beams.

The present experimental demonstration is based on results from the shape optimization approach of Dou et al.¹⁹ applied to altering the nonlinearity in a set of clamped-clamped microbeams. Beams with geometries corresponding to different stages in the optimization process were selected and fabricated in single-crystal silicon using standard SOI (silicon-on-insulator) processing techniques with DRIE (deep reactive ion etching) and released with HF vapor. Figure 1 shows SEM images and COMSOL models of three of representative designs of the clamped-clamped beams under study: a beam with uniform thickness (Beam_initial) and beams designed to minimize/maximize the cubic nonlinearity (Beam_min/Beam_max) under certain thickness constraints. The resonant frequencies of the beams range from 80 to 200 kHz.

In this work, we demonstrate manipulation of the Duffing coefficient of the fundamental mode with one variable parameter, the beam thickness, which can vary along the length, while maintaining all other parameters constant. Specifically, the beam lengths are all $L = 500 μ m$ (length) by $w = 20 μ m$ (width, dimension perpendicular to the beam deflection), and the transverse thickness, h, in the direction of the deflection, is allowed to take on values between 2 and $6 μ m$⁠, according to the optimization scheme described by Dou et al.¹⁹ As a result, the resonant frequency of the beams varies. While future work will allow for more practical constraints on the mode frequency, this work serves as initial demonstration that simple shape alterations can lead to significant changes in the Duffing nonlinearity of MEMS resonators.

The dynamic response of a clamped-clamped microbeam performing flexural vibrations in its fundamental mode when driven by a harmonic field can be modeled by the standard Duffing equation,^20,21 as follows:

where x is the modal amplitude, ω₀ is the natural frequency of the beam, Q is the quality factor (⁠ $Q = ω 0 / 2 Γ$ and Γ being the linear damping rate), γ is the coefficient of the cubic nonlinear stiffness term (Duffing nonlinearity), and g and ω are the drive amplitude and frequency. For a uniform doubly clamped thin elastic beam, the linear natural frequency of the fundamental mode is given by $ω 0 2 = 4 π 4 E h 2 / ( 9 ρ L 4 )$⁠, where E is Young's modulus and ρ is the material density. The nonlinear effect results from axial stretching due to bending as the beam is deflected, and the corresponding Duffing coefficient is given by $γ = 3 / ( 4 h 2 )$⁠.^22,23 For a beam of a given length, the ability to independently tune the natural frequency, which scales with h, and the nonlinearity, which scales with $h − 2$⁠, requires the use of non-uniform beams. One can systematically alter the Duffing nonlinearity of such non-uniform MEMS beams using a shape optimization process such as that described by Dou et al.¹⁹ 

In free vibration, the relationship between response amplitude (a_p) and the response frequency (ω_p) can be approximated as $ω p ( a p ) = ω 0 ( 1 + 3 8 γ a p 2 )$⁠, which provides the so-called Duffing backbone curve^22,23 that characterizes the system nonlinearity in the form of amplitude-frequency dependence. Thus, in order to assess the results of the proposed method on the system Duffing nonlinearity, we use the effective Duffing coefficient, $γ e f f = 3 γ / 8$⁠, as the objective function to be manipulated. When the system is subjected to harmonic drive, the critical vibration amplitude (a_nl) at which nonlinear effects result in a qualitative change of the resonator response, specifically, at which the frequency response becomes a multi-valued function of the drive frequency, is given by $a n l = ( 64 27 ) 1 / 4 1 Q γ$⁠, thus demonstrating that γ directly affects the vibration amplitude at which nonlinear effects come into play.

Clearly, one can alter the resonator critical vibration amplitude by changing the Duffing nonlinearity γ, and this can be achieved by simply changing the beam thickness h. In contrast, our goal is to systematically change the nonlinear coefficient by allowing the thickness to vary along the beam length, as described by Dou et al.¹⁹ Such an approach, when generalized to allow variation in the length and width, will provide an approach for manipulating the linear dynamic range for a given frequency and packaging space.

In the experiment, the beams were externally actuated by a shear piezoelectric stack actuator in order to achieve in-plane excitation via base acceleration. The in-plane motion of each beam was detected using a Laser Doppler Vibrometer (LDV) by tilting the microbeams at a $45 °$ angle with the laser focused on the sidewalls. The actual in-plane velocity of the microbeam can be calculated from simple geometry. This arrangement allowed characterization of the mechanical nonlinearity isolated from electrostatic effects. Amplitude-frequency responses for each beam were obtained at multiple drive amplitudes in atmospheric pressure from a spectrum analyzer. The schematic in Fig. 2 depicts the experimental setup.

The amplitude-frequency responses for the three representative beams shown in Fig. 1 are plotted in Fig. 3. As the driving amplitude increases, the amplitude-frequency response of each microbeam exhibits hardening behavior, as expected for a positive Duffing nonlinearity. Once the system response enters its bistable regime, since dissipation is small, the peaks of the steady-state response amplitude closely follow the backbone curve, as described by $ω p ( a p )$⁠. Using this result, we extracted the peak amplitudes from Fig. 3 and plotted the frequency as a function of the peak amplitude (a_p). Curve fitting is applied to the data to match the form $ω ( a p ) = ω 0 ( 1 + γ e f f a p 2 )$⁠, from which one can obtain estimates for the natural frequency and cubic nonlinearity. After these coefficients are determined, we normalized the frequency of each beam design by its natural frequency for direct comparison of the three representative beams (Fig. 4). We further compared the experimental results to the numerically computed coefficients, obtained using the techniques described in Ref. 19. Table I summarizes the comparisons between the beam designs, including those intermediate to the initial and optimized designs. Table I shows a comparison of the values of γ_eff for beams at different stages of the optimization iteration process with the corresponding Duffing nonlinearity of the nominal (uniform) beam, $γ e f f , 0$⁠. These characterization results show that the resonator shape optimization results in the natural frequency varying from its nominal, uniform beam, value. They also demonstrate that the designs that minimize the nonlinear coefficient reduce the Duffing nonlinearity by a factor of 2.6 when compared to the initial uniform beam, which corresponds to a factor 1.27 increase in the critical vibration amplitude at which the resonator becomes nonlinear, when accounting for the differences in quality factor among the two beams.²² Similarly, the beam design resulting in the maximum nonlinear coefficient achieves a 3.3 times increase in the Duffing nonlinearity as compared to the uniform beam. Table I shows that the experimental results for minimizing beam nonlinearity match computational results quite well, while those for maximizing the beam nonlinearity are off by a factor of two or more. Note that simply changing the thickness of the uniform beam to the maximum/minimum allowable thickness (⁠ $6 μ m$/ $2 μ m$⁠) results in factors of 2.25 (reduction)/4.0 (increase) compared to the uniform beam of thickness $4 μ m$⁠. While these non-uniform beams are not significantly different in this regard, these changes are achieved with less significant effects on the natural frequency. Specifically, a uniform beam of the type under consideration with thickness $6 μ m$ (⁠ $2 μ m$⁠) has a natural frequency of about 202 kHz (67 kHz), whereas the corresponding non-uniform beams have frequencies of 192 kHz (80 kHz), demonstrating some level of independent control over frequency and nonlinearity for beams of fixed length.

While the proposed designs have variable cross-sectional area, we account for the overall effect of h by considering the average thickness along the beam for the optimized designs. To illustrate the applicability of shape optimization methods for altering the nonlinearity, Fig. 5 shows a comparison of the nonlinear coefficients for the shaped beams with uniform beams with the same average thickness. These results show that to achieve the same amount of change in nonlinearity, the optimized beams required less average thickness changes when compared with uniform beams, thus allowing one to achieve a given nonlinearity with less shift in the natural frequency. We also performed numerical computations on the same beam geometries using the techniques described by Dou et al.¹⁹ for comparison. Comparison of the experimental results with the numerical predictions (Table I) shows that the minimized γ_eff agrees with the computed value within 8%. A larger discrepancy occurs in the maximization of γ_eff; however, the trend qualitatively agrees with the numerical results. We believe that the difference between numerical and experimental results for the beams designed to maximize γ_eff is due to the abrupt changes in the beam thickness profile, which causes errors in the numerical algorithm that utilizes simple beam elements. In contrast, beams designed to minimize γ_eff result in relatively small errors due to their relatively smooth beam profiles.

The ability to adjust the Duffing nonlinearity through design is especially useful in increasing the linear dynamic range of resonators, which plays a pivotal role in reducing phase noise in MEMS oscillators and for increasing the signal-to-noise ratio in resonant sensors.²⁴ The experimental results shown in this work are in good agreement with numerical predictions¹⁹ for the important case of minimizing the nonlinearity, and they demonstrate that one can utilize shape optimization methods for adjusting the resonator nonlinearity in a well-controlled manner. For more practical designs, it is desirable to include more variable parameters and constraints in the optimization process, for example, so that one can maintain the same resonant frequency or satisfy additional dimension constraints. Another important consideration is optimization for both mechanical and electrostatic effects and their relative contributions to the resonator parameters, such as stiffness coefficients and motional impedance. The latter is of significant importance for applications requiring low power consumption. These are possible and we are currently working on systematic designs of resonators with more freedom in the tuning of these parameters. Finally, the results presented provide strong confidence that nonlinear shape optimization methods can be applied in the design of MEMS resonators with more complicated geometries, with multi-physics effects, e.g., piezo-electric resonators, and for different nonlinear parameters, including difficult-to-control modal coupling coefficients.¹⁹ Additionally, it is possible to combine shape optimization with topology optimization algorithms to achieve the full potential of finite element modeling with considerations of fabrication constraints for targeted dynamical responses.

This work was supported by NSF Grant Nos. 1561829, 1561934, 1234067, 1234645, 442550-22324, and 442550-22031, ERC Starting Grant No. 279529 (INNODYN).