Steady state dynamics of pre‐stressed, piezoelectrically excited circular plates – a harmonic balance approach

We study the dynamics of radially symmetric circular multi‐layer plates. Geometric non‐linearity is taken into account through the Föppel‐von Karman strain‐displacement relationship. Employing the harmonic balance method we obtain a boundary value problem, which is numerically solved to obtain the non‐linear normal modes of the system.


Introduction
Typical MEMS (micro electro mechanical systems) for sound generation (speakers) are circular, plate-like structures. As actuation mechanism, piezoelectric layers are used, and the actuators are often operated in the large displacement regime. Additionally, the fabrication process may lead to static pre-stress in the structure. The aim of this contribution is to present a modeling strategy which incorporates above effects.

Modeling
We consider the speaker as a rotational symmetric plate consisting of multiple layers. In order to account for large displacements we apply the Föppel-Von Karman assumption. We adopt Voigt's linear description of piezoelectricity, and assume all materials are transversally isotropic, with the isotropy plane parallel to the reference surface of the plate. For thin piezoelectric layers, Gauss law can be integrated over the thickness of the layer to directly obtain an relationship between applied voltage and dielectric displacement [1,2], yielding a modified constitutive relation. Integrating over the complete thickness of the multi-layer plate, one obtains the plate level constitutive relations (relating reference surface displacements and section forces and moments). Writing the differential operators in cylindrical coordinates one obtains where () denotes derivatives with respect to the radial coordinate r. The unknown reference surface displacements are denoted by u i , where the indices i = r, ϕ are used to denote the radial and circumferential direction, respectively. Section forces are denoted by n i and section moments by m i . The parameters A ij , B ij and D ij are the relevant entries of the membrane, coupling and bending stiffness matrix of the plate, which may be functions of r. The know functionsn andm originate from the actuation of piezoelectic layers, and are dependent on time and the radial coordinate. The balance of momentum for plate level section forces and section moments (see e.g. [3]) in cylindrical coordinates reads where() denotes derivatives with respect to time, M the mass per area, and c i are damping coefficients. Inserting eqs. (1) to (4) into eq. (6) and the radial derivative of eq. (5) one obtains a system of non-linear partial differential equations (PDE) of order 4 in u z , and order 3 in u r , which is closed by the boundary conditions u r (0) = u r (0) = u z (0) = u z (0) = u z (R) = 0, n r (R) = −k n u r (R), m r (R) = k m u z (R), which allow for arbitrary linear stiffness k i at the outer edge r = R.
This transforms the PDE system into a system of orninary differential equations, which can be written in first order form x = f (x) with 7(2N + 1) unknowns, and solved by a collocation based boundary value problem solver [4,5]. This allows all coefficients to remain (known) functions of the radial coordinate

Results
All results were obtained for a typical MEMS actuator with a radius of 350 µm and a total thickness of h = 3.13 µm consisting of 2.33 µm thick silicon substrate, and an active piezoelectric Aluminium nitride layer with a thickness of 400 nm sandwiched between 160 nm thick Gold electrodes together with 40 nm thick Copper adhesion layers. We assume a static compressive pre-stress, e.g. from the production process.Considering three harmonics (N = 3) and neglecting damping, the non-linear normal modes were computed by a path following procedure starting from the linear oscillation modes. Figure 1 shows the typical stiffening non-linear response. The first mode shape (blue curves) appears three times, once at the base harmonic (starting from ≈ 0.5ω 0 due to pre-stress), and once for harmonic 2 and 3 starting at 1/2 and 1/3 of the base harmonic's frequency, respectively. Furthermore, we encounter a 1:3 interaction between mode 1 and mode 2.
Further modal interactions can be recognized by plotting the average kinetic energy in one oscillation period, versus the fundamental frequency of the response ( fig. 2). Each line represents a possible response, and whenever two lines intersect, modal interactions occur.