Modeling and Simulating MEMS Devices Using Finite Element Analysis

Introduction

In this paper, we discuss how Finite Element Analysis (FEA) can be used to simulate Micro Electro Mechanical Systems (MEMS). The mechanical behavior of many MEMS devices is a direct result of non-mechanical physical phenomena, such as electrostatics, electromagnetics or heat transfer. This coupling of mechanical behavior with another physical effect, combined with considerations associated with their small dimensions, poses significant challenges to the analysis of MEMS. These challenges can be met by FEA software capable of performing multiphysics analysis with no limitations on the dimensions of the objects in a model. With such software at their disposal, engineers can use the results of their simulations to design MEMS devices.

In this paper, we apply FEA to the analysis of a MEMS multiplexer used in fiberoptics. These micro-switches operate by orienting, and thus focusing, optical light rays. For proper operation, the switches must remain aligned. We concentrate on how thermal effects influence this alignment. Specifically, we examine how changes in the temperature of the MEMS structure influence its maximum and dynamic displacement. We also examine how the displacement is dependent on the thermal expansion coefficient of the material as a secondary dependence of temperature.

FEA Modeling

MEMS optical micro-switches take the form of long thin plates with the mirror located along the width direction. Typical lengths and widths are in the order of hundreds of microns, whereas a typical thickness may be below ten microns. Such a device can be idealized as a beam-like structure. In this study, we consider a multi-layered beam of dimensions shown in Figure 1. Table 1 lists the thickness of each layer and the material properties corresponding to each material.

The effect that temperature has on the mechanical behavior of this beam can be determined using FEA. The first step in the FEA modeling process is to consider only what is significant in the analysis. With this in mind, we consider whether a 2-D representation of the structure is adequate. Certainly, the primary motion of the beam is in the thickness direction. This primary motion can be modeled using a 2-D analysis in the length and thickness directions. The appropriateness of a 2-D analysis can best be established by comparing it with the results of equivalent 3-D analyses. Such comparisons give us limits on the applicability of the 2-D analysis. Specifically, we can ascertain at what excitation frequency levels secondary twisting or planar modes become present, thus, invalidating the 2-D analysis. ALGOR’s Linear Natural Frequency (Modal) application software was used to determine the modal behavior of the 2-D model of the structure (Figure 2). The same application was utilized to determine the modal behavior of several 3-D models of the beam, each with a different width. We considered other widths besides that of the beam under study in order to gain a better understanding of the dynamic behavior of such structures. Figure 3 depicts the first three mode shapes and frequencies of a 3-D model of the beam under study. Note how its first mode closely resembles that of the 2-D beam. Its second mode involves twist, which cannot be described by the 2-D model, and its third mode corresponds to the second mode of the 2-D beam. We can conclude that a 2-D analysis is only appropriate when modeling primary mode behavior for a beam with a width of 236 mm. Therefore, because we chose to consider a 2-D model for the remainder of this study, we will avoid secondary and higher mode behavior. Figure 4 summarizes some of the modal frequencies obtained for beams with other widths. As an aside, note how the frequencies of the secondary mode only coincide for the 100 mm beam. In this case, a 2-D analysis will capture the secondary mode behavior of such beams.

The eventual goal of the study is to examine the mechanical response of the beam when a harmonic load is applied to its free end in conjunction with a thermal load along its length. This response can be obtained in either the frequency or time domain. In this study, we consider only the time domain. Solutions in the time domain have the advantage that they can be considered virtual experiments. We employ ALGOR’s Mechanical Event Simulation (MES) to conduct these numerical experiments. The MES model of the beam consists of three layers of 2-D solid elements, with each layer composed of an isotropic material. Mesh refinement studies demonstrated that a convergent solution, capable of accurately describing bending behavior, could be achieved with a single element across the thickness of each layer as long as higher order elements with mid-side nodes were utilized.

The possibility of using a geometrically linear element formulation was considered by comparing simulations that utilized linear and nonlinear elements. These simulations consisted of applying a point force at the free end of the beam in the thickness direction. This load was increased linearly with respect to time from 0 to 0.0001 seconds, and was then maintained at a constant value. Several magnitudes of this constant value were considered, ranging from 1´ 10^-4 N to 8´10^-4 N. The beam exhibited sinusoidal, oscillatory motion once the load became constant; this is a consequence of simulating an energy-conserving, non-damped system. Figure 5 shows how the average value of this oscillatory motion differed for different load magnitudes and for simulations using linear and nonlinear formulations. In this figure, the base force was 1´ 10^-4 N, with the multiplier used to denote the higher loads. The two formulations give similar results for the lower force multipliers, but as Figure 6 clearly shows, the higher force multipliers do result in a significant difference between the two formulations – between 1% and 5%. Because we hoped to consider deflections possibly greater than 50 mm, which lies near the range where the formulations begin to differ significantly, we chose a nonlinear formulation for the remainder of the study. Actually, the extra computational effort involved in a nonlinear element formulation is practically insignificant when utilizing MES; run-times differed by less than 25% between the two formulations.

As mentioned above, the eventual goal of the study requires the application of a thermal load along the length of the beam. Once applied, the temperature of the beam should increase in a linear manner from the stress free reference temperature at its base to a value D T above this temperature at the free end. Such a temperature distribution, combined with the different thermal expansion coefficients of each layer, should induce bending in the beam. It is the combination of this thermally-induced bending with vibrational effects that must be considered in a proper design of such a MEMS structure. In a linear static stress solution, a thermal load is applied without temporal considerations. Because we aim to simulate a transient event, care has to be taken regarding how the thermal load is applied. If this load is imposed too abruptly, it can cause unwanted vibrations. Numerical experiments showed that increasing the temperature using a sinusoidal half period of duration 0.0004 seconds resulted in insignificant mechanical vibrations.

Our goal also required the application of a harmonic load. We chose a load composed of eleven harmonics whose amplitude decreases with increasing frequency. The frequency range was chosen so as to promote a vibrational response near the first natural frequency. The form of the load is given by Equation (1) in conjunction with the parameters listed in Table 2.

(1)

 

No force load is applied for times before t₀, which corresponds to a time 0.0001 seconds after the temperatures along the length reach constant values. The magnitude of the force load, F(t), was scaled using the parameter C so that it never exceeded 2´ 10^-4 N. This ensured that the displacement remains near 50 mm. Such a displacement would certainly significantly affect the optical alignment of the actual MEMS device.

Results

Simulations were run using five values for DT: 0, 50, 100, 150 and 200 °C. Figure 7 shows the displacement versus time curve for DT equal to 100 °C. This figure also shows that this temperature difference results in a displacement of 5.91 mm. A Fast Fourier Transform (FFT) of the curve after t₀ (0.0005 seconds) is given in Figure 8. The FFT clearly shows two peaks at 2441 Hz and 23440 Hz. It should be noted that all five simulations resulted in the same two frequency peaks. Thus, changing the value of the thermal load does not noticeably affect the frequency response of the beam. The displacement curve, however, is affected by the different values of DT. Figure 9 summarizes how the maximum displacement, as well as, the purely thermally-induced displacement depend on DT. A linear relationship is clearly visible from this figure. The fact that the difference between these displacements remains constant demonstrates that the variations prompted by different values of DT are primarily due to changes in the thermally-induced displacement. This is vital information to the designer of a MEMS multiplexer, especially since the maximum displacement changes dramatically with temperature.

This dramatic temperature effect was further studied by considering how changes in the thermal expansion coefficient of the middle layer affect the mechanical behavior of the beam. This layer contains the material with the highest expansion coefficient, thus increasing this value will result in the greatest change in the response. We considered only the case of DT equal to 100 °C, and simulated situations with this expansion coefficient multiplied by a factor of 2 and 10, respectively. The resulting values for the coefficient are high, but nevertheless are physical. Figure 10 summarizes how the maximum displacement as well as the thermally-induced displacement depend on the multiplier of the thermal expansion coefficient. Again, a linear response is obtained, as demonstrated by the nearly constant value for the difference of the displacements. The large sensitivity of the overall displacements to the thermal expansion coefficient is another vital piece of information for the designer. It should be mentioned that FFT plots of the curves obtained for the different expansion coefficients again showed prominent peaks at 2441 Hz and 23440 Hz. Thus the frequency of the beam seems constant, and even insensitive to a material property otherwise critically linked with the maximum displacement.

Conclusions

In this paper, we showed how to construct an FEA model used to analyze the mechanical behavior of a MEMS multiplexer. The model, though somewhat simplified, included thermal effects as well as accounting for geometrically nonlinear motion. This motion was coupled to the thermal load since the model’s three layers were composed of materials with different thermal expansion coefficients. ALGOR’s MES software was used to perform virtual experiments using this model. From these experiments, we found that the magnitude of the displacement of the beam is highly dependent on the thermal loads, and even more so on the expansion coefficient of the middle layer. Nevertheless, the frequency response of the model was insensitive to the temperature load or to the thermal expansion coefficient of the middle layer. This type of information is of great use to the designer of such a multiplexer, who must develop a design that maintains optical alignment under a varying range of thermal and vibrational operating conditions.

Table 1: Thickness and material properties of layers composing beam structure.