# Energy Localization in Periodic Structures: Application to Centrifugal Pendulum Vibration Absorber

Aurelien Grolet, Alexandre Renault, Olivier Thomas

► To cite this version:

Aurelien Grolet, Alexandre Renault, Olivier Thomas. Energy Localization in Periodic Structures: Application to Centrifugal Pendulum Vibration Absorber. 17th International Symposium on Transport Phenomena and Dynamics of Rotating Machinery (ISROMAC2017), Dec 2017, Maui, United States. hal-02387666

HAL Id: hal-02387666

https://hal.science/hal-02387666

Submitted on 29 Nov 2019

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
# Energy Localization in Periodic Structures: Application to Centrifugal Pendulum Vibration Absorber

Aurélien GROLET¹*, Alexandre RENAULT¹, Olivier THOMAS¹

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

ISROMAC 2017

International Symposium on Transport Phenomena and Dynamics of Rotating Machinery

Maui, Hawaii

December 16-21, 2017

**Abstract**

In this paper we study the non-linear dynamic of centrifugal pendulum vibration absorbers (CPVA), and we pay a special attention to localized state solutions. The prediction of such states of vibration, and their stability, is of particular importance because they can lead to inefficient behavior for the CPVA and/or unforeseen stress levels.

Using an approximated equation for the pendulums dynamics, we derive initial conditions to put the system in localized states. Following an analytic study of the non-linear normal mode of the system, the resolution of the approximated equation is carried out in the frequency domain using the Harmonic Balance Method and the Asymptotic Numeric Method.

For sufficiently low inertia ratio, we show that the system can possess stable localized states.

**Keywords**

Nonlinear Dynamic — Centrifugal Pendulum Vibration Absorbers — Localization

Aurélien GROLET, Alexandre RENAULT, Olivier THOMAS
Arts et Métiers Paritech, LSIS UMR CNRS 7296, Lille, FRANCE

*Corresponding author: aurelien.grolet@ensam.eu

## INTRODUCTION

In this work we study the non-linear dynamic of centrifugal pendulum vibration absorbers (CPVA). Those components are often found in rotating machinery to suppress undesired rotational vibrations in the main shaft. For example, in the automotive industry they are used to suppress engine order vibration at the entrance of a gearbox [1]. They can also be used to reduce translational vibration of bladed rotor in helicopters [2]. This kind of absorber is very convenient because, due to centrifugal force, the vibration absorber stay tuned to a particular engine order, even while ramping up.

However, due to the different choices for the trajectory of the pendulums and due to the large amplitude of motion, the system is non-linear and the resonance frequencies depend on the amplitude. The behaviour can be hardening or softening depending on the choice of the trajectory for the pendulums [1]. Such non-linearity results in a completely new set of phenomenon as compared to the linear case. In particular, it can happen that the energy localizes to a particular subsystem, leaving the others virtually motionless [3, 4]. Such localization phenomenon can happen in perturbed (mis tuned) linear system (see e.g. [5]). But here, the structure is perfectly symmetric, and it is the nonlinearly that induces the localization through bifurcation of homogeneous solution. The prediction of such states of vibration, and their stability, is of particular importance because they can lead to unexpected and or inefficient behaviours for the CPVA and/or unforeseen stress levels. Such has already been done using analytical technique such as the (complex) averaging

method [4, 6].

Here, the computations are realized in a fully numerical framework in the frequency domain. The Harmonic Balance Method (HBM) and the Asymptotic Numeric Method (ANM) will be used for the continuation of periodic solutions of free (Nonlinear Normal Modes, NNM) and/or forced systems (forced response). Moreover the stability and bifurcation analysis is performed using Floquet's theory (in spectral domain). Everything is linked together and controlled with the GUI of the MANLAB tool developed under the direction of B.Cochelin [7].

In the application, we consider only two pendulums and we study the effect of the inertia ratio $\alpha$ on the bifurcation diagram of the free and forced response. We conduct a stability/bifurcation analysis of the free and forced response on a simplified version of the equation of motion. Finally those results are compared to the results of the full equation using time integration.

## 1. SYSTEM DEFINITION

We consider a structure composed of $N$ identical pendulums attached on a main rotating structure which is excited by a periodic torque $C(t)$. The rotating motion of the main structure is measured by the angle $\theta$ and the pendulum motions are measured by their curvilinear abscissa $s_i$, $1 \le i \le N$ (see Fig.1). In the following, we suppose that all pendulums are identical.
<figure><img src="image_3.png" /><figcaption>Figure 1. schematic view of CPVA and rotating structure</figcaption></figure>

## 1.1 Equation of motion

The equations of motion are obtained using Lagrange's equations. We denotes $Z_i = [X(s_i), Y(s_i)]^T$ the coordinates of the pendulum in the rotating frame. Note that the function $X$ and $Y$ are given by the trajectory of the pendulum.

### Kinetic energy

In order to express the kinetic energy, we define the following scalar quantities related to the pendulums trajectory (primes denote derivation w.r.t. the curvilinear abscissa):

$$F(s_i) = \dot{Z}_i^T Z_i = X(s_i)^2 + Y(s_i)^2 \quad (1)$$

$$H(s_i) = \dot{Z}_i^T Z_i' = X'(s_i)^2 + Y'(s_i)^2 \quad (2)$$

and

$$G(s_i) = Z_i^T \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix} Z_i = X(s_i)Y'(s_i) - X'(s_i)Y(s_i) \\ = \sqrt{FH} - \frac{1}{4}(F')^2 \quad (3)$$

Since the parameter $s_i$ is chosen to be the curvilinear abscissa, one has the simplification $H(s_i) = 1$, so that the kinetic energy of the i-th pendulum can be rewritten simply as follows:

$$T_i = \frac{1}{2} m [\dot{s}_i^2 + \dot{\theta}^2 F(s_i) + 2\dot{\theta}\dot{s}_i G(s_i)] \quad (4)$$

Summing the kinetic energies for all pendulums and adding also the kinetic energy of the rotor $T_r = \frac{1}{2} J \dot{\theta}^2$, one obtain the kinetic energy for the whole system:

$$T = \frac{1}{2} J \dot{\theta}^2 + \sum_{i=1}^{N} T_i \quad (5)$$

### Damping and External forces

We consider that each pendulum is subjected to linear damping with coefficient $\mu_p$. Finally, we consider that the rotor is subjected to linear damping with coefficient $\mu_r$, and to an external torque $C(t) = C_0 + C_1(t)$, where $C_0$ is constant and $C_1$ is periodic.

### Equations of motion

Applying Lagrange's equations leads to the following equation of motion for the n-th pendulum (for $n \in [1, N]$):

$$m\ddot{s}_n + mG_n\ddot{\theta} - \frac{1}{2}mF'_n\dot{\theta}^2 = -\mu_p\dot{s}_n \quad (6)$$

For the rotating structure, the equation of motion is:

$$\left( I + m \sum_{i} F_{i} \right) \ddot{\theta} + m \sum_{i} G_{i} \ddot{s}_{i} + m \sum_{i} (F'_{i} \dot{s}_{i} \dot{\theta} + G'_{i} \dot{s}_{i}^{2}) = -\mu_{r} \dot{\theta} + C(t) \quad (7)$$

## 1.2 Non-dimensional form of EoM

We introduce the following non-dimensional variables

$$\tau = \omega t, \quad \hat{s} = \frac{s}{\sqrt{F_0}}, \quad \hat{F} = \frac{F}{F_0}, \quad \hat{G} = \frac{G}{\sqrt{F_0}} \quad (8)$$

where $\sqrt{F_0}$ is a reference length (for example $F_0 = F(s=0)$). For now, $\omega$ is just a scaling parameter having the dimension of a rotation speed (later on, we will see that $\omega$ can be chosen as $\omega = C_0/\mu_r$).

We note that $\frac{d}{ds} = \frac{1}{\sqrt{F_0}} \frac{d}{d\hat{s}}$ and therefore we have

$$F' = \frac{dF}{ds} = \sqrt{F_0} \frac{d\hat{F}}{d\hat{s}} = \sqrt{F_0} \hat{F}' \quad (9)$$

and

$$G' = \frac{dG}{ds} = \frac{d\hat{G}}{d\hat{s}} = \hat{G}' \quad (10)$$

In addition, we introduce the inertia ratio $\alpha = \frac{mF_0}{I}$ which is supposed to be small. Substituting the non dimensional variables into the equation of motion leads to non dimensional form of the equation.

For the pendulums this leads to following non-dimensional equations (for $n \in [1, N]$):

$$\frac{d^2 \hat{s}_n}{d\tau^2} + \hat{G}_n \frac{d^2 \theta}{d\tau^2} - \frac{1}{2} \hat{F}'_n \left(\frac{d\theta}{d\tau}\right)^2 = -\hat{\mu}_p \frac{d \hat{s}_n}{d\tau} \quad (11)$$

where $\hat{\mu}_p = \frac{\mu_p}{m\omega}$.

For the rotor the non dimensional form of the equation is given by:

$$\left( 1 + \alpha \sum_{i} \hat{F}_{i} \right) \frac{d^2 \theta}{d\tau^2} + \alpha \sum_{i} \hat{G}_{i} \frac{d^2 \hat{s}_{i}}{d\tau^2} + (\hat{F}'_{i} \frac{d \hat{s}_{i}}{d\tau} \frac{d\theta}{d\tau} + \hat{G}'_{i} \frac{d \hat{s}_{i}}{d\tau})^2 = -\hat{\mu}_{r} \frac{d\theta}{d\tau} + \hat{C}(\tau) \quad (12)$$

where $\hat{\mu}_r = \frac{\mu_r}{\omega J}$ and $\hat{C}(\tau) = \frac{C(\tau)}{J\omega^2}$
### 1.3 Using $\theta$ as the independent variable

In the previous equations (Eq. (11), (12)), the independent variable is a non-dimensional time $\tau$. However, the fluctuation $C_1$ in the excitation torque is most of the time expressed in term of $\theta$. So it is useful to change the independent variable from $\tau$ to $\theta$. For this, we introduce a new variable $y$ (corresponding to a non-dimensional rotation speed) and we set the following notation:

$$ \dot{\theta} = \omega \frac{d\theta}{d\tau} = \omega y \quad (13) $$

We want to substitute $\theta$ in place of $\tau$ as the independent variable. Applying the chain rule for the derivation leads to the following:

$$ \frac{d}{d\tau} = \frac{d\theta}{d\tau} \frac{d}{d\theta} = \frac{dy}{d\theta} \frac{d}{d\theta} \text{ and } \frac{d^2}{d\tau^2} = y^2 \frac{d^2}{d\theta^2} + \frac{dy}{d\theta} \frac{d}{d\theta} \quad (14) $$

Substituting (14) into the non-dimensional equations (11) and (12) results in the following:

$$ y^2 \frac{d^2 \hat{s}_n}{d\theta^2} + y \frac{dy}{d\theta} \frac{d\hat{s}_n}{d\theta} + \hat{G}_n y \frac{dy}{d\theta} - \frac{1}{2} \hat{F}'_n y^2 = - \hat{\mu}_p y \frac{d\hat{s}_n}{d\theta}, \quad (15) $$

and

$$ \left( 1 + \alpha \sum_i \hat{F}_i \right) y \frac{dy}{d\theta} + \alpha \sum_i \hat{G}_i (y^2 \frac{d^2 \hat{s}_i}{d\theta^2} + y \frac{dy}{d\theta} \frac{d\hat{s}_i}{d\theta}) + y^2 \hat{F}'_i \frac{d\hat{s}_i}{d\theta} + \hat{G}'_i y^2 (\frac{d\hat{s}_i}{d\theta})^2 = - \hat{\mu}_r y + \hat{C}(\theta) \quad (16) $$

Finally, in order to simplify the notation, the derivative with respect to $\theta$ will be denoted by ()'. This leads to the following final form for the equations of motions:

$$ y \hat{s}''_n + y' \hat{s}'_n + \hat{G}_n y' - \frac{1}{2} \frac{d \hat{F}_n}{d \hat{s}_i} y = - \hat{\mu}_p \hat{s}'_n \quad (17) $$

$$ \left( 1 + \alpha \sum_i \hat{F}_i + \hat{G}_i \hat{s}'_i \right) y y' + \alpha y^2 \sum_i \hat{G}_i \hat{s}''_i + \frac{d \hat{F}_i}{d \hat{s}_i} \hat{s}'_i + \frac{d \hat{G}_i}{d \hat{s}_i} \hat{s}''_i = - \hat{\mu}_r y + \hat{C}(\theta) \quad (18) $$

Equations (17) and (18) constitute the starting point for the analysis carried out in this paper.

#### Equilibrium

Taking only a constant torque $C(t) = C_0$ leads to an equilibrium $s_i = 0, y' = 0$. Substituting into Eq. (18) equation leads to $\hat{\mu}_r y = \hat{C}$, so we can take $\omega = \frac{C_0}{\hat{\mu}_r}$ for the scaling parameter in the equation defining the non-dimensional time. At equilibrium, this will set the non-dimensional rotation speed to one, i.e. $y_{eq} = \hat{C}/\hat{\mu}_r = 1$.

#### Excitation

As already mentioned, the excitation term in Eq. (18) is composed of a constant torque $\hat{C}_0$ and a fluctuation $\hat{C}_1(\theta)$. In the remaining of the paper, we will consider that the fluctuation is harmonic, and that it pulsates $n_e$ times per round, so it can be expressed as $\hat{C}_1(\theta) = \hat{C}_1 \sin(n_e \theta)$.

### 1.4 Trajectory

Here we consider three types of trajectories for the pendulums, namely circular, cycloidal and epycycloidal. Here, one only gives the definition of the function $\hat{F}(s)$ in non-dimensional form (every needed quantity can be reconstructed from it, e.g. using Eq(??)). In the following $n_p$ is called the order of the trajectory and is related to geometric parameters only.

#### Circular trajectory

$$ \hat{F}(s) = \frac{1 + n_p^4 + 2n_p^2 \cos(s(1 + n_p^2))}{(1 + n_p^2)^2} \quad (19) $$

#### Epycycloidal trajectory

$$ \hat{F}(s) = 1 - n_p^2 s^2 \quad (20) $$

#### Cycloidal trajectory

$$ \hat{F}(s) = 1 - s^2 \left(\frac{3}{4} + n_p^2\right) + \frac{1}{4(1 + n_p^2)^2} \operatorname{asin}\left(s(1 + n_p^2)\right) \times \left[ \operatorname{asin}\left(s(1 + n_p^2)\right) + 2s(1 + n_p^2) \sqrt{1 - s^2(1 + n_p^2)^2} \right] \quad (21) $$

#### Taylor expansion of the trajectories

For linear and perturbation analysis, it will be useful to expand the trajectory in a Taylor series. For the three previously defined trajectories, the Taylor expansion can be written under the following form [8]:

$$ F(s) = 1 - n_p^2 s^2 + \gamma s^4 + O(s^6) \quad (22) $$

The coefficient $\gamma$ depends on the type of path:

* Circular path: $\gamma = \frac{1}{12} n_p^2 (n_p^2 + 1)^2$

* Epycycloidal path: $\gamma = 0$

* Cycloidal path: $\gamma = -\frac{1}{12} (n_p^2 + 1)^2$

### 1.5 Reduction

Here we present the reduction first introduced in [8] and reused in particular in [4]. The main idea is to use Taylor expansion, multiple scales and variable scaling to eliminate the variable $y$ from the equation of motion at first order in $\alpha$ (supposing that the inertia ratio $\alpha$ is small). The result is a set of N differential equations (one for each pendulums) that are independent of the rotor dynamic (at first order in $\alpha$). The resulting equations are particularly useful for an analytic analysis.

Following Shaw and al [8] the resulting equation of motion for the i-th pendulum is of the following form (for $i \in [1, N]$):

$$ z_i'' + \alpha \tilde{\mu}_p z_i' + n_p^2 z_i - 2\alpha\gamma z_i^3 + \alpha n_p^2 \sum_{j=1}^{N} z_j = -\alpha \tilde{C}_1(\theta) \quad (23) $$

We recall that ()' denotes derivative with respect to $\theta$. $z_i$ is the rescaled displacement of the i-th pendulum: $s_i = \alpha^{1/2} z_i$,
$\tilde{\mu}_p$ is the rescaled damping coefficient: $\hat{\mu}_p = \alpha \tilde{\mu}_p$, and $\tilde{C}_1$ is the rescaled periodic force amplitude: $\tilde{C}_1 = \alpha^{3/2} \tilde{C}_1$. The constant $\gamma$ is related to the trajectory (see Eq. (22)).

The resulting equation shows that each pendulum is weakly coupled to every other pendulums. At the order retain in the multiple scale expansion, the non-linearity comes only from the trajectory. One have $\gamma > 0$ for a circular trajectory (softening path) (resp. $\gamma < 0$ for an cycloidal trajectory, hardening path). For an epicycloidal trajectory, one have $\gamma = 0$ and Eq. (23) reduces to a linear equation.

## 2. NORMAL MODES OF THE SYSTEM

In this section we compute the non-linear modes of the system. Using the reduced system in Eq. (23) and a single harmonic HBM approach, we derive analytic conditions for bifurcation to occur. We compare those results to the non-linear normal modes of the full system in Eq. (17) and (18).

### 2.1 Linear modes

We start by computing the linear modes of the reduced equations (23). A simple analysis indicates that there is only two kind of modes:

* One in-phase mode, with shape $\Phi = [1, 1, \dots, 1]$ and angular frequency $\omega^2 = n_p^2(1 + N\alpha)$

* $N-1$ degenerated modes with shape of the type $\Phi = [0, \dots, 1, -1, \dots, 0]$ and angular frequency $\omega^2 = n_p^2$

We will focus on those modes, and we will show that bifurcation to localization is possible under special condition.

### 2.2 In phase non-linear mode and its bifurcation

Here we will search for an approximation of the backbone curve of the in-phase mdoe (and for possible bifurcation points) using the harmonic balance method [3] with only one harmonic. In order to simplify, we will retain only the cosine terms in the HBM development (i.e. we search for monophase non-linear normal modes). The i-th pendulum displacement is searched for under the following form:

$$z_i = A_i \cos \omega t \quad (24)$$

$A_i$ and $\omega$ are the unknown of the problem. Substituting this approximation into Eq. (23) and applying the harmonic balance method results in the following sets of algebraic equations ($\forall i \in [1, N]$):

$$(n_p^2 - \omega^2)A_i - \frac{3}{2}\alpha\gamma A_i^3 + \alpha n_p^2 \sum_{j=1}^{N} A_j = 0 \quad (25)$$

The in phase mode is defined by the fact that all pendulums are vibrating with the same amplitude in a “unisson” way. Thus we assume $A_i = A$ and we substitute into Eq. (25)

to obtain an approximation of the frequency-amplitude relation of the in-phase mode:

$$\omega^2 = n_p^2(1 + N\alpha) - \frac{3}{2}\alpha\gamma A^2 \quad (26)$$

The bifurcation (branching) points are obtained when the determinant of the jacobian matrix of the algebraic system is equal to zero. Here the jacobian matrix along the in-phase mode ($A_i = A$ for $i \in [1, N]$) is given by the following:

$$J_{\text{in-phase}} = \begin{bmatrix} \xi & \mu & \dots & \dots & \mu \\ \mu & \xi & \ddots & (\mu) & \vdots \\ \vdots & \vdots & \ddots & \ddots & \vdots \\ \vdots & (\mu) & \ddots & \xi & \mu \\ \mu & \dots & \dots & \mu & \xi \end{bmatrix}, \quad (27)$$

with $\xi = n_p^2(1 + \alpha) - \omega^2 - \frac{9}{2}\alpha\gamma A^2$, and $\mu = n_p^2\alpha$. The eigenvalues of $J$ are $\xi + (N-1)\mu$ with multiplicity one (resp. $\xi - \mu$ with multiplicity $N-1$) and the associated eigenvector is the vector full of ones (resp. a series of $N-1$ vectors of the type $[0, \dots, 0, 1, -1, 0, \dots, 0]$). Therefore, the determinant of the jacobian matrix along the first mode is easily obtained as (product of the eigenvalues):

$$\det(J_{\text{in-phase}}) = (\xi + (N-1)\mu) \cdot (\xi - \mu)^{N-1} \quad (28)$$

The determinant is null if one of its factor is equal to zero. Combining Eq. (28) with the expression of frequency-amplitude relation for the in phase mode (Eq. (26) leads to two sets of equation defining the bifurcation points of the in-phase mode.

The first set of equations gives the bifurcation point from the trivial solution: $A = 0$, $\omega^2 = n_p^2(1 + N\alpha)$. The second set of equations gives a bifurcation point along the in-phase mode, it is characterized by the following amplitude and angular frequency:

$$A^2 = -n_p^2 \frac{N}{3\gamma}, \quad \omega^2 = n_p^2 \left(1 + \frac{3}{2}N\alpha\right) \quad (29)$$

There is several things to notice from this simple analysis. First the bifurcation point of the in-phase mode exists only for $\gamma < 0$, i.e for a hardening path (cycloidal trajectory). Second the position of the bifurcation point is directly linked to inertia ratio $\alpha$ and to the number of pendulums. The higher $\alpha$, the further the bifurcation point, and the more pendulums the further the bifurcation point. Interestingly enough, the amplitude at the bifurcation point does not depend on $\alpha$, but depends on the number of pendulum.

### 2.3 Degenerated non-linear modes and its bifurcation

As stated before, there exists $N-1$ degenerated linear modes. Lets consider the particular modes with shape $\Phi = [1, -1, 0, \dots, 0]$, (in other words $A = A_1 = -A_2$ and $A_i = 0$ for $i \neq 1, 2$). Substituting into Eq. (25) leads to an approximation of the
frequency amplitude relationship for the degenerated modes:

$$ \omega^2 = n_p^2 - \frac{3}{2}\gamma\alpha A^2 \quad (30) $$

The jacobian matrix along the degenerated mode ($A = A_1 = -A_2$ and $A_i = 0$ for $i \neq 1, 2$) is given by the following:

$$ J_{\text{deg.}} = \begin{bmatrix} \xi & \mu & \dots & \dots & \dots & \mu \\ \mu & \xi & \mu & & & \vdots \\ \vdots & \mu & \xi_l & \mu & & \vdots \\ \vdots & & \ddots & \ddots & \ddots & \vdots \\ \vdots & & & \mu & \xi_l & \mu \\ \mu & \dots & \dots & \dots & \mu & \xi_l \end{bmatrix}, \quad (31) $$

with $\xi = n_p^2(1 + \alpha) - \omega^2 - \frac{2}{3}\alpha\gamma A^2$, $\xi_l = n_p^2(1 + \alpha) - \omega^2$ and $\mu = n_p^2\alpha$.

We suppose that $N > 2$. One notes that $J$ has an eigenvalue given by $\xi - \mu$ (associated with the eigenvector $\Phi = [1, -1, 0, \dots, 0]$). Another eigenvalue is $\xi_l - \mu$, associated with the eigenvector $\Phi = [0, 0, 1, -1, 0, \dots, 0]$. There is $N - 3$ such eigenvectors, so this eigenvalue is with multiplicity $N - 3$.

The last two eigenvalues are associated with eigenvectors of the form $[1, 1, X, \dots, X]$ where $X$ is solution of the following second degree equation:

$$ (N - 2)\mu X^2 + (\xi - \xi_l - (N - 2)\mu)X - 2\mu = 0 \quad (32) $$

The corresponding eigenvalues are given by:

$$ \lambda_{\pm} = \xi + \mu + (N - 2)\mu X_{\pm} \quad (33) $$

Finally, one can express the determinant of $J_{\text{deg.}}$ as follows (product of the eigenvalues):

$$ \det(J_{\text{deg.}}) = (\xi - \mu) \cdot (\xi_l - \mu)^{N-3} \cdot (\xi + \mu + (N-2)\mu X_+) \cdot (\xi + \mu + (N-2)\mu X_-) \quad (34) $$

Following the same procedure as in the previous subsection, one is able to derive the position of the bifurcation points along the degenerated modes. Discarding the bifurcations from the trivial solution, one find that the condition for a non-trivial bifurcation is:

$$ A_{\pm}^2 = n_p^2 \frac{2 + (N - 2)X_{\pm}}{3\gamma} \quad (35) $$

This analysis suggest that there can exist two non trivial bifurcation points for the degenerated modes. The existence of the bifurcation points is linked to the sign of the left hand side. In particular if both solutions $X_{\pm}$ are positives, then those bifurcation points occur only for $\gamma > 0$, i.e for softening trajectories (circular). If the bifurcation point exists, the associated angular frequency is given by $\omega^2 = n_p^2[1 - \frac{\alpha}{2}(2 - (N-2)X_{\pm})]$. Note that for the particular case $N = 2$, the formula is still valid if we set $X_+ = X_- = 0$.

## 2.4 Detailed application to a system of two pendulums

In the case of only two pendulums, Eq. (25) reduces to the following:

$$ (n_p^2 - \omega^2)A_1 - \frac{3}{2}\gamma\alpha A_1^3 + \alpha(A_1 + A_2) = 0 \quad (36) $$

$$ (n_p^2 - \omega^2)A_2 - \frac{3}{2}\gamma\alpha A_2^3 + \alpha(A_1 + A_2) = 0 \quad (37) $$

Taking the sum and difference of the two previous equation and factoring we have:

$$ (A_1 - A_2) \left[ n_p^2 - \omega^2 - \frac{3}{2}\gamma\alpha(A_1^2 + A_2^2 + A_1 A_2) \right] = 0 \quad (38) $$

$$ (A_1 + A_2) \left[ n_p^2 + 2n_p^2\alpha - \omega^2 - \frac{3}{2}\gamma\alpha(A_1^2 + A_2^2 - A_1 A_2) \right] = 0 \quad (39) $$

From there we derive four solutions:

* $A_1 = A_2$ and $A_1 = -A_2$, so that $A_1 = A_2 = 0$ and $\omega \in \mathbb{R}$, this is the trivial solution

* $A_1 = -A_2$ and $\omega^2 = n_p^2 - \frac{3}{2}\gamma\alpha A_1^2$. This is the out of phase mode (not degenerated in this particular case of two pendulums).

* $A_1 = A_2$ and $\omega^2 = n_p^2(1 + 2\alpha) - \frac{3}{2}\gamma\alpha A_1^2$. This is the in phase mode.

* $A_1 A_2 = -\frac{n_p^2}{3\gamma}$ and $\omega^2 = n_p^2 - \frac{3}{2}\gamma\alpha(A_1^2 + A_2^2 + A_1 A_2)$. There is two cases to consider:

- $\gamma < 0$ (cycloid path, hardening): This solution is a bifurcation of the in-phase mode that corresponds to a localized solution. The bifurcation point is given by and $A_1^2 = A_2^2 = -\frac{n_p^2}{3\gamma}$ with angular frequency $\omega^2 = n_p^2(1 + 3\alpha)$.

- $\gamma > 0$ (circular path, softening): This solution is a bifurcation of the out-of-phase mode that corresponds to a localized solution. The bifurcation point is given by and $A_1^2 = A_2^2 = \frac{n_p^2}{3\gamma}$ with angular frequency $\omega^2 = n_p^2(1 - \alpha)$.

## 3. SOLUTION METHODS

In the next section we will search for localized states of the system by first solving the approximated equation in Eq. (23). The resolution of the equation is done using the harmonic balance method (HBM) (see e.g. [3]). This method consists in searching for periodic solutions under the form of a truncated Fourier series. This allows to transform the set of differential equation into a set of algebraic equation that can then be solved for example using newton-like algorithm.

Here the continuation of the periodic solution will be done using the Numeric Assymptotic Method (MAN) in the framework of the MANLAB software [7]. In this method each
variable is searched for under a power series of a pseudo-curvilinear abscissa.

The stability of the periodic solutions is computed using Floquet theory in spectral domain [9], and bifurcated branches are computed by adding a small perturbation in the vicinity of a branching point, detected my monitoring the Floquet multiplier w.r.t the unit circle.

Finally verification is carried out using classic time integration scheme on the complete equations in Eq. (17) and (18)

## 4. APPLICATION TO AN EXAMPLE

In the remaining of the paper, we will focus on a particular example to show the validity of the previous analysis. The numerical values to be used are given by the following:

$$n_p = 1, \mu_p = 0.02, \mu_r = 0.005 \qquad (40)$$

Only a cycloidal path will be considered, but the same analysis could be done for circular trajectories. In a first subsection, we will set the inertia ratio to $\alpha = 1/20$. Then we will study the influence of $\alpha$ on the results.

### 4.1 Particular case of $N = 2$ pendulums with inertia ratio $\alpha = 1/20$

**Non-linear normal mode**

In this section we compute the non-linear normal modes of a system with two pendulums from the undamped and unforced version of Eq. (23). Ten harmonics are retained for the HBM approximation, and the threshold on the norm of Floquet multipliers is set to 0.001 (i.e a solution is stable if the norm of all Floquet multipliers is less than $1 + 0.001$). The results are depicted on Fig.2, along with the results of the analytic analysis of section 2.4.

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Backbone curves for the non-linear modes computed from Eq. (23), □: stable, ×: unstable. Comparison with analytic results from section 2.4 (·).</figcaption></figure>

The numerical results indicate that the in-phase mode loses its stability through a pitchfork bifurcation and the stability is transferred to a solution corresponding to a localized

state. As expected, the stability change happens exactly at the bifurcation point. Finally, one can see that we have an excellent agreement between the analytic and numeric solutions for both modes and for the bifurcated solution.

**Forced response and bifurcation**

In this section we compute the forced response of a system with two pendulums for different forcing amplitude from Eq. (23), and we compare the results with time integration of the complete equations in Eq. (17) and Eq. (18).

The parameter $n_e$ will be varied around $n_p$. Strictly speaking, in the application, $n_e$ is given by the system (for example $n_e = 2$ for a 4 cylinders engine). However, it may happen that the pendulums are not perfectly tuned to the excitation. This mistuning may be intentional or may results from material imperfections. Therefore varying the excitation frequency is equivalent to introducing mistuning (provided that all pendulums have the same level of mistuning) and is relevant for studying the effect of mistuning on the system response.

For low amplitude forcing ($\hat{C}_1 = 0.005$) the response is similar to the response of a linear system, where the in-phase solution is stable for all excitation frequencies (Fig.3).

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Forced response for $\hat{C}_1 = 0.005$ computed from Eq. (23), □: stable, ×:unstable. Comparison with time integration of Eq. (17) and (18) (·)</figcaption></figure>

When the forcing increases ($\hat{C}_1 = 0.01$), the in phase solution loses its stability and a stable localized state appears (Fig.4). The loss of stability and the localization is confirmed by time integration of the complete equations, see Fig.4 and Fig.5. According to the literature, some authors consider localized motion only when the energy of the subsystem is at least one order of magnitude smaller than the energy of the remaining part. Here this definition will be relaxed a bit, and we consider localization when the response is not homogeneous anymore, i.e. when the response lies on a bifurcated branch of solution (bifurcation from the homogeneous solution). Note that both definitions coincide when the amplitude (along the bifurcated branch) increases.

When the forcing is even larger ($\hat{C}_1 = 0.02$) the localized
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Forced response for $\tilde{C}_1 = 0.01$ computed from Eq. (23), □: stable, ×:unstable. Comparison with time integration of Eq. (17) and (18) (○)</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Time evolution of a localized state computed by time integration of Eq. (17) and (18) for $\tilde{C}_1 = 0.01$ and $\omega = 1.07$</figcaption></figure>

solution can lose its stability through Neimark-Sacker (Hopf-2) bifurcation, leading to localized quasi-periodic states (Fig.6). An example of such quasi periodic solution obtained by time integration of the complete equation is given on Fig.7.

In this case of low inertia ratio ($\alpha = \frac{1}{20}$), the simple equation in Eq. (23) is sufficient to have a very good approximation of the forced solution as demonstrated by the previous simulation. Moreover, the bifurcation points and the loss of stability are also computed with good agreement as compared to time integration of the complete system. Finally it is interesting to note that Eq. (23) is again sufficient to detect the Neimark-Sacker bifurcation range with a very good level of accuracy.

As a partial conclusion, one can see that adding too much detuning to the pendulums can lead to a drastic change in the topology of the solution. Indeed, if the level of detuning

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Forced response for $\tilde{C}_1 = 0.02$ computed from Eq. (23), □: stable, ×:unstable. Comparison with time integration of Eq. (17) and (18) (○)</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Quasi periodic solution obtained from time integration of Eq. (17) and (18) for $\tilde{C}_1 = 0.02$ and $\omega = 1.08$ (indicated by a black square on Fig.6)</figcaption></figure>

is such that the solution amplitude overpass the bifurcation point of the in-phase mode, then localized states can occur in the system.

## 4.2 Influence of the inertia ratio $\alpha$

In this section we study the influence of the inertia ratio on the previous results. As indicated by the analytic study in section 2.4, the position of the bifurcation point of the in-phase mode is directly linked to the ratio $\alpha$ through the following formula:

$$A_1^2 = A_2^2 = -n_p^2 \frac{2}{3\gamma}, \quad \text{and} \quad \omega^2 = n_p^2(1 + 3\alpha) \quad (41)$$

As already mentioned, only the frequency at which the bifurcation occur depends on $\alpha$, the amplitude does not.

Four inertia ratio will be considered: $\alpha = \frac{1}{20}, \frac{1}{10}, \frac{1}{5}$ and $\frac{1}{3}$.
the forcing amplitude is fixed to $\hat{C}_1 = 0.02$. Results for the forced response are depicted on Fig. 6, 8,9, and 10.

One can observe that when $\alpha$ increases the resonance is shifted to the right as expected from the non-linear mode analysis. Moreover, the branch corresponding to the localized solution (the “croissant”) seems to shrink as $\alpha$ increases, until a point where no bifurcation occurs (Fig. 10). Actually, it seems that the amplitude of the solution (normalized with respect to $\alpha$) seem to decrease when $\alpha$ increases, thus one does not attain the amplitude of the bifurcation point any more, thus no localized solution. In order to check the results, time integration were performed for different frequencies. For $\alpha = \frac{1}{10}$, it is still possible to observe localized solution (just between the branching point and the Neimark-Sacker bifurcation), see Fig. 8. For higher $\alpha$ we were not able to confirm the localized branch by time integration. Actually, the amplitude of the solution obtained through the simplified equations (Eq. (23)) is too high to be considered as physical amplitudes. Indeed, for a cycloidal trajectory of the pendulums, there is an upper limit for the amplitude of the curvilinear abscissa (otherwise some complex quantities appear due to square root of negative numbers). Here the limit is $s^2 < \frac{1}{(1+n_p^2)^2}$, this also corresponds to $z^2 < \frac{1}{\alpha(1+n_p^2)^2}$. In our example, for $\alpha = \frac{1}{5}$ we should have $z < \sqrt{\frac{5}{4}} \approx 1.12$, which is why one cannot obtain time integration results (black circle) above that amplitude on Fig. 9.

Finally, one can see on Fig. 10 that Eq. (23) is valid only for $\alpha < \frac{1}{3}$. This is in agreement with the hypothesis used to derive this equation (e.g. $\alpha$ was considered small). In this case, the maximum physical amplitude is given by $z < \sqrt{\frac{3}{4}} \approx 0.87$, which is why no time integration results can be obtained over this amplitude.

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Forced response for $\alpha = \frac{1}{10}$ and $\hat{C}_1 = 0.02$ computed from Eq. (23), □: stable, ×:unstable. Comparison with time integration of Eq. (17) and (18) (○)</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 9.** Forced response for $\alpha = \frac{1}{5}$ and $\hat{C}_1 = 0.02$ computed from Eq. (23), □: stable, ×:unstable. Comparison with time integration of Eq. (17) and (18) (○)</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Forced response for $\alpha = \frac{1}{3}$ and $\hat{C}_1 = 0.02$ computed from Eq. (23), □: stable, ×:unstable. Comparison with time integration of Eq. (17) and (18) (○)</figcaption></figure>

## 5. CONCLUSION

In this work we study the non-linear dynamic of centrifugal pendulum vibration absorbers (CPVA), and we paid a special attention to localized solutions. The prediction of such states of vibration, and their stability, is of particular importance because they can lead to inefficient behavior for the CPVA and/or unforeseen stress levels.

Using an approximated equation for the pendulums dynamic, one were able to derive initial conditions to put the system in localized states. Following an analytical study of the non-linear normal mode of the system, the resolution of the approximated equation were realized in the frequency domain using the Harmonic Balance Method (HBM) and the Asymptotic Numeric Method (ANM). stability was computed using Floquet theory. Verification of the HBM solution was carried out using time integration of the complete equation.
For sufficiently low inertia ratio, we show that the system
can possess stable localized states. An interesting results is
that localized quasi periodic solutions can also exist through
Neimark-Sacker bifurcation of the periodic localized states
for sufficiently high forcing amplitude.

Further work will address the effect of adding more pendulums to the system.

REFERENCES

[1] A. Renault, O. Thomas, H. Mahé, and Y. Lefebvre. Hardening softening behavior of antiresonance for nonlinear torsional vibration absorbers. *proceedings of XXIV ICTAM, Montreal, 2016*, 2016.

[2] C. Shi, R. Parker, and S. Shaw. Tuning of centrifugal pendulum vibration absorbers for translational and rotational vibration reduction. *Mechanism and Machine Theory*, 66:56-65, 2013.

[3] A. Grolet and F. Thouverez. Computing multiple periodic solutions of nonlinear vibration problems using the harmonic balance method and groebner bases. *Mechanical Systems and Signal Processing*, 52-53:529-547, 2015.

[4] J. Issa and S. Shaw. Synchronous and non-synchronous responses of systems with multiple identical nonlinear vibration absorbers. *Journal of Sound and Vibration*, 348:105-125, 2015.

[5] A. Alsuwaiyan and S. Shaw. Steady-state responses of systems of nearly-identical torsional vibration absorbers. *Journal of Vibration and Acoustics*, 125:80-87, 2003.

[6] A. Alsuwaiyan and S. Shaw. Non-synchronous and localized responses of systems of identical centrifugal pendulum vibration absorbers. *Arabian Journal for Science and Engineering*, 39:9205-9217, 2014.

[7] R. Arquier, S. Karkar, A. Lazarus, O. Thomas, C. Vergez, and B. Cochelin. Manlab 2.0 : an interactive path-following and bifurcation analys is software. *Tech. rep.*, Laboratoire de Mécanique et d'Acoustique.

[8] A. Alsuwaiyan and S. Shaw. Perfomance and dynamic stability of general path centrifugal pendulum absorbers. *Journal of Sound and Vibration*, 252:791-815, 2002.

[9] A. Nayfeh and B. Balachandran. *Applied nonlinear dynamic*. Wiley series in nonlinear science, 1995.