# Modal Analysis of a Nonlinear Periodic Structure with Cyclic Symmetry

Fostios Georgiades, Maxime Peeters, Gaëtan Kerschen, Jean-Claude Golinval, Massimo Ruzzene

► To cite this version:

Fostios Georgiades, Maxime Peeters, Gaëtan Kerschen, Jean-Claude Golinval, Massimo Ruzzene. Modal Analysis of a Nonlinear Periodic Structure with Cyclic Symmetry. AIAA Journal, 2009, 47, pp.1014 - 1025. 10.2514/1.40461. hal-01385735

HAL Id: hal-01385735

https://hal.science/hal-01385735

Submitted on 22 Oct 2016

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
# Modal Analysis of a Nonlinear Periodic Structure with Cyclic Symmetry

F. Georgiades,\* M. Peeters,† G. Kerschen,‡ and J. C. Golinval§

University of Liège, B-4000 Liège, Belgium

and

M. Ruzzene¶

Georgia Institute of Technology, Atlanta, Georgia 30332

The objective of this paper is to examine nonlinear normal modes and their bifurcations in cyclic periodic structures. The nonlinear normal modes are computed using a numerical technique that combines shooting and pseudoarclength continuation. Unlike perturbation techniques, the resulting algorithm can investigate strongly nonlinear regimes of motion. This study reveals that modal interactions may occur without necessarily having commensurate natural frequencies in the underlying linear system. In addition, a countable infinity of such modal interactions are shown to exist in the system.

## Nomenclature

**F** = augmented two point boundary value problem

**f**<sub>nl</sub> = nonlinear restoring force

**g** = vector field

**H** = shooting function

**h** = phase condition

**K** = stiffness matrix

**k** = linear spring

**k**<sub>nl</sub> = nonlinear spring

**M** = mass matrix

**M** = disk mass

**m** = blade mass

**p** = tangent vector

**s** = predictor step size

**T** = period

**t** = time

**X**<sub>i</sub> = ith disk mass displacement

**x**<sub>i</sub> = ith blade mass displacement

**x**(**t**) = displacement

**ẋ**(**t**) = velocity

**ẋ̇**(**t**) = acceleration

**z**(**t**) = state space vector

**Δ****T** = correction to the period during shooting

**Δ****z**<sub>p0</sub> = correction to the initial conditions during shooting

**ε** = prescribed relative precision

## I. Introduction

MODAL analysis refers to the analysis of the dynamics of a vibrating structure in terms of its modal parameters. For linear systems, these parameters, termed the mode shapes, natural frequencies, and damping ratios, can be computed from the

governing equations of motion merely by solving an eigenvalue problem. A number of methods that can extract them directly from the system response have also been developed [1]. Modal analysis and testing of linear mechanical structures has been developed over the past 40–50 years, and the techniques available today are really quite sophisticated and advanced.

Modal analysis proved useful for analyzing the dynamics of linear periodic or nearly periodic structures (i.e., structures composed of a number of similar or identical elements coupled together, practical examples of which are bladed disk assemblies and space antennas). In particular, energy localization phenomena occurring in slightly mistuned periodic structures were shown to exist due the presence of strongly localized modes [2–4]. Realizing that nonlinearity is a frequent occurrence in real life applications, the dynamics of nonlinear periodic structures were also investigated in several studies. Wei and Pierre examined the effects of dry friction on a nearly cyclic structure using the harmonic balance method [5]. In a series of papers, Vakakis et al. demonstrated that, in contrast to the findings of linear theory, *nonlinear mode localization may occur in perfectly cyclic nonlinear systems* [6–9]. Other studies dealing with mode localization in nonlinear cyclic systems are [10–12].

Because any attempt to apply traditional linear analysis to nonlinear systems is bound to fail, the nonlinear normal mode (NNM) theory is used herein. NNMs and asymptotic analysis were used in [6–9] to examine the free and forced vibrations of nonlinear periodic systems with cyclic symmetry. In the present study, we revisit the free dynamics of such systems using the numerical algorithm for NNM computation proposed in [13,14]. This computational framework allows the relaxation of the assumptions of weak nonlinearity and/or weak coupling used in [6–9]. In addition, it is shown that modal interactions can take place in nonlinear systems without necessarily having commensurate natural frequencies in the underlying linear system. This dynamical phenomenon is rarely discussed in the literature.

This paper is organized as follows. In the next section, the two main definitions of an NNM are provided. The fundamental properties of NNMs are also described, and a numerical algorithm for their computation is introduced. In Sec. III, the NNM theory is applied to a nonlinear cyclic system, and the resulting dynamics are analyzed in detail. The conclusions are summarized in Sec. IV.

## II. Nonlinear Normal Modes

The free response of discrete mechanical systems is considered, assuming that continuous systems have been spatially discretized using, for example, the finite element method. The equations of motion are

<footnote>\*Postdoctoral Researcher; fotios.georgiades@ulg.ac.be.</footnote>

<footnote>†Ph.D. Student; m.peeters@ulg.ac.be.</footnote>

<footnote>‡Assistant Professor, Space Structures and Systems Lab, Structural Dynamics Research Group Department of Aerospace and Mechanical Engineering, 1 Chemin des Chevreuil; g.kerschen@ulg.ac.be.</footnote>

<footnote>§Professor; jc.golinval@ulg.ac.be.</footnote>

<footnote>¶Associate Professor, School of Aerospace Engineering; massimo.ruzzene@ae.gatech.edu.</footnote>
$$M\ddot{x}(t) + Kx(t) + f_{nl}\{x(t)\} = 0 \quad (1)$$

where M is the mass matrix; K is the stiffness matrix; x and $\ddot{x}$ are the displacement and acceleration vectors, respectively; and $f_{nl}$ is the nonlinear restoring force vector.

### A. Definitions

NNMs were pioneered in the 1960s thanks to the seminal work of Rosenberg [15]. They were regarded as a theoretical curiosity until the beginning of the 1990s, when they were given a new impetus through the efforts of Vakakis et al. [16 18] and Shaw and Pierre [19,20].

Two main definitions of an NNM exist in the literature. During the normal mode motion of a linear conservative system, each system component moves with the same frequency and with a fixed ratio amongst the displacements of the components. Targeting a straightforward nonlinear extension of the linear normal mode (LNM) concept, Rosenberg defined an NNM as a vibration in unison of the system (i.e., a synchronous oscillation). This definition requires that all material points of the system reach their extreme values and pass through zero simultaneously and allows all displacements to be expressed in terms of a single reference displacement. One fundamental difference with linear systems is that the modal lines of nonlinear systems are generally curves, resulting from the nonlinear relationship between the coordinates during the periodic motion.

Shaw and Pierre proposed a generalization of Rosenberg's definition that provides a direct and elegant extension of the NNM concept to damped systems [20]. Based on geometric arguments and inspired by the center manifold technique, they defined an NNM as a two dimensional invariant manifold in phase space. Such a manifold is invariant under the flow (i.e., orbits that start out in the manifold remain in it for all time), which extends the invariance property of LNM to nonlinear systems.

At first glance, Rosenberg's definition may appear restrictive in two cases:

1) It cannot be easily extended to nonconservative systems.

2) In the presence of internal resonances (i.e., when two or more NNMs interact), some coordinates may have a dominant frequency component that is different than that of the other coordinates (e.g., some coordinates may vibrate faster than others). In this case, the system no longer vibrates in unison.

However, these two limitations can be circumvented. First, as shown in [13,21], the damped dynamics can often be interpreted based on the topological structure and bifurcations of the NNMs of the underlying undamped system. Secondly, realizing that the motion is still periodic in the presence of internal resonances, Rosenberg's definition of an NNM can be extended to a (nonnecessarily synchronous) periodic motion of the system. This extended definition is particularly attractive when targeting a numerical computation of the NNMs. As will be shown in Sec. II.C, the definition also enables the nonlinear modes to be effectively computed using advanced algorithms for the continuation of periodic solutions. This NNM definition is considered throughout the present study.

### B. Fundamental Properties

A detailed description of NNMs and their fundamental properties is beyond the scope of this paper. However, for completeness, the frequency energy dependence, as well as bifurcations and stability of NNMs, are briefly reviewed in this section. For further details, the reader can refer to [13,17,18].

A typical dynamical feature of nonlinear systems is the frequency energy dependence of their oscillations. One important consequence is that the frequency response functions (FRFs) of nonlinear systems are no longer invariant. The modal curves and frequencies of oscillation of NNMs also depend on the total energy in the system. In contrast to linear theory, this energy dependence prevents the direct separation of space and time in the governing equations of motion, which complicates the analytical calculation of the NNMs.

Because of frequency energy dependence, an appropriate graphical depiction of the NNMs is key to their exploitation. In this study, the representation of the NNMs of the underlying Hamiltonian system in a frequency energy plot (FEP) is adopted. An NNM is represented by a point in the FEP, which is drawn at a frequency corresponding to the minimal period of the periodic motion and at an energy equal to the conserved total energy during the motion. A branch, represented by a solid line, corresponds to a family of NNM motions possessing the same qualitative features (e.g., the in phase NNM motions of a 2 degree of freedom (DOF) system). For illustration, a 2 DOF system, the motion of which is governed by the equations

$$\ddot{x}_1 + (2x_1 - x_2) + 0.5x_1^3 = 0 \qquad \ddot{x}_2 + (2x_2 - x_1) = 0 \quad (2)$$

is considered, and its FEP is represented in Fig. 1. The backbone of the plot is formed by two branches, which represent in phase ($S11+$) and out of phase ($S11-$) synchronous NNMs. The FEP clearly shows that the nonlinear modal parameters have a strong dependence on the total energy in the system:

1) The frequency of both the in phase and out of phase NNMs increases with the energy level, which reveals the hardening characteristic of the system.

2) For increasing energies, the in phase NNM tends to localize to the second DOF (i.e., it resembles a vertical curve), whereas the out of phase NNM localizes to the first DOF (i.e., it resembles an horizontal curve).

A second salient feature of nonlinear systems is that NNMs may interact during a general motion of the system. An energy exchange between the different modes involved may therefore be observed during the internal resonance. For instance, exciting a high frequency mode may produce a large amplitude response in a low frequency mode. Internally resonant NNMs have no counterpart in linear systems and are generated through bifurcations. The FEP in Fig. 1 does not seem to feature internally resonant NNMs. At higher energies, however, another branch of periodic solutions, termed a tongue, emanates from the backbone branch $S11+$ (Fig. 2). On this tongue, denoted $S31$, there is a 3:1 internal resonance between the in phase and out of phase NNMs. We note that other resonance scenarios (e.g., 2:1, 4:1, 5:1, etc.) exist in this seemingly simple system and are described in [13].

A third fundamental property of NNMs is that their number may exceed the number of DOFs of the system. Because of mode bifurcations, not all NNMs can be regarded as a nonlinear continuation of normal modes of linear systems. Internally resonant NNMs are one possible example. Moreover, NNMs can be stable or

<figure><img src="image_2.png" /><figcaption>Fig. 1 Frequency-energy plot of system (2). NNM motions depicted in the configuration space are inset. The horizontal and vertical axes in these plots are the displacements of the first and second DOFs, respectively; the aspect ratio is set so that increments on the horizontal and vertical axes are equal in size.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>Fig. 2 Frequency-energy plot of system (2) featuring a 3:1 internal resonance between the in-phase and out-of-phase NNMs.</figcaption></figure>

unstable, which is in contrast to linear theory, in which all modes are neutrally stable. In this context, instability means that small perturbations of the initial conditions that generate the NNM motion lead to the elimination of the mode oscillation. Therefore, unstable NNMs are not physically realizable. In the present study, the stability analysis is performed numerically through the eigenvalues of the monodromy matrix. Bifurcations and stability are interrelated concepts, because a change in stability necessarily occurs through a bifurcation.

NNMs have a fundamental limitation compared with their linear counterpart, that is, the general motion of a nonlinear system cannot be expressed as a superposition of individual NNM motions. The practical utility of NNMs might appear, at first, questionable. A first motivation to compute and exploit the NNMs is that forced resonances in nonlinear systems occur in their neighborhoods. The knowledge of the NNMs can therefore provide valuable insight into the structure of the resonances, a feature of considerable engineering importance [17]. A second motivation is that NNMs are attractors for the damped dynamics. For example, during the targeted energy transfer phenomenon [22], the NNMs attract the free dynamics, even though the motion is initiated using arbitrary initial conditions. A third motivation is that an experimental structure can be excited in one of its NNMs provided an appropriate exciter force pattern is applied. Nonlinear normal modes can therefore be isolated using force appropriation techniques, similar to what is achieved for linear structures (e.g., ground vibration testing of aircraft).

### C. Numerical Computation

There have been very few attempts to compute NNMs using numerical methods (see, for example, [21,23 25]). Interestingly, algorithms for the continuation of periodic solutions are really quite sophisticated and advanced [26,27], yet they have not been fully exploited for the computation of nonlinear modes.

The numerical method considered for the NNM computation relies on two main techniques, namely, a shooting technique and the pseudoarclength continuation method. A detailed description of the algorithm is given in [14].

#### 1. Shooting Method

The equations of motion of system (1) can be recast into state space form:

$$ \dot{\mathbf{z}} = \mathbf{g}(\mathbf{z}) \quad (3) $$

where $\mathbf{z} = [\mathbf{x}^* \ \dot{\mathbf{x}}^*]^*$ is the 2n dimensional state vector, with * denoting the transpose operation, and

$$ \mathbf{g}(\mathbf{z}) = \begin{pmatrix} \dot{\mathbf{x}} \\ -\mathbf{M}^{-1}[\mathbf{K}\mathbf{x} + \mathbf{f}_{nl}(\mathbf{x})] \end{pmatrix} \quad (4) $$

is the vector field. The solution of this dynamical system for initial conditions $\mathbf{z}(0) = \mathbf{z}_0 = [\mathbf{x}_0^* \ \dot{\mathbf{x}}_0^*]^*$ is written as $\mathbf{z}(t) = \mathbf{z}(t, \mathbf{z}_0)$ to exhibit the dependence on the initial conditions, $\mathbf{z}(0, \mathbf{z}_0) = \mathbf{z}_0$. A solution $\mathbf{z}_p(t, \mathbf{z}_{p0})$ is a periodic solution of the autonomous system (3) if $\mathbf{z}_p(t, \mathbf{z}_{p0}) = \mathbf{z}_p(t + T, \mathbf{z}_{p0})$, where T is the minimal period.

The NNM computation is carried out by finding the periodic solutions of the governing nonlinear equations of motion (3). In this context, the shooting method, probably one of the most popular numerical techniques, solves numerically the two point boundary value problem defined by the periodicity condition:

$$ \mathbf{H}(\mathbf{z}_{p0}, T) = \mathbf{z}_p(T, \mathbf{z}_{p0}) - \mathbf{z}_{p0} = \mathbf{0} \quad (5) $$

where $\mathbf{H}(\mathbf{z}_0, T) = \mathbf{z}(T, \mathbf{z}_0) - \mathbf{z}_0$ is called the shooting function and represents the difference between the initial conditions and the system response at time T. Unlike forced motion, the period T of the free response is not known a priori.

The shooting method consists of finding, in an iterative way, the initial conditions $\mathbf{z}_{p0}$ and the period T that realize a periodic motion. To this end, the method relies on direct numerical time integration and on the Newton Raphson algorithm.

Starting from some assumed initial conditions, $\mathbf{z}_{p0}^{(0)}$, the motion $\mathbf{z}_p^{(0)}(t, \mathbf{z}_{p0}^{(0)})$ at the assumed period $T^{(0)}$ can be obtained by numerical time integration methods (e.g., Runge Kutta or Newmark schemes). In general, the initial guess $(\mathbf{z}_{p0}^{(0)}, T^{(0)})$ does not satisfy the periodicity condition (5). A Newton Raphson iteration scheme is therefore used to correct the initial guess and to converge to the actual solution. The corrections $\Delta\mathbf{z}_{p0}^{(0)}$ and $\Delta T^{(0)}$ are found by expanding the nonlinear function

$$ \mathbf{H}(\mathbf{z}_{p0}^{(0)} + \Delta\mathbf{z}_{p0}^{(0)}, T^{(0)} + \Delta T^{(0)}) = 0 \quad (6) $$

in a Taylor series and neglecting higher order terms.

The phase of the periodic solutions is not fixed. If $\mathbf{z}(t)$ is a solution of the autonomous system (3), then $\mathbf{z}(t + \Delta t)$ is geometrically the same solution in phase space for any $\Delta t$. Hence, an additional condition, termed the phase condition, has to be specified to remove the arbitrariness of the initial conditions. This is discussed in detail in [14].

In summary, an isolated NNM is computed by solving the augmented two point boundary value problem defined by

$$ \mathbf{F}(\mathbf{z}_{p0}, T) = \begin{cases} \mathbf{H}(\mathbf{z}_{p0}, T) = 0 \\ h(\mathbf{z}_{p0}) = 0 \end{cases} \quad (7) $$

where $h(\mathbf{z}_{p0}) = 0$ is the phase condition.

#### 2. Continuation of Periodic Solutions

Because of the frequency energy dependence, the modal parameters of an NNM vary with the total energy. An NNM family, governed by Eqs. (7), therefore traces a curve, termed an NNM branch, in the $(2n + 1)$ dimensional space of initial conditions and period $(\mathbf{z}_{p0}, T)$. Starting from the corresponding LNM at low energy, the computation is carried out by finding successive points $(\mathbf{z}_{p0}, T)$ of the NNM branch using methods for the numerical continuation of periodic motions (also called path following methods) [27]. The space $(\mathbf{z}_{p0}, T)$ is termed the continuation space.

Among the different methods for numerical continuation found in the literature, the so called pseudoarclength continuation method is used herein.

Starting from a known solution $(\mathbf{z}_{p0,(j)}, T_{(j)})$, the method consists of computing the next periodic solution $(\mathbf{z}_{p0,(j+1)}, T_{(j+1)})$ on the branch using a predictor step and a corrector step.

1) Predictor step: At step j, a prediction $(\tilde{\mathbf{z}}_{p0,(j+1)}, \tilde{T}_{(j+1)})$ of the next solution $(\mathbf{z}_{p0,(j+1)}, T_{(j+1)})$ is generated along the tangent vector to the branch at the current point $\mathbf{z}_{p0,(j)}$:

$$ \begin{bmatrix} \tilde{\mathbf{z}}_{p0,(j+1)} \\ \tilde{T}_{(j+1)} \end{bmatrix} = \begin{bmatrix} \mathbf{z}_{p0,(j)} \\ T_{(j)} \end{bmatrix} + s_{(j)} \begin{bmatrix} \mathbf{p}_{z,(j)} \\ p_{T,(j)} \end{bmatrix} \quad (8) $$
where $s_{(j)}$ is the predictor step size. The tangent vector $\mathbf{p}_{(j)} = [\mathbf{p}_{z,(j)}^T \quad p_T_{(j)}]^T$ to the branch defined by Eqs. (7) is the solution of the system:

$$ \begin{bmatrix} \frac{\partial \mathbf{H}}{\partial \mathbf{z}_{p0}} |(\mathbf{z}_{p0,(j)}, T_{(j)}) & \frac{\partial \mathbf{H}}{\partial T} |(\mathbf{z}_{p0,(j)}, T_{(j)}) \\ \frac{\partial h^*}{\partial \mathbf{z}_{p0}} |(\mathbf{z}_{p0,(j)}) & 0 \end{bmatrix} \begin{bmatrix} \mathbf{p}_{z,(j)} \\ p_T_{(j)} \end{bmatrix} = \begin{bmatrix} \mathbf{0} \\ \mathbf{0} \end{bmatrix} \quad (9) $$

with the condition $||\mathbf{p}_{(j)}|| = 1$. This normalization can be taken into account by fixing one component of the tangent vector and solving the resulting overdetermined system using the Moore Penrose matrix inverse; the tangent vector is then normalized to 1.

2) Corrector step: The prediction is corrected by a shooting procedure to solve Eqs. (7), in which the variations of the initial conditions and the period are forced to be orthogonal to the predictor step. At iteration $k$, the corrections

$$ \mathbf{z}_{p0,(j+1)}^{(k+1)} = \mathbf{z}_{p0}^{(k)}_{(j+1)} + \Delta \mathbf{z}_{p0}^{(k)}_{(j+1)} \qquad T_{(j+1)}^{(k+1)} = T_{(j+1)}^{(k)} + \Delta T_{(j+1)}^{(k)} \quad (10) $$

are computed by solving the overdetermined linear system using the Moore Penrose matrix inverse:

$$ \begin{bmatrix} \frac{\partial \mathbf{H}}{\partial \mathbf{z}_{p0}} |(\mathbf{z}_{p0,(j+1)}^{(k)}, T_{(j+1)}^{(k)}) & \frac{\partial \mathbf{H}}{\partial T} |(\mathbf{z}_{p0,(j+1)}^{(k)}, T_{(j+1)}^{(k)}) \\ \frac{\partial h^*}{\partial \mathbf{z}_{p0}} |(\mathbf{z}_{p0,(j+1)}^{(k)}) & 0 \end{bmatrix} \begin{bmatrix} \Delta \mathbf{z}_{p0,(j+1)}^{(k)} \\ \Delta T_{(j+1)}^{(k)} \end{bmatrix} \mathbf{p}_{z,(j)}^T = \begin{bmatrix} -\mathbf{H}(\mathbf{z}_{p0,(j+1)}^{(k)}, T_{(j+1)}^{(k)}) \\ -h(\mathbf{z}_{p0,(j+1)}^{(k)}) \\ 0 \end{bmatrix} \quad (11) $$

where the prediction is used as initial guess, that is, $\mathbf{z}_{p0,(j+1)}^{(0)} = \tilde{\mathbf{z}}_{p0,(j+1)}$ and $T_{(j+1)}^{(0)} = \tilde{T}_{(j+1)}$. The last equation in Eq. (11) corresponds to the orthogonality condition for the corrector step.

This iterative process is carried out until convergence is achieved. The convergence test is based on the relative error of the periodicity condition:

$$ \frac{\|\mathbf{H}(\mathbf{z}_{p0}, T)\|}{\|\mathbf{z}_{p0}\|} = \frac{\|\mathbf{z}_p(T, \mathbf{z}_{p0}) - \mathbf{z}_{p0}\|}{\|\mathbf{z}_{p0}\|} < \epsilon \quad (12) $$

where $\epsilon$ is the prescribed relative precision.

### 3. Algorithm for NNM computation

The algorithm proposed for the computation of NNM motions is a combination of the shooting and pseudoarclength continuation methods, as shown in Fig. 3. It has been implemented in the MATLAB® environment. Other features of the algorithm, such as the step control, the reduction of the computational burden, and the method used for the numerical integration of the equations of motion, are discussed in [14].

Thus far, the NNMs have been considered as branches in the continuation space $(\mathbf{z}_{p0}, T)$. As explained in Sec. II.B, the NNMs are represented in an FEP in this study. This plot can be computed in a straightforward manner: 1) the conserved total energy is computed from the initial conditions realizing the NNM motion, and 2) the frequency of the NNM motion is calculated directly from the period.

## III. Modal Analysis of a Nonlinear Cyclic Structure

The objective of this section is to determine the essential features of the modal parameters of nonlinear periodic structures with cyclic symmetry. To this end, the modal shapes and frequencies of oscillation of the different NNM branches are computed from the free response of the undamped system.

The system considered throughout this paper is a simplified mathematical model of a bladed disk assembly. This model, though simplified, allows the investigation of interesting and complex

<figure><img src="image_4.png" /><figcaption>Fig. 3 Algorithm for NNM computation.</figcaption></figure>

dynamic phenomena related to the presence of nonlinearity. The lumped parameter model admits a single degree of freedom for each blade and includes a similarly simplified representation of the flexible disk. The bladed disk is composed of 30 sectors assembled with cyclic periodicity; a single sector is represented in Fig. 4. Each sector is modeled using disk (M) and blade (m) lumped masses coupled by linear (k) and cubic (knl) springs. The nonlinear springs can, for instance, be representative of geometrically nonlinear effects in slender blades. The disk masses are connected together by linear springs, K. The equations of motion of this 60 DOF system are

<figure><img src="image_5.png" /><figcaption>Fig. 4 One sector of the nonlinear bladed-disk assembly: a) continuous structure, and b) discrete model.</figcaption></figure>
<table><caption>Table 1 Natural frequencies of the underlying linear bladed assembly</caption><thead><tr><th>Mode</th><th>Nodal circles</th><th>Nodal diameters</th><th>Freq., rad/s</th><th>Mode</th><th>Nodal circles</th><th>Nodal diameters</th><th>Freq., rad/s</th></tr></thead><tbody><tr><td>1</td><td>0</td><td>0</td><td>0.000</td><td>31</td><td>1</td><td>0</td><td>2.082</td></tr><tr><td>2,3</td><td>0</td><td>1</td><td>0.183</td><td>32,33</td><td>1</td><td>1</td><td>2.084</td></tr><tr><td>4,5</td><td>0</td><td>2</td><td>0.363</td><td>34,35</td><td>1</td><td>2</td><td>2.092</td></tr><tr><td>6,7</td><td>0</td><td>3</td><td>0.536</td><td>36,37</td><td>1</td><td>3</td><td>2.104</td></tr><tr><td>8,9</td><td>0</td><td>4</td><td>0.700</td><td>38,39</td><td>1</td><td>4</td><td>2.123</td></tr><tr><td>10,11</td><td>0</td><td>5</td><td>0.850</td><td>40,41</td><td>1</td><td>5</td><td>2.147</td></tr><tr><td>12,13</td><td>0</td><td>6</td><td>0.985</td><td>42,43</td><td>1</td><td>6</td><td>2.178</td></tr><tr><td>14,15</td><td>0</td><td>7</td><td>1.103</td><td>44,45</td><td>1</td><td>7</td><td>2.215</td></tr><tr><td>16,17</td><td>0</td><td>8</td><td>1.202</td><td>46,47</td><td>1</td><td>8</td><td>2.258</td></tr><tr><td>18,19</td><td>0</td><td>9</td><td>1.282</td><td>48,49</td><td>1</td><td>9</td><td>2.304</td></tr><tr><td>20,21</td><td>0</td><td>10</td><td>1.346</td><td>50,51</td><td>1</td><td>10</td><td>2.350</td></tr><tr><td>22,23</td><td>0</td><td>11</td><td>1.394</td><td>52,53</td><td>1</td><td>11</td><td>2.394</td></tr><tr><td>24,25</td><td>0</td><td>12</td><td>1.428</td><td>54,55</td><td>1</td><td>12</td><td>2.431</td></tr><tr><td>26,27</td><td>0</td><td>13</td><td>1.452</td><td>56,57</td><td>1</td><td>13</td><td>2.460</td></tr><tr><td>28,29</td><td>0</td><td>14</td><td>1.465</td><td>58,59</td><td>1</td><td>14</td><td>2.478</td></tr><tr><td>30</td><td>0</td><td>15</td><td>1.470</td><td>60</td><td>1</td><td>15</td><td>2.485</td></tr></tbody></table>

$$m\ddot{x}_i + k(x_i - X_i) + k_{nl}(x_i - X_i)^3 = 0$$

$$M\ddot{X}_i + K(X_i - X_{i+1}) + K(X_i - X_{i-1}) + k(X_i - x_i) + k_{nl}(X_i - x_i)^3 = 0 \quad (13)$$

for $i = 1, \dots, 30$, where $X_{31} = X_1$ and $X_0 = X_{30}$ (conditions of cyclic periodicity). $X_i$ and $x_i$ are the displacements of the disk and blade masses of the $i$th sector, respectively. The values $M = 1$, $m = 0.3$, $K = 1$, $k = 1$, and $k_{nl} = 0.1$ are used in this study.

### A. Modal Analysis of the Underlying Linear System

Before studying the nonlinear bladed disk assembly, the natural frequencies and mode shapes of the underlying linear system are first discussed. All bladed assemblies with circumferential symmetry exhibit single and double modes [28].

Double modes represent the majority. They have the same natural frequency and similar mode shapes. In fact, no unique mode shapes can be specified for these modes. Rather, it is sufficient to specify two suitably orthogonal shapes and to note that, when vibrating freely at

that natural frequency, the structure can assume any form given by a linear combination of the two specified shapes. Further, at the corresponding natural frequency, the assembly can vibrate in any combination of $\cos n\theta$ and $\sin n\theta$ circumferential distributions of displacement around the assembly, that is, in a shape of the form $\cos n\theta + \phi$. The mode shape is characterized by $n$ nodal diameters because the displacement is constrained to be zero along $n$ equally spaced diametral lines. The mode shapes of a mode pair have mutually orthogonal nodal diameters.

Single modes correspond to motion with all the blades undergoing either in phase (zero nodal diameter) or out of phase ($N/2$ nodal diameters) motion.

The natural frequencies of the underlying linear bladed assembly are listed in Table 1, in which the modes are denoted by the integer pair $(n, m)$, which corresponds to the number of nodal circles $n$ and nodal diameters $m$ for the considered mode. In the model (13), the nodal circle parameter $n$ can only take the values of $n = 0$ or $n = 1$ according to whether the blade and disk masses undergo in phase or out of phase motion, respectively. One observes the existence of 28 pairs of double modes and 4 single modes. Figure 5 depicts four

<figure><img src="image_6.png" /><figcaption>Fig. 5 Representative LNM of the bladed assembly in which the blade and disk masses are shown in black and gray, respectively: a) mode (0,15) and one mode of the mode pair, b) (0,1), c) (1,5), and d) (1,14).</figcaption></figure>
representative LNMs of the bladed assembly, namely, mode (0,15) and one mode of the mode pairs (0,1), (1,5) and (1,14).

## B. Nonlinear Normal Modes

Modal analysis of the nonlinear bladed assembly is carried out in this section using the algorithm described earlier. NNM branches are computed by starting from the corresponding LNMAs at low energy and gradually increasing the total energy in the system. These branches, termed backbone branches, are represented in Fig. 6 and form the skeleton of the FEP. As we shall see, other NNM branches bifurcate from and coalesce into these backbone branches.

The first noticeable feature in Fig. 6 is the frequency energy dependence of the NNMs. The oscillation frequency of the modes with one nodal circle is strongly affected by the nonlinearities in the system. For these modes, the blade and disk masses vibrate in an out of phase fashion, which enhances nonlinear effects. On the other hand, the oscillation frequency of the modes with zero nodal circles is much less affected. This is because the blade and disk masses vibrate in an in phase fashion for these modes.

### 1. Similar and Nonsimilar Nonlinear Normal Modes

In addition to the dependence of their oscillation frequency, the NNMs may also have their modal shapes vary with the total energy in the system. According to Rosenberg's terminology [15], a similar NNM corresponds to an (energy independent) straight modal line in the configuration space and occurs only in systems presenting certain spatial symmetries. A nonsimilar NNM corresponds to a curve in the configuration space, the shape of which varies with the total energy. Because of its symmetry properties, the system possesses both similar and nonsimilar NNMs. Two examples of similar NNMs in the bladed disk are the nonlinear extension of the LNMAs with zero nodal diameters, namely, modes (0,0) and (1,0). Mode (0,0) is a rigid body mode, which is obviously unaffected by nonlinearity. The FEP of mode (1,0) in Fig. 7 clearly depicts that, although the NNM frequency is altered by the nonlinearities in the system, the modal shape remains unchanged.

Nonsimilar NNMs resemble the corresponding LNMAs at low energy. The structure (i.e., the number of nodal circles and diameters) is preserved and, as for the modes of the linear system, they mostly appear in pairs. Nonsimilar NNMAs in this system are either weakly, moderately, or strongly affected by nonlinearity for increasing energy levels. Figure 8 represents a mode of the mode pair (0,2), whose shape is almost energy independent. Figure 9 shows that the NNM motions of mode pair (0,14) have a marked energy dependence.

A remarkable property of the NNM motions of mode (1,14) is that the vibrational energy localizes to a limited number of sectors (four in

<figure><img src="image_7.png" /><figcaption>Fig. 7 FEP of mode (1,0). NNMAs represented by bar graphs are inset; they are given in terms of the initial displacements that realize the periodic motion (with zero initial velocities assumed). The blade and disk masses are shown in black and gray, respectively.</figcaption></figure>

this case), the rest of the system being virtually motionless (see Fig. 10). The resulting spatial confinement of the energy causes the responses of some blades to become dangerously high and might lead to premature failure of the blades. For illustration, the time series corresponding to such an NNM motion is displayed in Fig. 11. This localization phenomenon was also observed in linear mistuned bladed assemblies [29], but, here, it occurs even in the absence of structural disorder and direct interblade coupling. Localization is in fact the result of the frequency energy dependence inherent to nonlinear oscillations, as discussed in [8].

All these NNM motions correspond to standing wave motions in the sense that the system coordinates vibrate in a synchronous manner. They are represented by lines or curves in the configuration space. The phase condition used for their computation assumes that all initial velocities are zero. One therefore starts the motion from a maximum of the potential energy (see, for example, Fig. 11).

### 2. Modal Interaction: Internally Resonant Nonlinear Normal Modes

A first example of modal interaction is the 1:1 internal resonance between the two modes of a mode pair. This resonance scenario results in NNM motions, which take the form of traveling waves and are represented by ellipses in the configuration space. A detailed

<figure><img src="image_8.png" /><figcaption>Fig. 8 FEP of one mode of the mode pair (0,2). NNMAs represented by bar graphs are inset; they are given in terms of the initial displacements that realize the periodic motion (with zero initial velocities assumed). The blade and disk masses are shown in black and gray, respectively.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 6 Evolution of the NNM frequencies with the total energy in the system. The solid lines represent an NNM with zero nodal circles and the dashed lines represent an NNM with 1 nodal circle.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>Fig. 9 FEP of one mode of the mode pair (0,14). NNMs represented by bar graphs are inset; they are given in terms of the initial displacements that realize the periodic motion (with zero initial velocities assumed). The blade and disk masses are shown in black and gray, respectively.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 10 FEP of one mode of the mode pair (1,14). NNMs represented by bar graphs are inset; they are given in terms of the initial displacements that realize the periodic motion (with zero initial velocities assumed). The blade and disk masses are shown in black and gray, respectively.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 11 Time series corresponding to the localized NNM motion of mode (1,14) (see Fig. 10).</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Fig. 12 FEP of the traveling-wave NNM corresponding to a 1:1 internal resonance between the modes of mode pair (1,3) (solid line). For comparison, the dashed line represents the backbone of one standing-wave NNM of the mode pair (1,3). NNMs represented by bar graphs are inset.</figcaption></figure>

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

a)

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

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

b)

<caption>Fig. 13 Representative motion of the bladed assembly: a) time series of the first six blades during traveling-wave NNM motion of mode pair (1,3) (see Fig. 12), and b) contour plot.</caption>
<figure><img src="image_17.png" /><figcaption>Fig. 14 A 3:1 internal resonance between modes of the mode pairs (0,6) and (1,12). The solid line corresponds to the backbone of one mode of the mode pair (0,6), which is continued by a tongue of internally resonant NNM s. The dashed line corresponds to the backbone of one mode of the mode pair (1,12) represented at the third of its dominant frequency.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Fig. 15 Time series corresponding to NNM motions on the tongue of a 3:1 internal resonance (solid line: blade 1; dashed line: disk 10; dotted line: disk 14): a) beginning of the tongue (in the vicinity of the branch of mode (0,6)), b) middle of the tongue, and c) extremity of the tongue (in the vicinity of the branch of mode (1,12)).</figcaption></figure>
<figure><img src="image_19.png" /><figcaption>Fig. 16 Intricate succession of modal interactions between mode (0,6) and other modes of the system.</figcaption></figure>

analytical study of these modes is given in [6]. Because of the existence of a phase difference between the coordinates, a different phase condition is considered for the NNM computation: only one initial velocity is set to zero, which is compatible with a traveling wave motion. For instance, Fig. 12 depicts the NNMs corresponding to a 1:1 internal resonance between the modes of mode pair (1,3). The mode structure is preserved in the sense that this traveling wave

<figure><img src="image_20.png" /><figcaption>Fig. 17 Splitting of one mode of the mode pair (1,5). NNMs represented by bar graphs are inset; they are given in terms of the initial displacements that realize the periodic motion (with zero initial velocities assumed). The blade and disk masses are shown in black and gray, respectively.</figcaption></figure>

motion also features one nodal circle and three nodal diameters. Representative time series are shown in Fig. 13 and clearly highlight that the motion is no longer synchronous. In this particular case, the traveling wave is propagating in the counterclockwise direction, but its companion propagating in the clockwise direction also exists.

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

<figure><img src="image_22.png" /><figcaption>Fig. 18 Bifurcation of the mode pair (1,14).</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>Fig. 19 Bifurcation of mode (1,15).</figcaption></figure>

These modes have an important practical significance because they can be excited with an appropriate engine order excitation.

Other resonance scenarios can be observed in this system with the occurrence of tongues of internally resonant NNMs (see Sec. II.B). Unlike backbone branches, tongues are localized to a specific region of the FEP. They bifurcate from the backbone branch of a specific mode and coalesce into the backbone branch of another mode, thereby realizing an internal resonance between the two modes. For instance, Fig. 14 depicts a 3:1 internal resonance between modes (0,6) and (1,12) in the FEP. To better understand the resonance mechanism, the backbone of mode (1,12) is represented at the third of its characteristic frequency (this is relevant because a periodic solution of period *T* is also periodic with period 3T). This demonstrates that a smooth transition from mode (0,6) to mode (1,12) occurs on the tongue. A further illustration is that modes M1 and M2, which are the modes right after and before the coalescence of the two NNM branches, are almost identical.

During this 3:1 internal resonance, the system vibrates along a subharmonic NNM, that is, an NNM motion characterized by more than one dominant frequency component. On the branch of mode (0,6), the motion is characterized by one dominant frequency component, say, $ω$. As we move along the tongue from this branch, a third harmonic progressively appears, and the system vibrates with two dominant frequency components, $ω$ and $3ω$. As we progress further on the tongue, the third harmonic tends to dominate the component at the fundamental frequency until this latter completely disappears. At this precise moment, a transition to mode (1,12) is realized. This transition is illustrated in Fig. 15 using a time series representative of the NNM motion at three different locations on the tongue.

Surprisingly, the ratio of the linear natural frequencies of modes (0,6) and (1,12) is far from 3; it is equal to 2.47. A 3:1 internal resonance between the two modes can still be realized,

because the frequency of mode (0,6) increases much less rapidly than that of mode (1,12), as shown in Fig. 6. From this discussion, it turns out that a 3:1 internal resonance is not the only possible interaction between modes (0,6) and (1,12). Depending on the relative evolution of the frequencies on the backbones of these modes, other *n*:*m* resonances with *n* and *m* as relatively prime integers can exist. This highlights that NNMs can be internally resonant without necessarily having commensurate linear natural frequencies, a feature that is rarely discussed in the literature.

Another interesting finding is that there is a countable infinity of branches of internally resonant NNMs in this system, similar to what was reported for a 2 DOF system in [13]. Figure 16 depicts the same FEP as in Fig. 14, but the algorithm is not stopped after the tongue. Clearly, there is an intricate succession of modal interactions, each one a different realization of an internal resonance between mode (0,6) and another mode of the system.

### 3. Additional Mode Bifurcations

According to the degree of detuning for linear mistuned bladed assemblies, double modes with identical frequencies may split into two different modes with distinct natural frequencies [28]. Because of mode bifurcations, mode splitting may still occur in nonlinear systems with cyclic symmetry. A direct consequence is that a mode pair can bifurcate into two mode pairs, and the number of NNMs exceeds the number of DOFs of the system, as mentioned in Sec. II.B. For illustration, the splitting of one mode of the mode pair (1,5) is depicted in Fig. 17. Clearly, after the bifurcation, two NNM branches exist and are characterized by different oscillation frequencies and modal curves.

Other examples of mode bifurcations are shown in Figs. 18 and 19 for the mode pair (1,14) and mode (1,15), respectively. Mode pair
<figure><img src="image_24.png" /><figcaption>Fig. 20 Stability analysis of mode (1,15): a) FEP, in which a solid (dashed) line indicates stability (instability); and b) unit circle and Floquet multipliers, which are represented by squares.</figcaption></figure>

(1,14) undergoes a bifurcation during which a new NNM branch is generated. Interestingly, this branch is characterized by modes that are localized to two sectors only. Even more interesting is the bifurcation of mode (1,15), which generates mode shapes that are localized to only one sector. Depending on their stability properties, these NNMs may be excited in practice, giving rise to a potentially harmful motion that must be accounted for.

### 4. Mode Stability

In the present study, a stability analysis is performed numerically using the eigenvalues of the monodromy matrix (i.e., the Floquet multipliers), which are a by product of the proposed algorithm. NNM stability is an important concept because it dictates whether or not a mode is physically realizable. For instance, Fig. 20 shows the stability properties of mode (1,15). From very low energies to energies slightly above $1J$, the Floquet multipliers lie on the unit circle. The NNM motions are stable and, hence, physically realizable. From this latter energy, the Floquet multipliers leave the unit circle, and the NNM motions become unstable. This stability change occurs through a bifurcation, which coincides exactly with the generation of the branch of NNMs localized to one sector in Fig. 19.

## IV. Conclusions

In this study, the free vibrations of a nonlinear periodic structure with cyclic symmetry were examined using the NNM theory. The

NNMs were computed numerically by combining a shooting technique with pseudoarclength continuation. One advantage of this approach is that it provides the capability for analysis of strongly nonlinear regimes of motion.

A very complicated structure of NNMs, including similar and nonsimilar NNMs, nonlocalized and localized NNMs, and bifurcating and internally resonant NNMs, was observed. One important finding of this study is that modal interactions can occur without necessarily having commensurate natural frequencies in the underlying linear system. Furthermore, a countable infinity of such modal interactions were shown to exist in this system.

Further research will investigate how forced resonances of nonlinear periodic structures relate to the topological structure of their NNMs. This will provide useful insight into the forced dynamics, which is of particular practical significance.

## References

[1] Ewins, D. J., *Modal Testing: Theory, Practice and Application*, 2nd ed., Research Studies Press, Hertfordshire, England, U.K., 2000.

[2] Pierre, C., and Dowell, E. H., "Localization of Vibrations by Structural Irregularity," *Journal of Sound and Vibration*, Vol. 114, 1987, pp. 549-564.
doi:10.1016/S0022 460X(87)80023 8

[3] Bendiksen, O. O., "Localization Phenomena in Structural Dynamics," *Chaos, Solitons and Fractals*, Vol. 11, 2000, pp. 1621-1660.
doi:10.1016/S0960 0779(00)00013 8

[4] Castanier, M. P., and Pierre, C., "Using Intentional Mistuning in the
Design of Turbomachinery Rotors," *AIAA Journal*, Vol. 40, 2002, pp. 2077-2086.
doi:10.2514/2.1542

[5] Wei, S. T., and Pierre, C., "Effects of Dry Friction Damping on the Occurrence of Localized Forced Vibrations in Nearly Cyclic Structures," *Journal of Sound and Vibration*, Vol. 129, 1989, pp. 397-416.
doi:10.1016/0022-460X(89)90432 X

[6] Vakakis, A. F., "Dynamics of a Nonlinear Periodic Structure with Cyclic Symmetry," *Acta Mechanica*, Vol. 95, 1992, pp. 197-226.
doi:10.1007/BF01170813

[7] Vakakis, A. F., and Cetinkaya, C., "Mode Localization in a Class of Multidegree of Freedom Nonlinear Systems with Cyclic Symmetry," *SIAM Journal on Applied Mathematics*, Vol. 53, 1993, pp. 265-282.
doi:10.1137/0153016

[8] Vakakis, A. F., Nayfeh, T., and King, M. E., "A Multiple Scales Analysis of Nonlinear, Localized Modes in a Cyclic Periodic System," *Journal of Applied Mechanics*, Vol. 60, 1993, pp. 388-397.
doi:10.1115/1.2900806

[9] King, M. E., and Vakakis, A. F., "A Very Complicated Structure of Resonances in a Nonlinear System with Cyclic Symmetry: Nonlinear Forced Localization," *Nonlinear Dynamics*, Vol. 7, No. 1, Jan. 1995, pp. 85-104.
doi:10.1007/BF00045127

[10] King, M. E., and Layne, P. A., "Dynamics of Nonlinear Cyclic Systems with Structural Irregularity," *Nonlinear Dynamics*, Vol. 15, 1998, pp. 225-244.
doi:10.1023/A:1008291628528

[11] Folley, C. N., "Bifurcation and Symmetry in Cyclic Structures," Ph.D. Thesis, Purdue University, West Lafayette, IN, 1999.

[12] Samaranayake, S., Samaranayake, G., and Bajaj, A. K., "Resonant Vibrations in Harmonically Excited Weakly Coupled Mechanical Systems with Cyclic Symmetry," *Chaos, Solitons and Fractals*, Vol. 11, 2000, pp. 1519-1534.
doi:10.1016/S0960-0779(99)00075-2

[13] Kerschen, G., Peeters, M., Golinval, J. C., and Vakakis, A. F., "Nonlinear Normal Modes, Part I: A Useful Framework for the Structural Dynamicist," *Mechanical Systems and Signal Processing*, Vol. 23, 2009, pp. 170-194.
doi:10.1016/j.ymssp.2008.04.002

[14] Peeters, M., Viguié, R., Sérandour, G., Kerschen, G., and Golinval, J. C., "Nonlinear Normal Modes, Part II: Practical Computation Using Numerical Continuation Techniques," *Mechanical Systems and Signal Processing*, Vol. 23, 2009, pp. 195-216.
doi:10.1016/j.ymssp.2008.04.003

[15] Rosenberg, R. M., "On Nonlinear Vibrations of Systems with Many Degrees Of Freedom," *Advances in Applied Mechanics*, Vol. 9, Academic Press, New York/London/Orlando, FL, 1966, pp. 155-242.

[16] Vakakis, A. F., "Analysis and Identification of Linear and Nonlinear Normal Modes in Vibrating Systems," Ph.D. Dissertation, California

Inst. of Technology, Pasadena, CA, 1990.

[17] Vakakis, A. F., Manevitch, L. I., Mikhlin, Y. V., Pilipchuk, V. N., and Zevin, A. A., *Normal Modes and Localization in Nonlinear Systems*, Wiley, New York, 1996.

[18] Vakakis, A. F., "Non Linear Normal Modes and Their Applications in Vibration Theory: An Overview," *Mechanical Systems and Signal Processing*, Vol. 11, 1997, pp. 3-22.
doi:10.1006/mssp.1996.9999

[19] Shaw, S. W., and Pierre, C., "Non Linear Normal Modes and Invariant Manifolds," *Journal of Sound and Vibration*, Vol. 150, 1991, pp. 170-173.
doi:10.1016/0022-460X(91)90412 D

[20] Shaw, S. W., and Pierre, C., "Normal Modes for Non Linear Vibratory Systems," *Journal of Sound and Vibration*, Vol. 164, 1993, pp. 85-124.
doi:10.1006/jsvi.1993.1198

[21] Lee, Y. S., Kerschen, G., Vakakis, A. F., Panagopoulos, P. N., Bergman, L. A., and McFarland, D. M., "Complicated Dynamics of a Linear Oscillator with a Light, Essentially Nonlinear Attachment," *Physica D*, Vol. 204, 2005, pp. 41-69.
doi:10.1016/j.physd.2005.03.014

[22] Vakakis, A. F., Gendelman, O., Bergman, L. A., McFarland, D. M., Kerschen, G., and Lee, Y. S., *Nonlinear Targeted Energy Transfer in Mechanical and Structural Systems*, Solid Mechanics and Its Applications, Springer, Vol. 156, 2009.

[23] Slater, J. C., "A Numerical Method for Determining Nonlinear Normal Modes," *Nonlinear Dynamics*, Vol. 10, 1996, pp. 19-30.
doi:10.1007/BF00114796

[24] Pescheck, E., "Reduced Order Modeling of Nonlinear Structural Systems Using Nonlinear Normal Modes and Invariant Manifolds," Ph.D. Thesis, Univ. of Michigan, Ann Arbor, MI, 2000.

[25] Arquier, R., Bellizzi, S., Bouc, R., and Cochelin, B., "Two Methods for the Computation of Nonlinear Modes of Vibrating Systems at Large Amplitudes," *Computers and Structures*, Vol. 84, 2006, pp. 1565-1576.
doi:10.1016/j.compstruc.2006.01.011

[26] Doedel, E., AUTO, Software for Continuation and Bifurcation Problems in Ordinary Differential Equations, 2007.

[27] Seydel, R., *Practical Bifurcation and Stability Analysis, from Equilibrium to Chaos*, 2nd ed., Springer Verlag, Berlin, 1994.

[28] Ewins, D. J., "Structural Dynamics Characteristics of Bladed Assemblies," Lecture Series 1992 06 on Vibration and Rotor Dynamics, Von Karman Institute for Fluid Dynamics, 1992.

[29] Castanier, M. P., and Pierre, C., "Modeling and Analysis of Mistuned Bladed Disk Vibration: Status and Emerging Directions," *Journal of Propulsion and Power*, Vol. 22, 2006, pp. 384-396.
doi:10.2514/1.16345