**4. MEMS coupled system-level model**

The modeling of MEMS provides a very challenging task in modern engineering. This field of research is inherently multiphysics of nature, since different physical phenomena are

Macromodels of Micro-Electro-Mechanical Systems (MEMS) 183

According to the assumption that conductors are equipotential, all the nodes connected to a certain conductor are subjected to the same voltage boundary conditions. The total current

> [ ( )] *k kk n <sup>d</sup> i CV V dt* .

The nodes used for electrostatic force application can be also used for monitoring the distance between the movable structure and the electrode. When this is equal to the transduction gap, the contact condition is reached and contact forces with the stiffness of the

Fcont = Kn(d − gapmin), (47)

Computation of the electrostatic forces adds some complexity to the development of the MEMS system-level model, but this is largely compensated by the speed-up simulation of

There is the special ANSYS' transducer element *TRANS126* which has the possibility to calculate the capacitance of a parallel plate capacitor model ( at particular nodes or as a

flowing in the conductor is simply given by the sum of the currents at those nodes:

Fig. 16. MEMS system-level simulation approach [33]

contact *Kn* , given by equation

can be applied to nodes.

the full model.

whole) as [47]

tightly intertwined at microscale. Typically, up to four different physical domains are usually considered in the analysis of microsystems: *mechanical, electrical, thermal* and *fluidic*. For each of these separate domains, well-established reduced order modeling and analysis techniques are available. However, one of the main challenges in the field of microsystems engineering is to connect models for the behavior of the device in each of these domains to equivalent lumped or reduced-order models without making unacceptably inaccurate assumptions and simplifications and to couple these domains correctly and efficiently.

Micromechanical membrane devices (capacitive pressure transducers, ultrasonic transducers), surface micro machined devices (RF switches, micro optical devices) as well as bulk micro machined devices (accelerometers, inclinometers, laser scanning mirrors) are driven or sensed by nearly parallel electrode pairs in many cases. The motion of these electrodes is strictly normal to the surfaces.

It means that MEMS electrical parts of these MEMS had to be combined with mechanical ones (fig.16). Typically, block-diagram descriptions and lumped-element circuit models for components are connected into a full system. Mostly, this description is used for functional analysis of a design concept.

The coefficients and electrostatic nodal force are obtained from the capacitancedisplacement function *C(w)* of the associated electrode portions and gap space. The function can be input by one of three means [44]:


The electrostatic forces acting on the movable conductor of the device are included in the model as nonlinear input forces, which are applied on nodes distributed over the conductor surface. Nodes divide the surface into *N* smaller portions. A lumped force is applied to *k*- th node at the center of each portion, in its preferential direction of movement *xi*. The capacitance *Ck* between the *k*-th portion and the fixed electrode of the device is computed, for its undeformed configuration, using an electrostatic analysis

The entity of each force *fk* is then approximated as:

$$f\_k = \frac{1}{2} \frac{\varepsilon A\_k}{\left(d\_0 + w\_i^k\right)^2} (V\_p - V\_n)^2 \tag{45}$$

where *d0* is the initial distance between the conductor and the electrode, *ε* is the relative dielectric permittivity, *<sup>k</sup> wi* is the displacement of the *k*- th node along the direction *xi* and (*Vp-Vn)* is the voltage difference between the structure and the fixed electrode. The capacitance *Ck* can be calculated from

$$\mathbf{C} \mathbf{\bar{e}} \mathbf{e} \mathbf{A} \boldsymbol{\psi} (\mathbf{d} \mathbf{\bar{u}} \, \mathbf{w}\_i^k \,) \,. \tag{46}$$

If *Ck* varies considerably with the deformation of the structure, then a series of electrostatic computations for different device deflections in its operation range can be performed. The results can be used for extracting the dependency *Ck* ( *<sup>k</sup> wi* ) and calculating the electrostatic force.

tightly intertwined at microscale. Typically, up to four different physical domains are usually considered in the analysis of microsystems: *mechanical, electrical, thermal* and *fluidic*. For each of these separate domains, well-established reduced order modeling and analysis techniques are available. However, one of the main challenges in the field of microsystems engineering is to connect models for the behavior of the device in each of these domains to equivalent lumped or reduced-order models without making unacceptably inaccurate assumptions and simplifications and to couple these domains correctly and efficiently. Micromechanical membrane devices (capacitive pressure transducers, ultrasonic transducers), surface micro machined devices (RF switches, micro optical devices) as well as bulk micro machined devices (accelerometers, inclinometers, laser scanning mirrors) are driven or sensed by nearly parallel electrode pairs in many cases. The motion of these

It means that MEMS electrical parts of these MEMS had to be combined with mechanical ones (fig.16). Typically, block-diagram descriptions and lumped-element circuit models for components are connected into a full system. Mostly, this description is used for functional

The coefficients and electrostatic nodal force are obtained from the capacitancedisplacement function *C(w)* of the associated electrode portions and gap space. The function

as analytical function if the electrode portions make up a plate capacitor geometry with

2

<sup>1</sup> ( ) <sup>2</sup> ( ) *k k k p n i <sup>A</sup> <sup>f</sup> V V d w* 

where *d0* is the initial distance between the conductor and the electrode, *ε* is the relative dielectric permittivity, *<sup>k</sup> wi* is the displacement of the *k*- th node along the direction *xi* and

*Ck=εAk/(d- <sup>k</sup> w ) . i* (46)

If *Ck* varies considerably with the deformation of the structure, then a series of electrostatic computations for different device deflections in its operation range can be performed. The results can be used for extracting the dependency *Ck* ( *<sup>k</sup> wi* ) and calculating the electrostatic

0

(*Vp-Vn)* is the voltage difference between the structure and the fixed electrode.

2

(45)

 as data table wherein the element subroutine interpolates values during solution. The electrostatic forces acting on the movable conductor of the device are included in the model as nonlinear input forces, which are applied on nodes distributed over the conductor surface. Nodes divide the surface into *N* smaller portions. A lumped force is applied to *k*- th node at the center of each portion, in its preferential direction of movement *xi*. The capacitance *Ck* between the *k*-th portion and the fixed electrode of the device is computed,

as polynomial approximation of a function given by data points;

for its undeformed configuration, using an electrostatic analysis

The entity of each force *fk* is then approximated as:

The capacitance *Ck* can be calculated from

force.

electrodes is strictly normal to the surfaces.

can be input by one of three means [44]:

homogeneous intermediate field;

analysis of a design concept.

Fig. 16. MEMS system-level simulation approach [33]

According to the assumption that conductors are equipotential, all the nodes connected to a certain conductor are subjected to the same voltage boundary conditions. The total current flowing in the conductor is simply given by the sum of the currents at those nodes:

$$i\_k = \frac{d}{dt} [C\_k(V\_k - V\_n)]\,.$$

The nodes used for electrostatic force application can be also used for monitoring the distance between the movable structure and the electrode. When this is equal to the transduction gap, the contact condition is reached and contact forces with the stiffness of the contact *Kn* , given by equation

$$\mathbf{F}\_{\rm cont} = \mathbf{K}\_n (\mathbf{d} - \mathbf{g} \mathbf{a} \mathbf{p}\_{\rm min}) \tag{47}$$

can be applied to nodes.

Computation of the electrostatic forces adds some complexity to the development of the MEMS system-level model, but this is largely compensated by the speed-up simulation of the full model.

There is the special ANSYS' transducer element *TRANS126* which has the possibility to calculate the capacitance of a parallel plate capacitor model ( at particular nodes or as a whole) as [47]

Macromodels of Micro-Electro-Mechanical Systems (MEMS) 185

As the displacement reaches the value , the hard stop will restrict its further increase, but input voltage *Vn* and can be further increased. Assume that the input voltage *Vn* is a superposition of a constant voltage *VDC* and a time dependent signal *V t* and that*V Vt DC* . The elastic-plastic properties of the points of contact are simulated by a spring with rigidity *Kn* and a damper of damping factor *b0*. When the beam center moves past *wmax*, it starts interacting with the spring that represents the contact. The damper is introduced to take into account the energy dissipation at the contact. The following equations

*Felec* 

*<sup>E</sup> <sup>V</sup> <sup>C</sup> A V <sup>F</sup>*

*w t w d* max 0

2 2 0 <sup>2</sup> 2 2 ( ,) <sup>2</sup>

*e*

*in eq C in*

*w d <sup>L</sup> dw t*

max max

The process of interaction simulated by the force *R* sometimes can be highly sensitive to the values of the model parameters (rigidity and damping factor), especially if . The interaction can induce high-frequency motions and slow down the rate of convergence considerably. Therefore, a good deal of attention must be paid to accurately modeling and

If the input voltage is increased more, there will be no equilibrium and the plate collapse takes place. In this case, a hard stop or some other arrangement must be introduced to limit the plate motion. During the plate collapse, the difference between the electrostatic force and the elastic force of the spring will continue to increase. Therefore, when the plate drops down into the hard stop, it is not enough to reduce the input voltage below to release the plate. The input voltage should be reduced more to make the electrostatic force at least equal to the elastic force. Hence, the plate capacitance exhibits the *hysteretic* behavior with

.

max

*w d* max *<sup>e</sup>*

Fig. 17. Capacitive –voltage characteristic an ultrasonic transducer

*elec*

( ,) 2 *L*

0

representing this process using experimental data.

( ) *n o for w w <sup>R</sup> K w w b w for w w* ,

are used to represent this model:

where *R* is the interaction force.

respect to voltage change.

$$C(w) = \frac{C\_0 d\_0}{d\_0 - w} = \frac{C\_0}{d\_0} (1 + \frac{w}{d\_0} + \frac{w^2}{d\_0^2} + \frac{w^3}{d\_0^3} + \frac{w^4}{d\_0^4} + \dots) \tag{48}$$

where *do* and *w* are the initial distance and the displacement between the plate.

The element has two nodes; the gap distance is calculated as the sum of the initial displacement and the difference of the nodal displacements in the direction of the element. The force is calculated by equation being similar to (44) so that in the constant voltage case

$$F = \frac{1}{2} \frac{\partial C(w)}{\partial d} V^2 = \frac{C\_0 V^2}{2} (\frac{1}{d\_0} + 2 \frac{w}{d\_0^2} + 3 \frac{w^2}{d\_0^3} + 4 \frac{w^3}{d\_0^4} + \dots) \,. \tag{49}$$

The element has also contact capabilities (47). It is possible to specify a minimal gap and a spring stiffness *Kn* for the repelling force.

The drawback of this kind of element is that it is limited to the case where the electrodes are (almost) parallel plates, so that the stroke/capacitance function can be evaluated from a single degree of freedom. But the extensions for rotation plates and for 2D cases were developed. The last one with a triangular shape (element *TRANS109*) is useful for simulating structures such as comb drivers and optical MEMS, in which capacitance between the device parts is generally a function of a two-directional displacement. *TRANS126* and *TRANS109* elements enable a huge reduction of the complexity of the system-level simulation.

Let's consider for example the system-level macromodel of an ultrasonic transducer which has two plates with the bottom electrode area *Ac* and the plate dimension *L* and which can be presented by nonlinear capacitance:

$$C\_{aq} = C\_o + (C\_{L/2} - C\_0)(1 - e^{-rt}),\tag{50}$$

where *C0* is the smallest capacitance in the absent of voltage *V:* <sup>0</sup> <sup>0</sup> *<sup>c</sup> e <sup>A</sup> <sup>C</sup> d* 

dielectric permittivity of the insulator; *w* is a plate deflection .

*CL/2* is the largest capacitance when a plate center displacement is calculated from the ROM macromodel equations (26): ( ,) 2 *L w t*

$$C\_{L/2} = \frac{\varepsilon\_0 A\_c}{d\_c - w\left(\frac{L}{2}, t\right)'} $$

where *<sup>e</sup> <sup>d</sup>* is an equivalent gap ( <sup>1</sup> 0 1 *ins <sup>e</sup> ins <sup>d</sup> <sup>d</sup> d d* ); 0 is the absolute dielectric permittivity of the vacuum, <sup>1</sup> is the relative dielectric permittivity of the poly silicon , *ins* is the relative

The largest value of this capacitance corresponds to the value . The electrostatic force acting on the capacitor surfaces is the Coulomb force: *w d* max 0

0 00 000 ( ) (1 ...) *Cd C ww w w C w d wd d ddd*

The element has two nodes; the gap distance is calculated as the sum of the initial displacement and the difference of the nodal displacements in the direction of the element. The force is calculated by equation being similar to (44) so that in the constant

1 () 1 ( 2 3 4 ...) 2 2

*d d ddd* 

The element has also contact capabilities (47). It is possible to specify a minimal gap and a

The drawback of this kind of element is that it is limited to the case where the electrodes are (almost) parallel plates, so that the stroke/capacitance function can be evaluated from a single degree of freedom. But the extensions for rotation plates and for 2D cases were developed. The last one with a triangular shape (element *TRANS109*) is useful for simulating structures such as comb drivers and optical MEMS, in which capacitance between the device parts is generally a function of a two-directional displacement. *TRANS126* and *TRANS109* elements enable a huge reduction of the complexity of the

Let's consider for example the system-level macromodel of an ultrasonic transducer which has two plates with the bottom electrode area *Ac* and the plate dimension *L* and which can

/2 0 ( )(1 ), *<sup>t</sup> CCC C e eq o L*

*CL/2* is the largest capacitance when a plate center displacement is calculated from

0

*<sup>L</sup> dw t* 

*ins*

is the relative dielectric permittivity of the poly silicon , *ins*

The largest value of this capacitance corresponds to the value . The electrostatic

*<sup>c</sup> <sup>L</sup> e*

); 0

*<sup>A</sup> <sup>C</sup>*

1 *ins <sup>e</sup>*

, 2

,

/2

0

dielectric permittivity of the insulator; *w* is a plate deflection .

force acting on the capacitor surfaces is the Coulomb force:

*<sup>d</sup> <sup>d</sup> d d* 

*C w C V ww w F V*

00 0

2 0

where *C0* is the smallest capacitance in the absent of voltage *V:* <sup>0</sup>

where *do* and *w* are the initial distance and the displacement between the plate.

234

234

2 2 3

23 4 0 00 0

. (49)

(50)

is the absolute dielectric permittivity

*w d* max 0

is the relative

<sup>0</sup> *<sup>c</sup> e*

> ( ,) 2 *L w t*

*<sup>A</sup> <sup>C</sup> d* 

(48)

voltage case

spring stiffness *Kn* for the repelling force.

be presented by nonlinear capacitance:

the ROM macromodel equations (26):

where *<sup>e</sup> <sup>d</sup>* is an equivalent gap ( <sup>1</sup>

of the vacuum, <sup>1</sup>

system-level simulation.

Fig. 17. Capacitive –voltage characteristic an ultrasonic transducer

$$\tilde{F}\_{\text{elec}} = -\frac{\partial E}{\partial w} = \frac{V\_{\text{io}}^2}{2} \frac{\partial C\_{\text{eq}}}{\partial d} = \frac{\mathcal{E}\_0 A\_c V\_{\text{io}}^2}{2\left(d\_s - w(\frac{L}{2}, t)\right)^2}.$$

As the displacement reaches the value , the hard stop will restrict its further increase, but input voltage *Vn* and can be further increased. Assume that the input voltage *Vn* is a superposition of a constant voltage *VDC* and a time dependent signal *V t* and that*V Vt DC* . The elastic-plastic properties of the points of contact are simulated by a spring with rigidity *Kn* and a damper of damping factor *b0*. When the beam center moves past *wmax*, it starts interacting with the spring that represents the contact. The damper is introduced to take into account the energy dissipation at the contact. The following equations are used to represent this model: ( ,) 2 *L w t w d* max 0 *Felec* 

$$R = \begin{cases} 0 & \text{for } w \le w\_{\text{max}} \\ K\_n(w\_{\text{max}} - w) - b\_o w' & \text{for } w > w\_{\text{max}} \end{cases}\_{\text{max}}$$

where *R* is the interaction force.

The process of interaction simulated by the force *R* sometimes can be highly sensitive to the values of the model parameters (rigidity and damping factor), especially if . The interaction can induce high-frequency motions and slow down the rate of convergence considerably. Therefore, a good deal of attention must be paid to accurately modeling and representing this process using experimental data. *w d* max *<sup>e</sup>*

If the input voltage is increased more, there will be no equilibrium and the plate collapse takes place. In this case, a hard stop or some other arrangement must be introduced to limit the plate motion. During the plate collapse, the difference between the electrostatic force and the elastic force of the spring will continue to increase. Therefore, when the plate drops down into the hard stop, it is not enough to reduce the input voltage below to release the plate. The input voltage should be reduced more to make the electrostatic force at least equal to the elastic force. Hence, the plate capacitance exhibits the *hysteretic* behavior with respect to voltage change.

Macromodels of Micro-Electro-Mechanical Systems (MEMS) 187

Providing a single representation of a MEMS operating in multiple physical domains, the electrical circuit approach is very convenient. Moreover, powerful mathematical techniques and circuit simulation programs are available for solving design tasks. It is possible to develop a library of schematic model for different MEMS elements and then use their combinations to build a system-level macromodels for entire rather complicated MEMS

In this chapter, the methods and issues encountered in the development of MEMS macromodels at the system level have been presented. System level modeling is the highest and most abstract level of modeling. This level requires various devices` linking of MEMS *component level models* – both electronic and micromechanical – into a micro-electromechanical system. *System-level models* of MEMS components are needed to allow a fast and

 Starting point for the extraction of a reduced order model (ROM or a macromodel) is already its description with a large ODE system, which is typically derived using physical modeling techniques based on Finite Element Method (FEM) which is rather time consuming. Macromodels application allows the extraction of lower order ODE system that reproduces the input/output behavior with good accuracy. Particular attention has been

There are special methods for generating ROM for MEMS components and entire MEMS based on FEM descriptions. To derive macromodels of smaller sizes different approaches (*Modal decomposition, Moment matching, Equivalent circuit presentation)* were developed. Usage of the reduced MEMS components models allows applying successfully modern circuit

Three automatic procedures to generate device reduced order macromodels, being based on full FEM/FDM models, were demonstrated in this chapter. Two of them are suitable for simulators with possibilities to get input information in the equation forms (ODE or OAE). The third one in opposite produces macromodels in circuit presentation and so it is more suitable for circuit simulators. The Modal ROM approach is based on using natural (modal or resonant) frequencies of MEMS structure and it is spread mostly in the USA and Asia. The Moment matching ROM approach is based on using the Krylov subspace for transfer functions and it is popular in Western Europe and Asia. The Equivalent circuit ROM approach is based on using a capacitive-inductive-resistive circuit model for mass, damping and stiffness matrices and it is used mostly in Eastern Europe. It is worth to notice that the Modal ROM approach requires some full ANSYS runs to perform a proper orthogonal decomposition during basic functions determination in opposite to the Moment matching and the Equivalent circuit ROM approaches for which it is enough to use ANSYS only for

It seems to be interesting and perspective trying to combine mentioned approaches, for example, to start with Krylov/Arnoldi reduction of ODE dimension, then to build the proper equivalent circuit for obtained ODE systems and finally to apply Y/Δ transformation

[1] S. D. Senturia, "CAD challenges for microsensors, microactuators, and microsystems,"

sufficiently exact investigation of their behavior to simulate entire MEMS.

posed in the chapter on the possibility to get a macromodel circuit presentation.

simulators in workflow for MEMS design on system level.

or n-port transformation for further reducing macromodel order.

*Proc. IEEE*, vol. 86, pp. 1611–1626, 1998.

FEM model matrices building.

**6. References** 

constructions.

**5. Conclusion** 

But if the insulation layer is rather thick its restrictive effect should be taken into account. If the maximal center displacement is equal to initial thickness air gap the moving plate touches the insulator top surface of the electrode when the input *Vn* voltage reaches the value of (*V*max may be calculated from simulation). But as soon as the input voltage will be decreased under this value the plate will leave the hard stop. So, its capacitance does not demonstrate the hysteretic behavior (fig.17). *w*max *V*max

Parameter *τ* in (50) defines an ultrasonic transducer frequency band and can be calculated through the plate displacement and its velocity , which are defied from ROM equations in the following way: ( ,) 2 *L w t* ( ,) ( ,) 2 2 *L L v tw t*

$$\tau = \frac{1}{3} \mathbf{w}(\frac{L}{2}, t) / \mathbf{v}(\frac{L}{2}, t) \tag{51}$$

The coefficient 3 appears in (51) due to the fact that a capacitance recharge to 98% for a time value which is equal approximately 3τ.

It is possible to see two included procedures in according to fig.18: one for development of ROM for an ultrasonic transducer plate, where a deflection and speed of central point`s deflection of transducer plate is calculated for the value *Vn*, and second - for determine of MEMS system- equivalent capacity value (SLM), using values of plate central point coordinates .Then the cycle of calculations recurs whereupon.

Fig. 18. ROM- system-level model coupled simulations

Instead of using two sequence procedures mentioned above it is possible using functional possibilities of the circuit simulator NetALLTED to built a single system-level equivalent circuit model for an ultrasonic transducer by introducing directly into the equivalent- circuit ROM of mechanical MEMS part the additional arbitrarily connected element (a Depended Source) with an informative function which is determined by equation (50) [42]. Optimization procedures of NetALLTED allow getting the desirable values of this transducer capacity and through it to get a desirable value of output signal of an ultrasonic transducer system-level model by the changing ROM parameters, which, in turn, are depended upon an ultrasonic transducer construction sizes and used material properties.

Providing a single representation of a MEMS operating in multiple physical domains, the electrical circuit approach is very convenient. Moreover, powerful mathematical techniques and circuit simulation programs are available for solving design tasks. It is possible to develop a library of schematic model for different MEMS elements and then use their combinations to build a system-level macromodels for entire rather complicated MEMS constructions.
