# Theoretical and experimental investigation of a 1:3 internal resonance in a beam with piezoelectric patches

Vinciane Guillot, Arthur Givois, Mathieu Colin, Olivier Thomas, Alireza Ture

Savadkoohi, Claude-Henri Lamarque

► To cite this version:

Vinciane Guillot, Arthur Givois, Mathieu Colin, Olivier Thomas, Alireza Ture Savadkoohi, et al. Theoretical and experimental investigation of a 1:3 internal resonance in a beam with piezoelectric patches. Journal of Vibration and Control, 2020, pp.107754632091053. 10.1177/1077546320910536. hal-02500418

HAL Id: hal-02500418

https://hal.science/hal-02500418

Submitted on 4 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Theoretical and experimental investigation of a 1:3 internal resonance in a beam with piezoelectric patches

Vinciane Guillot<fnref n="1" />, Arthur Givois<fnref n="2" />,<fnref n="3" />, Mathieu Colin<fnref n="2" />, Olivier Thomas<fnref n="2" />, Alireza Ture Savadkoohi<fnref n="1" /> and Claude-Henri Lamarque<fnref n="1" />

## Abstract

Experimental and theoretical results on the nonlinear dynamics of a homogeneous thin beam equipped with piezoelectric patches, presenting internal resonances, are provided. Two configurations are considered: a unimorph configuration composed of a beam with a single piezoelectric patch and a bimorph configuration with two collocated piezoelectric patches symmetrically glued on the two faces of the beam. The natural frequencies and mode shapes are measured and compared with those obtained by theoretical developments. Ratios of frequencies highlight the realization of 1:2 and 1:3 internal resonances, for both configurations, depending on the position of the piezoelectric patches on the length of the beam. Focusing on the 1:3 internal resonance, the governing equations are solved via a numerical harmonic balance method to find the periodic solutions of the system under harmonic forcing. A homodyne detection method is used experimentally to extract the harmonics of the measured vibration signals, on both configurations, and exchanges of energy between the modes in the 1:3 internal resonance are observed. A qualitative agreement is obtained with the model.

## Keywords

Nonlinear beam, internal resonance, vibration control, piezoelectric, modal energy exchange, experiments

## I. Introduction

Piezoelectric (PZT) materials constitute an efficient mean of coupling mechanical vibrations to an electrical circuit. Several applications are usually targeted such as micro/nano electromechanical systems (Bhugra and Piazza, 2017; Brand et al., 2015), vibration control (Bricault et al., 2019; Collet et al., 2008; Preumont, 2011; Soltani and Kerschen, 2015), and energy harvesting (Erturk and Inman, 2011; Jacquelin et al., 2011; Mam et al., 2016). In most of those cases, even if linear systems are traditionally considered for their simplicity, taking advantage of nonlinearities is of the most interest (Cao et al., 2015). Among others, in vibration control and energy harvesting, the synchronized switch strategies (see, e.g. Richard et al. (1999), Ducarne et al. (2012) and Lallart (2016)) are intrinsically nonlinear because of the electronic switching between two electromechanical states synchronized with the structure oscillations. When a primary structure is nonlinear, fully passive PZT nonlinear tuned vibration absorbers, that follow its change of frequency with the increasing vibration amplitude, have been recently proposed (Lossouarn et al., 2018). Another concept, the nonlinear energy sink, for which the absorber is intrinsically nonlinear, has been realized with a PZT shunt

(Silva et al., 2018; Zhou et al., 2014), thanks to nonlinear electrical circuits.

PZT materials are usually glued on elastic structures and both can present nonlinear behaviors. The PZT material nonlinearities were studied to understand their dynamics (Abdelkefi et al., 2012; Guyomar et al., 1997, 2011; Parashar and Wagner, 2004; Von Wagner and Hagedorn, 2002; Wolf and Gottlieb, 2001), but if deriving a proper nonlinear PZT law, thermodynamically consistent, has been already addressed, obtaining the values of the coefficients for a practical application seems to be still an open field of investigation. Moreover, the aforementioned publications use a classical electric enthalpy function which is a high-

## Corresponding author:

Vinciane Guillot, Univ Lyon, ENTPE, LTDS UMR CNRS 5513, 3 Rue Maurice Audin, Vaulx-en-Velin Cedex 69518, France.

Email: vinciane.guillot@entpe.fr

<footnote n="1">¹Univ Lyon, ENTPE, LTDS UMR CNRS 5513, France</footnote>

<footnote n="2">²Arts et Métiers, LISPEN EA 7515, France</footnote>

<footnote n="3">³Laboratoire de Mécanique des Structures et des Systèmes Couplés, Conservatoire National des Arts et Métiers, (LMSSC EA 3196), France</footnote>
order (smooth) polynomial in the strain, whereas in Leadenham and Erturks (2015), a low-order function that includes the absolute value of the strain (it is, thus, non-smooth), is proposed and seems to be closer to what is measured in practice. On the other hand, thin structures are well known to present geometric nonlinearities, especially cantilever beams, that have been thoroughly studied in the past (Anderson et al., 1994, 1996; Crespo Da Silva, 1988; Pai and Nayfeh, 1990). In this work, only the geometrical nonlinearities of the beam are taken into account, whereas the PZT material constitutive law is considered linear.

Contrary to linear systems, nonlinear oscillations can present mode coupling(s) and energy exchange(s) between different modes due to internal resonances. This phenomenon can be used for passive control and/or energy harvesting (Mook et al., 1985). The creation of internal resonances in nonlinear systems is influenced by their degrees of non-linearities. Let us suppose $\omega_q$, $\omega_r$, and $\omega_k$ are resonant frequencies of a system, where $q, r$, and $k$ represent integers corresponding to a specific frequency. For system which presents cubic nonlinearities, the internal resonance can occur if we have $\omega_q \approx \omega_r$, $\omega_q = 3\omega_r$, and $\omega_k \approx \pm 2\omega_q \pm \omega_r$. In a system with additional quadratic nonlinearities, besides the aforementioned conditions, the internal resonance can emerge if $\omega_q \approx 2\omega_r$ and $\omega_k \approx \omega_q \pm \omega_r$ (Nayfeh, 1979; Nayfeh and Balachandran, 1989). When an internal resonance is reached, a strong energy exchange between the modes can arise, depending on the nonlinear coefficients of the system.

The 1:3 internal resonance has been studied for plates (Sun et al., 2018; Zhang and Guo, 2012), MEMS (Czaplewski et al., 2019; Houri et al., 2019; Ramini et al., 2016), buckled beams (Emam and Nayfeh, 2013), and clamped-clamped beams (Ghayesh et al., 2012; Özkaya et al., 2008), demonstrating bifurcations and chaotic behaviors. Garg and Dwivedy (2019) used a beam with PZT materials excited parametrically, where the 1:3 internal resonance being achieved by an added mass.

The aim of this work was to investigate experimentally the behavior of an elastic homogeneous beam equipped with PZT patches presenting a 1:3 internal resonance between modes. Two configurations are examined: a unimorph configuration with only one PZT patch on one side of the beam, and a bimorph configuration for which two PZT patches are placed symmetrically on each side of the beam, as shown in Figure 1. As it will be observed, both configurations present theoretically quadratic and cubic nonlinearities, thus, 1:2 or 1:3 internal resonances can allow energy exchanges between modes. From the theoretical developments, it will be deduced that the position of the PZT patches on the beam can be used to tune the natural frequencies of some modes to create resonances such as $\omega_3 \approx 2\omega_4$ and $\omega_3 \approx 3\omega_2$, to achieve 1:2 and 1:3 internal resonances. No other internal resonances were identified between at least the first four modes. The focus of this study

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Photograph of the different systems under study: a homogeneous beam, a unimorph configuration, and a bimorph configuration. Sketch of the unimorph and bimorph configurations.</figcaption></figure>

was to investigate the behavior of the systems, while the 1:3 internal resonance is happening between the second and the third mode of the structures. The additive mass and stiffness of PZT patches are used to tune the natural frequencies of the structure in a 1:3 internal resonance relation to study the exchanges of energy between the two involved modes. The PZT patches are not used as sensors or actuators because this study is a preliminary work to investigate the positive use of nonlinearities to favor the energy exchanges, before using the structure as an active beam for energy harvesting and vibration damping applications.

The article organization is as follows: the experimental setup and methodology are first presented. Then, modal analyses of the systems are carried out and experimental results are compared with those obtained from an analytical model. Periodic solutions under harmonic excitation of the system, in the case of 1:3 internal resonance, are obtained both theoretically and experimentally. Finally, conclusions are provided.

## 2. Experimental methodology

Experimentations were achieved on three different systems: a homogeneous beam without PZT patches, a unimorph configuration, and a bimorph configuration (Figure 1). For the two PZT configurations, the PZT materials are glued between $x = x_1$ and $x = x_2$ measured from the clamping. $L_b$ denotes the length of the beam, and $\Delta = L_b - x_2$ is the distance between the end of the PZT patches and the end of the beam. The mechanical properties of the homogeneous beam, the PZT materials (we used a PIC151 PZT material from PI Ceramic), and the epoxy glue (from Mam et al. (2016)) are defined in Table 1.

### 2.1. Experimental setup

The apparatus is shown in Figure 2. The vibrating beam was fixed on the head of a shaker (Brüel & Kjær 4808), thanks to a homemade clamping system. An accelerometer was glued on this clamping system, allowing to measure the base acceleration prescribed to the structures. A scanning
laser vibrometer (Polytec PSV-400) was used to measure the velocity of several points of the vibrating beam. The shaker was driven by a power amplifier and a computer with input/output acquisition cards enabled to generate the input driving signal and to record the velocity and acceleration signals.

### 2.1.1. Detection of the mode shapes and natural frequencies.
From the experimental setup, an experimental modal analysis was performed. First, we measured the frequency response functions (FRFs) between the velocity of several points on the beam and the base acceleration. We first investigated the influence of the position of the PZT patches on the natural frequencies and mode shapes, to verify whether it is possible to tune the natural frequencies to achieve 1:2 or 1:3 internal resonance, that are

$$ \omega_j \approx 2\omega_i \quad \text{or} \quad \omega_j \approx 3\omega_i \qquad (1) $$

The first four natural frequencies were detected for the different length of the homogeneous beam $L_b$. In practice, the PZT patches are glued at a fixed distance $\Delta$ of the free end of the beam and we changed the position of the beam in the clamping system, thus, modifying $L_b$ (see Figure 1). Because the length $L_p = x_2 - x_1$ of the PZT patches is fixed, $x_1$ and $x_2$ are changed according to $x_1 = L_b - \Delta - L_p$ and $x_2 = L_b - \Delta$. For the unimorph configuration, $\Delta = 60$ mm and the length $L_b$ of the beam has been increased from 110 mm to 170 mm with steps of 5 mm. For the bimorph configuration, $\Delta = 70$ mm and the length $L_b$ of the beam has been increased from 120 mm to 170 mm with steps of 5 mm. In both cases,

<table><caption>**Table 1.** Mechanical properties of the homogeneous beam (experimentally identified), the piezoelectric material PIC151 (from Thomas et al. (2009)), and epoxy glue (from Mam et al. (2016)).</caption><thead><tr><th></th><th>ρ (kg m<sup>-3</sup>)</th><th>Y (GPa)</th><th>h, h<sub>p</sub> (mm)</th><th>b, b<sub>p</sub> (mm)</th></tr></thead><tbody><tr><td>Beam</td><td>7810</td><td>179</td><td>0.5</td><td>13</td></tr><tr><td>PIC151</td><td>8500</td><td>66.7</td><td>0.5</td><td>10</td></tr><tr><td>Epoxy glue</td><td>8000</td><td>5</td><td>0.01</td><td>10</td></tr></tbody></table>

the length of the PZT patches is $L_p = 50$ mm. Around the particular lengths for which the 1:2 and 1:3 internal resonance frequency ratios of equation (1) were obtained, the step was reduced to 1 mm. To obtain the mode shapes, the FRFs were measured at several points along the length of the beam, on a line of points in the middle of the upper surface. Then, the operational deflection shapes at resonances were plotted as an estimation of the mode shapes.

### 2.1.2. Forced response measurements.
We were interested in the velocity response of the system when excited by a harmonic base acceleration of frequency $\Omega$. In our studies, the response of the structure was periodic and we used a homodyne detection to measure the amplitude and the phase of several harmonics of the signals at fixed excitation frequencies. The structure was excited at a certain amplitude of acceleration at several successive frequencies $\Omega$, stepped in a given frequency band of interest. Because of the retroaction of the vibrating structure on the shaker, a control loop was used to keep the excitation amplitude constant during operations. For each frequency measurement, a time delay was prescribed to wait for the end of the transient and to reach a steady state. Then, each amplitude of harmonic component was calculated by the computer from the velocity signal of the laser vibrometer with the homodyne detection strategy. It consists in multiplying the velocity signal by $\sin(h\Omega t)$ and $\cos(h\Omega t)$ functions and taking the average of the result to extract the amplitude of the $h$th harmonic (see, e.g. Monteil et al. (2015) and Denis et al. (2018) for details).

Because we focus in this article on a 1:3 internal resonance between the second and third modes of the structure, with a direct driving of the lower mode, the frequency bandwidth is defined around the second natural frequency $\omega_2$. Then, to measure precisely the response of each of the involved modes (we assume that the response of the structure is the modal superposition of only the two modes involved in the internal resonance), for each frequency measurement, two experiments were performed, one after the other. To obtain the response of the second mode, the laser was first pointed on a vibration node of the third mode

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Photograph of the unimorph piezoelectric beam clamped on the head of the shaker (left). Sketch of the experimental setup (right).</figcaption></figure>
shape, and, second, the response of the third mode was obtained by pointing the laser on a vibration node of the second mode shape.

### 3. Theoretical modeling

We consider the cantilever beams with PZT patches sketched in Figure 1. We denote by $b \times h$ the width and the thickness of the elastic layer cross section and by $b_p \times h_p$ the ones of the PZT patches. Because of the clamped/free boundary conditions, we assume that in large amplitude vibrations, the beam remains inextensible (Dowell and McHugh, 2016). This assumption may be broken by the PZT patches, in particular in the asymmetric configuration, because their mode of operation is extension/compression. However, we assume that this effect is small because of the axial stiffness of the beam, which is very large than that of the bending one, and also because we are interested in bending vibrations that are assumed to be uncoupled of axial vibrations due to the clamped/free boundary conditions. Parametric excitation is also not considered. Moreover, we assume a linear PZT constitutive law, mostly because no available numerical values for those parameters are available in the literature. As shown in the Appendix 1 section and in Guillot et al. (2019), the axial/bending coupling due to an asymmetry of the lamination of the beam (this is the case for the unimorph configuration (Ducarne et al., 2012)) is exactly canceled by the inextensibility property.

The governing electromechanical equations of the system, involving the transverse displacement $v(x, t)$, where $t$ is the time, and the voltage $V(t)$ across the PZT patches, can be written as

$$ m(x)\ddot{v} + [D(x)v^{\prime\prime}]^{\prime\prime} + (v^{\prime}[v^{\prime}D(x)v^{\prime\prime}]^{\prime})^{\prime} + \left[ v^{\prime} \int_{x}^{L} \frac{m(x)}{2} \left( \int_{x}^{0} v^{\prime 2} dx \right)^{\prime\prime} dx \right]^{\prime} + \left[ \delta_{x2}^{\prime}(x) - \delta_{x1}^{\prime}(x) \right] \Theta V \left( 1 + \frac{v^{\prime 2}}{2} \right) = p(x, t) \quad (2) $$

$$ Q = \Theta[v^{\prime}]_{x1}^{x2} \left(1 + \frac{v^{\prime 2}}{2}\right) + CV \quad (3) $$

where $[v^{\prime}]_{x1}^{x2} = v^{\prime}(x_2) - v^{\prime}(x_1)$ as observed in Ducarne et al. (2012) and Guillot et al. (2019).

Equation (2) represents the dynamical equilibrium of the beam and equation (3) the electrical state of the PZT material.

In the aforementioned equations, $m(x)$ is the mass per unit length, $D(x)$ is the bending stiffness, $p(x, t)$ is an external force per unit length, $Q(t)$ is the electric charge contained in one of the electrodes of the PZT patches, and $\delta(x)$ is the Dirac $\delta$ function. For the unimorph configuration, the PZT coupling coefficient is $\Theta = e_{31}b_p(h + h_p)/2$ and the

electric capacitance is $C = \epsilon_{33}b_p(x_2 - x_1)/h_p$, where $e_{31}$ is the PZT constant and $\epsilon_{33}$ is the dielectric permittivity of the PZT material (Ducarne et al., 2012), see Table 2. First equation predicts dynamical equilibrium of the beam and electrical state of the PZT. For the bimorph configuration, because the two collocated patches are in series, $\Theta$ is the same and the capacitance $C$ is half that of the unimorph configuration (Ducarne et al., 2012).

#### 3.1. Tuning of the natural frequencies for 1:2 and 1:3 internal resonance

The linear part of the governing equation (2) is now considered in short circuit ($V = 0$)

$$ m(x)\ddot{v}(s, t) + [D(x)v^{\prime\prime}]^{\prime\prime} = 0 \quad (4) $$

and the associated eigenmodes of the stepped beams (the mass and stiffness addition of the PZT patches are included) are computed, as explained in Ducarne et al. (2012) and Guillot et al. (2019). To be rigorous, in the case of the unimorph configuration, which has a nonsymmetric lamination, the linear axial/bending coupling, not present in equation (4) because of the inextensibility condition, is included in those computations, as explained in Ducarne et al. (2012).

The natural frequencies and mode shapes are computed for several values of the beam length $L_b$ and of the placement of the PZT patches, defined by $\Delta$ (see Figure 1). The length of the PZT patches is fixed at $L_p = x_2 - x_1 = 50$ mm. Because the bending inertia and stiffness depend on the location of the PZT patches on the beam, the natural frequencies depend on the PZT patches' positions, and one can find values of $(L, \Delta)$ for which the relations one are fulfilled to create a 1:2 internal resonance between mode 4 and mode 3 or a 1:3 internal resonance between mode 3 and mode 2. To illustrate these results, the maps of the ratios $\omega_4/\omega_3$ and $\omega_3/\omega_2$ are plotted as functions of $(L, \Delta)$ in Figure 3. They have been computed with the parameters of the beam gathered in Table 1 without considering the epoxy layer. Those plots have to be considered as a qualitative first insight into the tuning of the natural frequencies, before a precise experimental tuning, as presented in the next section. For the unimorph configuration and $\Delta = 60$ mm, the special lengths $L_b = 0.142$ mm and $L_b = 0.129$ mm allow to reach to the 1:2 and 1:3 internal resonances, respectively. As for the bimorph configuration, we have $\Delta = 70$ mm and the

<table><caption>**Table 2.** Piezoelectric and dielectric constant parameters of the piezoelectric material PIC151 (from Thomas et al. (2009)), with $\epsilon_0 = 8.854 \times 10^{-12}$ F/m.</caption><thead><tr><th>ε<sub>33</sub> (F/m)</th><th>e<sub>31</sub> (C/m<sup>2</sup>)</th></tr></thead><tbody><tr><td>2068 ε<sub>0</sub></td><td>-14</td></tr></tbody></table>
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>Figure 3. Theoretical map of the realization of the internal resonance 1:2 and 1:3. Plots of the frequency ratios $\omega_4/\omega_3$ and $\omega_3/\omega_2$ as functions of the beam length $L_b$ and the position of the piezoelectric patches $\Delta$, for the unimorph and bimorph configurations. The solid black lines give the 1:2 and 1:3 contour lines. In those computations, the piezoelectric patches are supposed to be perfectly glued on the elastic layer without the additional epoxy layer.</figcaption></figure>

special lengths are $L_b = 0.158$ mm and $L_b = 0.152$ mm to achieve the prescribed internal resonances.

## 3.2. Periodic solutions for the 1:3 internal resonance

We now discuss the full nonlinear equations of motions. Because we target a 1:3 internal resonance between the second and the third modes, we assume that the response of the system in the steady state, under harmonic forcing of frequency $\Omega$, is mainly governed by those two modes. This way, only the second and third modes are supposed to mainly exist, neglecting the effects of other possible modes, such as the first mode because it is not resonant. We, thus, expand $v(x, t)$ on the two corresponding mode shapes ($\phi_2(x)$ and $\phi_3(x)$)

$$v(x, t) = \phi_2(x)r_2(t) + \phi_3(x)r_3(t) \quad (5)$$

where ($r_2(t)$ and $r_3(t)$) are the corresponding unknown modal coordinates of the second and third mode, respectively. Assuming an open-circuit condition (the current is equal to zero, which leads to $Q = 0$, see Trindade and Benjeddou (2009) and Ducarne et al. (2012)), introducing equation (5) into (2) and (3), multiplying the result successively by $\phi_2(x)$

and $\phi_3(x)$, and finally using the orthogonality properties of the modes, one obtains

$$\begin{aligned} \ddot{r}_2 + 2\mu_1\dot{r}_2 + \omega_2^2 r_2 + N_2 r_2^3 + N_3 r_2^2 r_3 + N_4 r_2 r_3^2 + N_5 r_3^3 \\ + N_6 r_2 (r_2^2)^{\ddagger} + N_7 r_2 (r_2 r_3)^{\ddagger} + N_8 r_2 (r_3^2)^{\ddagger} + N_9 r_3 (r_2^2)^{\ddagger} \\ + N_{10} r_3 (r_2 r_3)^{\ddagger} + N_{11} r_3 (r_3^2)^{\ddagger} + N_{12} V = F_2 \cos(\Omega t) \end{aligned} \quad (6)$$

$$\begin{aligned} \ddot{r}_3 + 2\mu_2\dot{r}_3 + \omega_3^2 r_3 + M_2 r_2^3 + M_3 r_2^2 r_3 + M_4 r_2 r_3^2 + M_5 r_3^3 \\ + M_6 r_2 (r_2^2)^{\ddagger} + M_7 r_2 (r_2 r_3)^{\ddagger} + M_8 r_2 (r_3^2)^{\ddagger} + M_9 r_3 (r_2^2)^{\ddagger} \\ + M_{10} r_3 (r_2 r_3)^{\ddagger} + M_{11} r_3 (r_3^2)^{\ddagger} + M_{12} V = F_3 \cos(\Omega t) \end{aligned} \quad (7)$$

$$0 = L_1 V + L_2 r_2 + L_3 r_3 \quad (8)$$

Modal linear damping terms of coefficients $\mu_1$ and $\mu_2$ have been added in the aforementioned equations. Moreover, the analytical expressions of all the $N_i$, $M_i$, $i = 2, \dots, 12$ coefficients can be found in Guillot et al. (2019). They are here computed according to the geometrical and material parameters of the systems under study. Finally, $F_2$ and $F_3$ are the modal forcing amplitude, that are both nonzero because of the base acceleration driving, equivalent to a uniform force per unit length $p(x, t) = -m\gamma_e \cos(\Omega t)$ (where $\gamma_e$ is the amplitude of the base acceleration prescribed by the shaker)
and, thus, not orthogonal to the mode shapes. The modal damping coefficients are experimentally estimated with the modal damping factors $\xi_i$ ($\mu_i = \xi_i\omega_i$), computed with the -3 dB bandwidth $\Delta\omega_i$ at resonance ($\xi_i = 2\Delta\omega_i/\omega_i$), in the case of the lowest amplitude excitation for which no non-linear coupling occurs. The coefficients $N_i$ and $M_i$ from equation (6) were calculated from their analytical expressions deduced from equations (2) and (3) and the mechanical and electrical known values from Tables 1 and 2.

The periodic solutions of this system, in the steady state, are numerically computed by a higher order harmonic balance coupled to a continuation method (the asymptotic numerical method) implemented into software Manlab (Arquier et al., 2005–2019). Because of the targeted 1:3 internal resonance, we are mainly interested in the first and third harmonic of the second and third mode modal coordinates $r_2(t)$ and $r_3(t)$, respectively, for different forcing amplitudes $F$. Theoretically, because the nonlinear forces in equation (6) are odd, all even harmonics are zero in the periodic responses (they can be obtained after a symmetry-breaking bifurcation which has not been observed here). Thus, in the theoretical results, the second harmonic of the modal coordinates is equal to zero and was not displayed (see Figures 4 and 5). The displacement $v(x, t)$ can be constructed from equation (5).

To look at each mode independently, as explained in the experimental methodology section, the motion of the beam is measured at two positions successively. First, at the node of the third mode shape ($x = x_{f3}$), to measure the contribution of mode 2, and then, at the node of the second mode shape ($x = x_{f2}$), to investigate the contribution of mode 3. Thus, from equation (5), the beam’s displacements are defined by

$$ 
\begin{aligned}
\text{At } x = x_{f3} \quad v(x_{f3}, t) &= \phi_2(x_{f3})r_2(t) \\
\text{At } x = x_{f2} \quad v(x_{f2}, t) &= \phi_3(x_{f2})r_3(t)
\end{aligned}
\qquad (9) $$

For the unimorph configuration, the nodes of the second and third modes closest to the tip end are defined as $x_{f2} \approx 0.099$ mm and $x_{f3} = 0.113$ mm, respectively, as for the bimorph configuration, we have $x_{f2} \approx 0.114$ mm and $x_{f3} \approx 0.133$ mm.

In addition, because of the excitation bandwidth defined around the second mode and the internal resonance 1:3 between $\omega_3$ and $\omega_2$, we have

$$ \Omega \approx \omega_2, \quad \omega_3 \approx 3\omega_2 \approx 3\Omega \qquad (10) $$

The evolutions of the different harmonics amplitudes of $v(x = x_{f2}, t)$ and $v(x = x_{f3}, t)$ versus $\Omega$ are presented in Figures 4 and 5.

<figure><img src="image_8.png" /><figcaption>**Figure 4.** Theoretical results for the periodic forced response of the unimorph configuration in 1:3 internal resonance. Harmonics one and three of the displacements of the beam $v(x = x_{f3}, t)$ (low frequency—mode 2—response, left column) and $v(x = x_{f2}, t)$ (high frequency—mode 3—response, right column), as a function of the excitation frequency $\Omega$, for different amplitude of base acceleration (each color corresponds to a given base acceleration).</figcaption></figure>
<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>**Figure 5.** Theoretical results for the periodic forced response of the bimorph configuration in 1:3 internal resonance. Harmonics one and three of the displacements of the beam $v(x = x_3, t)$ (low frequency—mode 2—response, left column) and $v(x = x_2, t)$ (high frequency—mode 3—response, right column), as a function of the excitation frequency $\Omega$, for different amplitude of base acceleration (each color corresponds to a given base acceleration).</figcaption></figure>

Because of the 1:3 internal resonance, one would expect that the first harmonic of the low-frequency (LF) mode two and the third harmonic of the high-frequency (HF) mode 3 be the dominant harmonics, in the case of a strong transfer of energy from mode two to mode 3. On the contrary, those simulations show that the effect of the 1:3 internal resonance, in term of energy transfer between the modes, is not as strong as it is for other cases (see, e.g. Cao et al. (2015), Thomas et al. (2003, 2007) and Monteil et al. (2015)). It would tend to prove that for cantilever beams in direct bending driving, a 1:3 internal resonance between two modes (that have been tuned here by the structural effect of the PZT patches) does not lead to a strong energy exchange. To the knowledge of authors, this result is original because no previous work addressed the 1:3 internal resonance between two modes of a cantilever beam. However, further investigations are needed to fully understand what has to be done on the system to boost the exchanges of energy. The present case of a cantilever beam is here rendered complex by the presence of 20 resonant cubic terms in the oscillators (6), which a priori have all influence on the energy exchange. Finally, these results can probably be explained by considering that the clamped/free boundary conditions lead to very small geometrically nonlinear effects (the first mode,

for instance, is slightly hardening), as compared for instance with clamped-clamped beams or 2D structures, such as plates and shells, for which internal resonances have very strong effects (Amabili, 2008; Monteil et al., 2015; Thomas et al., 2003, 2007). No experiment was performed on the beam without PZT materials because for a homogeneous cantilever beam, no 1:3 ratio of the natural frequencies is observed in theory, which prevents any internal resonance to easily occur.

# 4. Experimental investigation and comparison to theory

## 4.1. Tuning of the natural frequencies

With the protocol explained in the previous sections, we first measured the natural frequencies of the homogeneous beam as a reference. By plotting the results as a function of the theoretical frequencies of a cantilever beam $\omega_k = \beta_k^2 h / L^2 \sqrt{Y/(12\rho)}$ (with $\beta_1 = 1.875$, $\beta_2 = 4.694$, and $\beta_3 = 7.855...$), it is possible to estimate the value of the Young's modulus $Y$ of the stainless steel because all the other parameters are known (the geometry has been measured and the mass density $\rho$ has been obtained by weighing the beam). We obtained $Y = 179$ GPa.
Then, we measured $\omega_2$, $\omega_3$, and $\omega_4$ for several values of the beam length, for the two PZT beam configurations (unimorph and bimorph). Those values are compared with the theoretical ones in Figure 6. A good agreement is obtained between experiments and theory, better for the unimorph configuration than for the bimorph one. The effect of the epoxy layer has also been investigated and is found to slightly change the natural frequency values. Those discrepancies can be explained by usual experimental characteristics not taking into account in the model: a nonperfect clamping (the clamping device, made in steel, is not infinitely rigid), some material values not experimentally identified (the material characteristics of the PZT and the epoxy layer, see Table 1), the Euler–Bernoulli assumptions of a rigid cross section without transverse shear, etc. All those characteristics could be verified with numerical computations (with a commercial finite-element code) and are out of the scope of this study because our interest is in the experimental tuning of the natural frequencies.

As predicted by Figure 3, it is shown that there exist particular values of the beam length $L_b$ for which $\omega_4/\omega_3 \approx 2$ and $\omega_3/\omega_2 \approx 3$, so that 1:2 and 1:3 internal resonances are at hand. Precisely, setting $L_b = 129$ mm for the unimorph configuration (Figure 6(c)) and $L_b = 145$ mm for the bimorph configuration (Figure 6(d)) will lead to the required

frequency ratio to favor a 1:3 internal resonance between the second and third modes of the beam. As for the internal resonance 1:2, the special lengths $L_b = 0.139$ mm and $L_b = 0.153$ mm allow the desired ratio to occur for the unimorph and bimorph configurations, respectively. The corresponding natural frequency values can be found in Tables 3 and 4. No other internal resonances were found, at special lengths $L_b$, between the other first four modes of the structure. Tables 3 and 4.

## 4.2. Mode shapes of the system

The homogeneous beam and both PZT configurations were set at a length $L_b = 170$ mm, and experiments were conducted to estimate the mode shapes of the three systems. The experimental results are shown in Figure 7 along with the theoretical mode shapes. An excellent agreement between the test and theory is obtained. One can also observe that the mode shapes show an overall modification (a change of curvature) in the area where the PZT patches are glued, as compared with the mode shape of the homogeneous beam (without PZT materials). Moreover, this kind of figure is useful to estimate the positions $x_1$ and $x_3$ of the nodes of mode 2 and 3 for the measurements of the 1:3 internal resonance presented in the following section. Note that the

<figure><img src="image_13.png" /><figcaption>**Figure 6.** Experimental verification of the creation of the internal resonances 1:2 and 1:3 compared with theoretical results. Blue points correspond to the experimental open-circuit condition and red points to the short circuit one. Theoretical frequencies are calculated in short-circuit condition. Solid line: computation without the epoxy layer and dotted line: computation with the epoxy layer. (a) and (c) unimorph configuration with $\Delta = 60$ mm; (b) and (d) bimorph configuration with $\Delta = 70$ mm.</figcaption></figure>
beam lengths $L_b$ used in the following section are not the length shown in Figure 7 because of the required particular tuning of $\omega_3/\omega_2 = 3$.

### 4.3. Periodic solutions with the 1:3 internal resonance

Both PZT configurations are settled at the special length $L_b$ so that $\omega_3/\omega_2 = 3$, to naturally favor the 1:3 internal resonance ($L_b = 129$ mm for the unimorph configuration and

<table><caption>**Table 3.** Comparison of theoretical and experimental results for the natural frequencies implicated in the internal resonance 1:2.</caption><thead><tr><th></th><th>ω<sub>3</sub></th><th>ω<sub>3,theo</sub></th><th>ω<sub>4</sub></th><th>ω<sub>4,theo</sub></th></tr></thead><tbody><tr><td>Unimorph configuration</td><td>2342</td><td>2280</td><td>4697</td><td>4584</td></tr><tr><td>Bimorph configuration</td><td>1857</td><td>1823</td><td>3709</td><td>3642</td></tr></tbody></table>

<table><caption>**Table 4.** Comparison of theoretical and experimental results for the natural frequencies implicated in the internal resonance 1:3.</caption><thead><tr><th></th><th>ω<sub>2</sub></th><th>ω<sub>2,theo</sub></th><th>ω<sub>3</sub></th><th>ω<sub>3,theo</sub></th></tr></thead><tbody><tr><td>Unimorph configuration</td><td>975.2</td><td>962.9</td><td>2944</td><td>2894</td></tr><tr><td>Bimorph configuration</td><td>711</td><td>686.4</td><td>2142</td><td>2894</td></tr></tbody></table>

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 7.** Comparison between the theoretically (lines) and experimentally (points) obtained second (top) and third (bottom) mode shapes for: the homogeneous beam (continuous line with cross points) without piezoelectric materials, the unimorph configuration (dashed line with x points), and the bimorph configuration (dotted line with star points).</figcaption></figure>

$L_b = 145$ mm for the bimorph configuration). The system is excited at different levels of amplitude of the base acceleration. We assume that the modes which are in internal resonance are those which dominate the response of the system, as shown by equation (5). Then, as shown by equation (9), pointing the laser at a node of the third mode ($x = x_f$) enables to measure the response of mode two only and pointing the laser at a node of the second mode ($x = x_f$) enables to measure the response of mode three only. A stepped frequency sweep is then performed in a narrow frequency band around $\omega_2$: we choose to directly drive the LF mode (mode 2) around its resonance and to observe a transfer of energy toward the HF mode (mode 3). The results are shown in Figure 8 for the unimorph configuration and Figure 9 for the bimorph configuration.

From the theoretical analysis, it was shown in Figures 4 and 5 that no strong transfer of energy was observed between the first harmonic (H1) of the LF mode and the third harmonic (H3) of the HF mode. Although, the present experimental results show the opposite. For the unimorph configuration (Figure 8), the amplitude of the third harmonic of the HF mode is about five times that of its first harmonic, whereas H3 for the LF mode is negligible with respect to H1, meaning that the dominant harmonics of the HF mode is H3. Moreover, the third harmonics of the HF mode is about 40 times larger than H1 of the LF mode, showing that this third harmonic is not directly excited by the third harmonic of the LF mode. Analogous results are shown in Figure 9 for the bimorph configuration. This confirms that a clear transfer of energy from the LF mode to the HF mode, associated to the 1:3 internal resonance, is at hand.

The fact that the energy transfer is observed experimentally, whereas it is not clearly predicted by the theoretical model is now addressed. In the theoretical model, only geometrical nonlinearities are considered, with a linear PZT constitutive law. Several works in the past showed that nonlinearities in the PZT constitutive law are experimentally observed (among others, see Guyomar et al. (1997, 2011) for a PZT stack with a 33 effect and Parashar and Wagner, (2004) for a sole PZT ceramics) In particular, in the case of cantilever beams with PZT ceramics in a 31 effect (see, e.g. Von Wagner and Hagedorn (2002) and Leadenham and Erturk (2015)), the nonlinear PZT material effect seems to be of the same order of magnitude than the geometrical one because the first mode of the cantilever beam is observed with a softening nonlinear effect, whereas it is predicted hardening in the case of geometrical nonlinearities only. From a theoretical point of view, it is also shown in Guillot et al. (2019) that the nonlinear PZT constitutive law adds more cubic terms in the model and also additional quadratic nonlinearities. However, in practice, some additional work, out of the scope of the present article, has to be carried out to correctly experimentally identify all the parameters of the constitutive law,
<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /><figcaption>**Figure 8.** Experimental results for the periodic forced response of the unimorph configuration in 1:3 internal resonance. First three harmonics of the signals measured at $x = x_f$ (low frequency—mode 2—response, left column) and $x = x_f$ (high frequency—mode 3—response, right column) on the beam, as a function of the excitation frequency $\Omega$, for different amplitude of base acceleration (each color corresponds to a given base acceleration).</figcaption></figure>

to obtain a usable theoretical model and compare its results to experiments. The asymmetry of the structure in the transverse direction, in the case of the unimorph configuration, canceled by the inextensibility constraint (see Appendix 1), can also be a factor that could change the nonlinear behavior of the structure. In shells and laminated plates, it adds quadratic nonlinearities (Lazarus et al., 2012; Thomas et al., 2005). We believe that the inextensibility constraint is very realistic and, thus, that this effect is probably negligible in practice.

Experimentally, some small amplitudes of the second harmonics are detected, which are naturally not present in our theoretical developments because no quadratic nonlinearities are included in the model. From Guillot et al.

(2019), some of the assumed nonlinear terms of the PZT materials could be responsible for the quadratic behavior, leading to the detection of second harmonic in the displacements. This can also be explained by unavoidable imperfections in the beam or by some aeroelastic damping effects Mam et al. (2016). One can also observe that the ratio between H2 and H1 for the LF mode is of about 1/20 for the unimorph configuration, whereas it is of about 1/75 for the bimorph configuration. This can be explained by the nonsymmetric lamination of the unimorph configuration, which breaks the transverse structural symmetry of the beam. This effect is not taken into account in the model because of the inextensibility condition. The Appendix 1 section gives more details about this point.
<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>

<figure><img src="image_24.png" /></figure>

<figure><img src="image_25.png" /></figure>

<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /><figcaption>**Figure 9.** Experimental results for the periodic forced response of the bimorph configuration in 1:3 internal resonance. First three harmonics of the signals measured at $x = x_f$ (low frequency—mode 2—response, left column) and $x = x_s$ (high frequency—mode 3—response, right column) on the beam, as a function of the excitation frequency $\Omega$, for different amplitude of base acceleration (each color corresponds to a given base acceleration).</figcaption></figure>

It is observed that the amplitude of the third harmonic of the second mode in Figure 8(e) shows different behavior for the different amplitudes of excitation. Because the amplitudes of the harmonic are really low (0.02 mm maximum), we assume that they are more related to noise and were not representative of the real amplitude of the third harmonic.

## 5. Conclusion

In this study, several results have been obtained about the design of a cantilever beam with PZT patches to favor internal resonance and energy exchanges between the modes. First, both theory and experiments allowed us to confirm that the position of the PZT patches on the beam can help to tune the natural frequencies of the beam to achieve 1:2 and 1:3 internal resonances, between the

fourth and third mode and the second and third mode, respectively. Experimentally, both 1:2 and 1:3 internal resonances have been tested, with a success only for the 1:3 internal resonance, on which this article has been focused. Experimentally, a clear exchange of energy between the lower mode to the higher one has been exhibited. However, the theoretical results did not predict such a clear energy exchange and explanations have been proposed to improve the model. In particular, PZT constitutive nonlinearities, not included in the model, could have an important effect on the nonlinear behavior and on the energy exchanges. In any cases, the measured energy exchange due to the 1:3 internal resonance could be improved by connecting the PZT patches to a proper nonlinear circuit, to enhance the nonlinear effects and achieve an efficient control of the beam.
## Acknowledgments

The Structural Mechanics and Coupled Systems Laboratory in Cnam, Paris, is warmly thanked for the bonding of the piezoelectric patches on the homogeneous beams, and especially F. Guillerm, B. Lossouarn, and J.-F. Deü. Moreover, the authors would like to thank the following organizations for supporting this research: (i) The "Ministère de la transition écologique et solidaire" and (ii) LABEX CELYA (ANR-10-LABX-0060) of the "Université de Lyon" within the program "Investissement d'Avenir" (ANR-11-IDEX-0007) operated by the French National Research Agency (ANR).

## Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

## Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

## ORCID iD

Vinciane Guillot https://orcid.org/0000-0002-6239-918X

## References

Abdelkefi A, Nayfeh AH and Hajj MR (2012) Effects of nonlinear piezoelectric coupling on energy harvesters under direct excitation. *Nonlinear Dynamics* 67(2): 1221–1232.

Amabili M (2008) *Nonlinear Vibrations and Stability of Shells and Plates*. USA: Cambridge University Press.

Anderson TJ, Balachandran B and Nayfeh AH (1994) Nonlinear resonances in a flexible cantilever beam. *Journal of Vibration and Acoustics* 116(4): 480–484.

Anderson TJ, Nayfeh AH and Balachandran B (1996) Coupling between high-frequency modes and a low-frequency mode: theory and experiment. *Nonlinear Dynamics* 11(1): 17–36.

Arquier R, Karkar S, Guillot L, et al. (2005–2019) *Manlab 4.0: An Interactive Path Following and Bifurcation Analysis Software*. France: CNRS. Available at: http://manlab.lma.cnrs-mrs.fr.

Bhugra H and Piazza G (2017) *Piezoelectric MEMS Resonators*. Switzerland: Springer.

Brand O, Dufour I, Heinrich S, et al. (2015) *Resonant MEMS: Fundamentals, Implementation, and Application*. Germany: Wiley-VCH.

Bricault C, Pézerat C, Collet M, et al. (2019) Multimodal reduction of acoustic radiation of thin plates by using a single piezoelectric patch with a negative capacitance shunt. *Applied Acoustics* 145: 320–327.

Cao DX, Leadenham S and Erturk A (2015) Internal resonance for nonlinear vibration energy harvesting. *The European Physical Journal Special Topics* 224(14–15): 2867–2880.

Collet M, Cunefare KA and Ichchou MN (2008) Wave motion optimization in periodically distributed shunted piezocomposite beam structures. *Journal of Intelligent Material Systems and Structures* 20(7): 787–808.

Crespo Da Silva MRM (1988) Non-linear flexural-flexural-torsional-extensional dynamics of beams-I. Formulation. *International Journal of Solids and Structures* 24(12): 1225–1234.

Czaplewski DA, Strachan S, Shoshani O, et al. (2019) Bifurcation diagram and dynamic response of a MEMS resonator with a 1:3 internal resonance. *Applied Physics Letters* 114(25): 254104.

Denis V, Jossic M, Giraud-Audine C, et al. (2018) Identification of nonlinear modes using phase-locked-loop experimental continuation and normal form. *Mechanical Systems and Signal Processing* 106: 430–452.

Dowell E and McHugh K (2016) Equations of motion for an inextensible beam undergoing large deflections. *Journal of Applied Mechanics* 83(5): 051007.

Ducarne J, Thomas O and Deü J-F (2012) Placement and dimension optimization of shunted piezoelectric patches for vibration reduction. *Journal of Sound and Vibration* 331(14): 3286–3303.

Emam SA and Nayfeh AH (2013) Non-linear response of buckled beams to 1:1 and 3:1 internal resonances. *International Journal of Non-Linear Mechanics* 52: 12–25.

Erturk A and Inman DJ (2011) *Piezoelectric Energy Harvesting*. UK: Wiley.

Garg A and Dwivedy SK (2019) Nonlinear dynamics of parametrically excited piezoelectric energy harvester with 1:3 internal resonance. *International Journal of Non-Linear Mechanics* 111: 82–94.

Ghayesh MH, Kazemirad S and Amabili M (2012) Coupled longitudinal-transverse dynamics of an axially moving beam with an internal resonance. *Mechanism and Machine Theory* 52: 18–34.

Guillot V, Ture Savadkoohi A and Lamarque C-H (2019) Analysis of a reduced-order nonlinear model of a multi-physics beam. *Nonlinear Dynamics* 97(2): 1371–1401.

Guyomar D, Aurelle N and Eyraud L (1997) Piezoelectric ceramics nonlinear behavior. Application to Langevin transducer. *Journal De Physique III* 7(6): 1197–1208.

Guyomar D, Ducharne B and Sebald G (2011) High nonlinearities in Langevin transducer: a comprehensive model. *Ultrasonics* 51(8): 1006–1013.

Houri S, Hatanaka D, Asano M, et al. (2019) Limit cycles and bifurcations in a nonlinear MEMS resonator with a 1:3 internal resonance. *Applied Physics Letters* 114(10): 103103.

Jacquelin E, Adhikari S and Friswell MI (2011) A piezoelectric device for impact energy harvesting. *Smart Materials and Structures* 20(10): 105008.

Lallart M (2017) Nonlinear technique and self-powered circuit for efficient piezoelectric energy harvesting under unloaded cases. *Energy Conversion and Management* 133: 444–457.

Lazarus A, Thomas O and Deü J-F (2012) Finite element reduced order models for nonlinear vibrations of piezoelectric layered beams with applications to NEMS. *Finite Elements in Analysis and Design* 49(1): 35–51.

Leadenham S and Erturk A (2015) Unified nonlinear electroelastic dynamics of a bimorph piezoelectric cantilever for energy harvesting, sensing, and actuation. *Nonlinear dynamics* 79(3): 1727–1743.

Lossouarn B, Deü J-F and Kerschen G (2018) A fully passive nonlinear piezoelectric vibration absorber. *Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences* 376(2127): 20170142.

Mam K, Peigney M and Siegert D (2017) Finite strain effects in piezoelectric energy harvesters under direct and parametric excitations. *Journal of Sound and Vibration* 389: 411–437.
Monteil M, Thomas O and Touzé C (2015) Identification of mode couplings in nonlinear vibrations of the steelpan. *Applied Acoustics* 89: 1–15.

Mook DT, Plaut RH and Haquang N (1985) The influence of an internal resonance on non-linear structural vibrations under subharmonic resonance conditions. *Journal of Sound and Vibration* 102(4): 473–492.

Nayfeh AH (1979) *Nonlinear Oscillations*. USA: Wiley-VCH.

Nayfeh AH and Balachandran B (1989) Modal interactions in dynamical and structural systems. *Applied Mechanics Reviews* 42(11S): S175–S201.

Özkaya E, Bağdatlı SM and Öz HR (2008) Nonlinear transverse vibrations and 3:1 internal resonances of a beam with multiple supports. *Journal of Vibration and Acoustics* 130(2): 021013.

Pai PF and Nayfeh AH (1990) Three-dimensional nonlinear vibrations of composite beams - I. Equations of motion. *Non-linear Dynamics* 1(6): 477–502.

Parashar SK and Von Wagner U (2004) Nonlinear longitudinal vibrations of transversally polarized piezoceramics: experiments and modeling. *Nonlinear Dynamics* 37: 51–73.

Preumont A (2011) *Vibration Control of Active Structures*. the Netherlands: Springer.

Ramini A, Hajjaj A and Younis MI (2016) Experimental investigation of 2:1 and 3:1 internal resonances in nonlinear MEMS arch resonators. In: Proceeding of the international design engineering technical conferences and computers and information in engineering conference, Charlotte, NC, Charlotte, NC: ASME.

Richard C, Guyomar D, Audigier D, et al. (1999) Semi-passive damping using continuous switching of a piezoelectric device. *Smart Structures and Materials 1999: Passive Damping and Isolation* 3672: 104–111.

Silva TMP, Clementino MA, De Marqui C Jr, et al. (2018) An experimentally validated piezoelectric nonlinear energy sink for wideband vibration attenuation. *Journal of Sound and Vibration* 437: 68–78.

Soltani P and Kerschen G (2015) The nonlinear piezoelectric tuned vibration absorber. *Smart Materials and Structures* 24(7): 075015.

Sun M, Quan T and Wang D (2018) Nonlinear oscillations of rectangular plate with 1:3 internal resonance between different modes. *Results in Physics* 11: 495–500.

Thomas O, Deü J-F and Ducarne J (2009) Vibrations of an elastic structure with shunted piezoelectric patches: efficient finite element formulation and electromechanical coupling coefficients. *International Journal for Numerical Methods in Engineering* 80(2): 235–268.

Thomas O, Sénéchal A and Deü J-F (2016) Hardening/softening behavior and reduced order modeling of nonlinear vibrations of rotating cantilever beams. *Nonlinear Dynamics* 86(2): 1293–1318.

Thomas O, Touzé C and Chaigne A (2003) Asymmetric non-linear forced vibrations of free-edge circular plates. part II: experiments. *Journal of Sound and Vibration* 265(5): 1075–1101.

Thomas O, Touzé C and Chaigne A (2005) Non-linear vibrations of free-edge thin spherical shells: modal interaction rules and 1: 1:2 internal resonance. *International Journal of Solids and Structures* 42(11–12): 3339–3373.

Thomas O, Touzé C and Luminais É (2007) Non-linear vibrations of free-edge thin spherical shells: experiments on a 1:1:2 internal resonance. *Nonlinear Dynamics* 49(1–2): 259–284.

Trindade M and Benjeddou A (2009) Effective electromechanical coupling coefficients of piezoelectric adaptive structures: critical evaluation and optimization. *Mechanics of Advanced Materials and Structures* 16(3): 210–223.

Von Wagner U and Hagedorn P (2002) Piezo-beam systems subjected to weak electric field: experiments and modelling of non-linearities. *Journal of Sound and Vibration* 256(5): 861–872.

Wolf K and Gottlieb O (2001) Nonlinear dynamics of a cantilever beam actuated by piezoelectric layers in symmetric and asymmetric configuration. In: Proceedings of the ASME design engineering technical conference, Pittsburgh, PA, 9–12 September 2001.

Zhang W and Guo X (2012) Periodic and chaotic oscillations of a composite laminated plate using the third-order shear deformation plate theory. *International Journal of Bifurcation and Chaos in Applied Sciences and Engineering* 22: 1202003.

Zhou B, Thouverez F and Lenoir D (2014) Essentially nonlinear piezoelectric shunt circuits applied to mistuned bladed disks. *Journal of Sound and Vibration* 333(9): 2520–2542.

# Appendix I

## Derivation of the theoretical model

This section is focused on the derivation of the model of equations (2) and (3) using the strong form of the equations of motions and a linear piezoelectric constitutive law. It is complementary to the derivation proposed in Guillot et al. (2019), where Lagrange equations are used. In essence, we combine Ducarne et al. (2012), in which linear models of piezoelectric (unimorph and bimorph) beams are proposed, and Thomas et al. (2016), in which a geometrically exact model for a cantilever beam is derived.

From Ducarne et al. (2012), we consider a laminated beam with K elastic (possibly piezoelectric, all built in the same material of constant $e_{31}$) layers whose axial direction is in the x-direction, with a transverse motion along the y-direction. There is no dependence on the z coordinate. The axial stress at a point $(x, y)$ of the beam in the k-th layer ($k = 1, \dots, K$) can be written

$$ \sigma_{xx} = Y_k(e - y\kappa) - e_{31}E \quad (11) $$

where $Y_k$ is the Young's modulus of the k-th layer and $E$ is the electric field. In any elastic layer, $e_{31} = 0$. From Thomas et al. (2016), with Euler-Bernoulli assumptions, the axial strain $e$, the curvature $\kappa$, the transverse displacement $v(x, t)$, the axial displacement $u(x, t)$, and the rotation $\theta(x, t)$ of the cross section are related by

$$ e = \sqrt{(1+u)^2 + v^2 - 1}, \quad \kappa = \theta' \quad (12) $$
$$ \sin \theta = \frac{v'}{1+e}, \quad \cos \theta = \frac{1+u'}{1+e} \qquad (13) $$

Then, the bending moment is, with equation (11)

$$ M = \int_S y\sigma_{xx}dS = Be - D\kappa + \sum_{p=1}^{P} \Theta_p V_p \qquad (14) $$

where S is the cross section and S is its area and $V_p$ is the voltage difference across the p-th piezoelectric layer and

$$ B = \sum_{k=1}^{K} b_k \frac{y_k^2 - y_{k-1}^2}{2} Y^{(k)}, \quad D = \sum_{k=1}^{K} b_k \frac{y_k^3 - y_{k-1}^3}{3} Y^{(k)}, \qquad (15) $$

$$ \Theta_p = b_p \frac{y_k + y_{k-1}}{2} e_{31} \qquad (16) $$

where $y_k$ is the vertical position of the interface between the $(k-1)$th and the $k$th layer. Those constants are, respectively, the axial/bending stiffness, the bending stiffness, and the piezoelectric coupling in bending.

Considering that the beam is inextensible, which is mandatory to obtain a simple nonlinear model as observed in the following, leads to $e = 0$. Consequently, the beam generalized constitutive law (14) reduces to

$$ M = -D\kappa + \sum_{p=1}^{P} \Theta_p V_p \qquad (17) $$

The nonlinear equations of motion of the beam, with no restriction on the cross-section rotation $\theta$ (geometrically exact beam model), with Euler-Bernoulli assumptions, and rotatory inertia neglected, are Thomas et al. (2016)

$$ \left\{ \begin{array}{l} (N \cos \theta - T \sin \theta)' = m\ddot{u}, \\ (N \sin \theta + T \cos \theta)' + p = m\ddot{v}, \\ T(1+e) + M' = 0. \end{array} \right. \qquad (18) $$

where $m$ is the mass per unit length of the beam. Then, the first of the aforementioned equations leads to

$$ N = \frac{1}{\cos \theta} \int_{L}^{x} m \ddot{u} dx + T \tan \theta \qquad (19) $$

By eliminating $N$ and $T$ using equation (19) and the last of equation (18) in the second one, and using the clamped/free boundary conditions, one obtains

$$ m\ddot{v} - \left( \frac{M'}{(1+e)\cos\theta} \right)' - \left( \tan\theta \int_{L}^{x} \rho S \ddot{u} dx \right)' = p \qquad (20) $$

Finally, using constitutive law (17), expanding all functions in Taylor series up to the third order in $v'$

($\theta = \arcsin v' \approx v' + v'^3/6$, $1/\cos \theta \approx 1 + v'^2/2$, and $\tan \theta \approx v'(1 - v'^2/2)$, using the inextensibility conditions (from (12) $e = 0 \Rightarrow u = \sqrt{1 - v'^2 - 1} \approx -v'^2/2$), and considering that $\Theta_p$ is zero out of $x \in [x_1, x_2]$, one obtains equation (2). The electric charge of equation (3) is obtained in the same manner, as explained in Ducarne et al. (2012).

A comment about these equations should be made. If a part of the beam has a nonsymmetric lamination (it is the case in the present study for the unimorph beam in the $x \in [x_1, x_2]$ area), $B \neq 0$ and equation (14) shows that an axial/bending coupling occurs. In the linear model of Ducarne et al. (2012), this effect is eliminated in the equations for a cantilever beam because $N=0$, which leads to a modified bending stiffness called $\hat{D}$. In the present nonlinear model, $N \neq 0$ because of the geometrical nonlinearities (see equation (19)), so that one cannot eliminate this axial/bending coupling of the equations without further assumption. It is the inextensibility condition ($e=0$), formulated to eliminate $u$ in the equations, which naturally cancels the linear axial/bending coupling in equation (14). However, the axial strain $e$ is clearly nonzero in the nonsymmetric area of the beam. This effect neglected in the present nonlinear model would probably be responsible of quadratic nonlinearities.

To give a first step to a more exact model, one has to compute $N$ by integrating equation (11)

$$ N = \int_S \sigma_{xx} dS = Ae - B\kappa + \sum_{p=1}^{P} \Xi_p V_p \qquad (21) $$

where $A = \sum_{k=1}^{K} b_k h_k Y_k$ and $\Xi_p = b_p e_{31}$. Eliminating the axial strain $e$ between (21) and (14) leads to

$$ M = \frac{B}{A} N - \hat{D}\kappa + \sum_{p=1}^{P} \hat{\Theta}_p V_p \qquad (22) $$

Then, eliminating $N$ between (22) and (19) leads to replacing (17) by

$$ M = \frac{B}{A} M' \tan \theta + \frac{B}{A \cos \theta} \int_{L}^{x} m \ddot{u} dx - \hat{D}\kappa + \sum_{p=1}^{P} \hat{\Theta}_p V_p \qquad (23) $$

where $\hat{D} = D - (B^2/A)$ and $\hat{\Theta}_p = \Theta_p - (\Xi B/A)$ are the modified constants that take into account the elimination of the axial strain, as introduced in Ducarne et al. (2012). A Taylor expansion of the functions of $\theta$ in $v'$ up to third order shows that the additional (underlined) terms are quadratic in $v'$. Then, even if the substitution of $M$ in (20) is not possible because (23) is now a differential equation in $M$, it shows that a nonsymmetric lamination in the beam adds quadratic terms in the transverse equation of motion.